A matrix-free, differentiable PyTorch solver for phase-field fracture: Formulation, benchmarks, and inverse analysis

2026-06-22

AI Summary

This work addresses the limitations of traditional phase-field fracture simulations, which rely on global sparse matrix assembly and thus hinder efficient GPU utilization and compatibility with automatic differentiation, impeding integration with optimization and machine learning. The authors propose the first purely PyTorch-based, matrix-free, and differentiable phase-field fracture solver, formulated along an explicit dynamics pathway. By recasting finite element operations as tensor contractions and scatter-add accumulations, the method avoids assembling the global stiffness matrix entirely. It supports multiple energy decompositions (e.g., spectral) and fracture models (AT1/AT2), operates seamlessly on both CPU and GPU, and exhibits memory consumption independent of conjugate gradient iterations. Custom backward propagation enables implicit differentiation through implicitly solved damage fields. The solver accurately reproduces six canonical fracture benchmarks on meshes with up to one million nodes and achieves end-to-end inversion of the critical fracture energy \( G_c \) with relative errors below \( 10^{-3} \).

Abstract

A matrix-free, open-source PyTorch solver is presented for phase-field fracture on central processing units (CPUs) and graphics processing units (GPUs) without custom compiled extensions. In the explicit dynamic pathway, finite-element operations are formulated as element-wise tensor contractions with scatter-based accumulation, removing global sparse mechanics-stiffness assembly from the core time-stepping loop. Both Ambrosio-Tortorelli regularisations (AT1 and AT2), multiple energy decompositions (spectral, volumetric-deviatoric, and star-convex), and plane strain or plane stress assumptions are supported. The explicit mechanics kernels are compatible with PyTorch's automatic differentiation engine (autograd), while the implicit, bound-constrained damage solve is wrapped in a custom backward rule. This rule implements implicit differentiation through the conjugate-gradient (CG) linear solve and keeps memory independent of the internal CG iteration count. The same implementation runs unmodified across macOS, Linux, and Windows, and has been run on meshes of order $10^6$ nodes on a single NVIDIA A100 GPU. The solver is compared against four dynamic fracture cases (straight crack propagation, shear-induced kinking, dynamic branching, and crack-hole interaction in perforated plates) and two quasi-static cases (single-edge notched tension and a notched-holed plate). As a differentiability demonstration, the scalar fracture energy $G_c$ is recovered from observed crack patterns using PyTorch gradients through the forward solve and limited-memory Broyden-Fletcher-Goldfarb-Shanno (L-BFGS) optimisation. Recovery of $G_c$ with relative error below $10^{-3}$ is achieved after three accepted L-BFGS states for glass and two for alumina. The implementation can be extended and combined with differentiable optimisation and machine-learning components.

Problem

• phase-field fracture
• matrix-free solver
• differentiable simulation
• inverse analysis
• PyTorch

Innovation

• matrix-free
• differentiable solver
• phase-field fracture
• implicit differentiation
• PyTorch autograd