Source-linked AI summary
A simple and robust Abaqus implementation of the phase field fracture method
Yousef Navidtehrani, Covadonga Betegón, Emilio Martínez-Pañeda
TL;DR
The paper addresses Abaqus implementation challenges for phase-field fracture, where ad hoc finite-element programming may be required. It presents an implementation that avoids user-defined-element coding, reaches convergence across cases, and favors single-pass staggered efficiency for demanding problems.
Problem
Phase-field fracture implementations in Abaqus can require programming an ad hoc finite element.
Method
The paper presents a phase-field fracture implementation in Abaqus that avoids coding user-defined elements and exploits built-in finite-element features.
Results
Standard Newton converges in all cases, while a single-pass staggered scheme appears more efficient for convergence-demanding problems.
Takeaways & Limitations
The implementation combines Abaqus’s built-in features with convergence across the reported cases and improved efficiency from single-pass staggering.
Takeaways & Limitations
The solution scheme is not unconditionally stable, and results can be sensitive to the number of time increments.
Abstract
from arXiv · showhide
The phase field fracture method is attracting significant interest. Phase field approaches have enabled predicting - on arbitrary geometries and dimensions - complex fracture phenomena such as crack branching, coalescence, deflection and nucleation. In this work, we present a simple and robust implementation of the phase field fracture method in the commercial finite element package Abaqus. The implementation exploits the analogy between the phase field evolution law and the heat transfer equation, enabling the use of Abaqus' in-built features and circumventing the need for defining user elements. The framework is general, and is shown to accommodate different solution schemes (staggered and monolithic), as well as various constitutive choices for preventing damage under compression. The robustness and applicability of the numerical framework presented is demonstrated by addressing several 2D and 3D boundary value problems of particular interest. Focus is on the solution of paradigmatic case studies that are known to be particularly demanding from a convergence perspective. The results reveal that our phase field fracture implementation can be readily combined with other advanced computational features, such as contact, and deliver robust and precise solutions. The code developed can be downloaded from www.empaneda.com/codes.
1. Introduction
Phase-field fracture methods avoid explicit crack-interface tracking while modelling complex fracture on arbitrary geometries and dimensions. This paper presents a robust Abaqus implementation that uses the heat-transfer analogy, supports multiple solution schemes, and is demonstrated on demanding 2D and 3D problems.
- Motivation: Explicitly tracking moving fracture boundaries requires interfacial boundary conditions and manual topology adjustments, creating implementation challenges.These difficulties motivate diffuse-interface formulations.
- Motivation: Phase-field fracture models complex crack trajectories and inertia-driven branching on the original finite-element mesh without arbitrary propagation criteria.The approach applies to arbitrary geometries and dimensions.
- Motivation: In Abaqus, solving the phase field as an additional degree of freedom traditionally requires user-element subroutines, limiting access to built-in features.Existing implementations can effectively use Abaqus only as a solver and lose much of its functionality.
- Contribution: The proposed implementation exploits the analogy between phase-field evolution and heat conduction to use Abaqus features and avoid coding user-defined elements.Coupled temperature-displacement elements support the implementation and reduce meshing and visualization complications.
2. Phase field fracture model
The phase-field fracture model regularizes Griffith fracture by replacing explicit crack tracking with a damage-like phase field coupled to displacement. Its governing equations can be solved using either monolithic or staggered minimization schemes.
- Variational formulation: The model minimizes elastic and fracture energies using displacement u and phase field φ as primary kinematic variables.The phase field is solved at finite-element nodes as an additional degree of freedom.
- Solution schemes: The discretized equations support monolithic solution, which solves u and φ simultaneously, and staggered solution, which alternates minimization steps.These are the two principal numerical solution strategies described for the model.
- Phase-field regularization: The phase field φ tracks the crack interface, varying from 0 in intact material to 1 inside the crack.A degradation function g(φ) = (1 −φ)2 reduces material stiffness as damage evolves.
- Phase-field regularization: The regularized formulation introduces a length scale ℓ governing the fracture process-zone size and, for finite ℓ, the material strength.Finite ℓ can capture crack nucleation and its transition from strength-driven to fracture-driven behavior.
- Governing equations: The coupled field equations follow from variations of the regularized energy with respect to displacement and phase field under equilibrium assumptions.The formulation restricts the analysis to linear elastic materials.
3. Finite element implementation
The implementation reformulates phase-field evolution as heat transfer, using Abaqus subroutines and displacement-temperature elements to solve coupled fracture problems. It supports monolithic and staggered schemes, trading unconditional stability against convergence robustness and efficiency.
- Damage irreversibility: A history field enforces damage irreversibility by preventing crack healing during loading and unloading.The field satisfies Kuhn–Tucker conditions and uses current and previous phase-field values.
- Heat-transfer analogy: The phase-field equation is mapped to a heat-transfer problem by identifying temperature with the phase field and defining an internal heat source.Under steady-state conditions, the rate term vanishes; transient heat capacity can instead provide viscous regularisation and stabilisation.
- Abaqus implementation: The heat-transfer analogy enables integration-point implementation with built-in displacement-temperature elements, avoiding user elements and additional preprocessing or postprocessing.The phase field is assigned as temperature and visualised through the nodal solution temperature NT11.
- Abaqus implementation: The Abaqus implementation combines UMAT and HETVAL subroutines, transferring history variables through solution-dependent state variables without external Fortran modules.UMAT computes degraded stresses and energy quantities, while HETVAL defines the internal heat flux and its derivative.
- Solution schemes: Monolithic and staggered schemes provide alternative solution strategies: monolithic solves displacement and phase field simultaneously, whereas staggered solves their subsystems sequentially.The staggered approach freezes the history field during iterations, facilitating convergence in demanding problems but sacrificing unconditional stability.
- Solution schemes: Abaqus lacks quasi-Newton thermo-mechanical solving, so the framework implements conventional Newton-based monolithic and single-pass staggered schemes.The framework is intended to balance robustness and efficiency according to the problem at hand.
4. Results
The implementation is evaluated on demanding two- and three-dimensional fracture problems, including unstable and stable crack growth, screw failure, and a contact-inclusive Brazilian test. Across these cases, it captures expected crack patterns and generally achieves convergence, while solution-scheme efficiency depends on the boundary-value problem.
- Test cases: The study evaluates fracture in notched plates under tension and shear, screws with and without initial cracks, and a 3D Brazilian test with contact.These cases address crack initiation, growth, nucleation, coalescence, and convergence in representative boundary-value problems.
- Notched square plate under shear: The shear case produces the expected deflected crack path and stable force reduction, but the monolithic implementation requires many iterations and is less efficient than the staggered scheme.The staggered implementation converges after two increments and predicts a similar maximum force, while remaining sensitive to the number of time increments.
- Notched square plate under shear: For the shear plate, the staggered implementation appears more robust and efficient than the monolithic one, although both are outperformed by the referenced quasi-Newton approach.The comparison is specific to this case study, and the staggered scheme is not unconditionally stable.
- Screw tension tests: Screw simulations converge readily and reproduce qualitatively reported crack trajectories, including nucleation near the head without an initial defect and branching from an initial crack.The defect-free sample carries a larger force, whereas the case with a smaller crack reaches a similar maximum-force magnitude.
- 3D Brazilian test: The implementation also resolves a 3D Brazilian test with contact, predicting a central crack that propagates toward the jaws and smaller cracks near the loading region without observed monolithic convergence issues.The Brazilian-test passages identify the test and report the observed crack pattern and convergence behavior.
5. Conclusions
The paper presents a simple Abaqus implementation of phase field fracture that uses built-in features, supports multiple schemes and constitutive choices, and handles demanding 2D and 3D problems. Both monolithic and staggered approaches converge, while the staggered scheme is more efficient for convergence-demanding cases and neither scheme matches quasi-Newton efficiency.
- Implementation: The implementation uses the analogy between heat transfer and phase field evolution to avoid user-defined elements while exploiting Abaqus features.It combines a UMAT and HETVAL subroutine and can be developed within Abaqus’ graphical user interface.
- Solution schemes: The monolithic Newton implementation reaches convergence in all examined cases, while a single-pass staggered scheme appears more efficient for convergence-demanding problems.The conclusions compare monolithic and staggered schemes across the reported boundary-value problems.
- Efficiency: Both schemes appear less efficient than quasi-Newton methods, and interpolation schemes might not improve phase field fracture efficiency.These observations qualify the computational efficiency of the provided implementations.
- Extensions: The framework can be extended to other material models and damage mechanisms, including plasticity and fatigue.The stated extension concerns the framework’s constitutive and damage-model scope.
Appendix A. Additional details of numerical implementation
The implementation can incorporate alternative constitutive choices and anisotropic phase field formulations through the UMAT, but the accompanying code omits these additional features for simplicity.
- Extensions: The framework can be extended with other constitutive choices and an anisotropic phase field formulation.The extensions include alternative strain-energy splits and anisotropic treatment in the governing formulation.
- Compression treatment: A tension-compression split of the fracture driving force should prevent damage from developing under compressive stresses.The text identifies this split as a constitutive choice for suppressing compression-driven damage.
- Implementation: The additional constitutive and anisotropic extensions are implemented in the UMAT subroutine.The implementation location is specified for these extensions.
- Scope boundary: The accompanying code does not include the additional features because it is provided in a simplified form.The omission is explicitly stated as a limitation of the distributed code.
Appendix A.1. Strain energy density decomposition
The appendix describes volumetric-deviatoric and spectral tension-compression decompositions of strain energy, with only the tensile contribution entering the history field.
- Decomposition choices: The paper considers the Miehe spectral tension-compression split and the Amor volumetric-deviatoric split.These are presented as the two most widely used strain-energy decompositions.
- History field: Both decompositions separate the strain energy density into tensile and compressive parts, but only ψ+ contributes to the history field.The history field is evaluated using the tensile energy contribution.
- Volumetric-deviatoric split: In the volumetric-deviatoric split, compressive volumetric strain energy does not contribute to fracture.The split retains the tensile volumetric contribution and deviatoric contribution while excluding the compressive volumetric part.
- Notation: The appendix defines the deviatoric strain from the strain tensor and its trace, and uses Macaulay brackets to separate positive and negative parts.The definitions specify ε′ and the positive and negative bracket operators.
- Spectral split: The spectral decomposition uses principal strains and directions to obtain tensile and compressive strain components.The tensile and compressive components are formed from the corresponding principal strains.
Appendix A.2. Anisotropic formulation
The anisotropic formulation extends the implementation by applying tension-compression decomposition in the displacement field equation, with the associated constitutive treatment handled in the implementation.
- Anisotropic extension: The implementation is extended with the anisotropic phase field approach using a hybrid formulation.This extension is attributed to the hybrid approach proposed by Ambati et al.
- Displacement formulation: The anisotropic formulation applies the tension-compression decomposition in the displacement field equation as well as the phase field treatment.The Cauchy stress is consequently modified in the formulation.
- Projection operators: The appendix introduces projection tensors for the tensile and compressive components, with P− defined from the fourth-order identity and P+.The text states P− = I − P+ and identifies I as the fourth-order identity tensor.
- Tensor construction: The projection tensor construction uses principal strain directions and Heaviside-function definitions.The appendix gives the principal-direction components and the Heaviside convention used in the formulation.
- Special case: The formulation encounters a special case in which the cited tensor expression cannot be evaluated.The text notes this issue without specifying a broader failure boundary.