# What Solvers Are Available in Phasefieldx and How Do They Differ: A Complete Guide

> Explore Phasefieldx's nine finite-element solvers, from staggered fracture to monolithic formulations. Understand their differences and find the best fit for your material models and coupling strategies.

- Repository: [Miguel Castillón/phasefieldx](https://github.com/castillonmiguel/phasefieldx)
- Tags: deep-dive
- Published: 2026-02-26

---

**Phasefieldx provides nine distinct finite-element solvers ranging from staggered phase-field fracture to monolithic energy-controlled formulations, each optimized for specific coupling strategies, constraint enforcement, and material models.**

The `castillonmiguel/phasefieldx` repository implements a modular family of finite-element solvers built on FEniCSx/dolfinx. While all solvers share common infrastructure for MPI communication, logging, and Paraview output, they differ fundamentally in how they couple physics, enforce constraints, and handle time integration. Understanding these distinctions is essential for selecting the right tool for fracture, elasticity, or phase-field diffusion problems.

## Phase-Field Fracture Solvers

The fracture solvers reside in `src/phasefieldx/Element/Phase_Field_Fracture/solver/` and represent the core capability of the library. They differ primarily in their coupling strategy between displacement `u` and phase-field `Φ`.

### Staggered Solver (solver.py)

The **staggered solver** ([`solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver.py)) treats the displacement and phase-field equations sequentially within each load step. First, it solves the nonlinear displacement problem using a Newton-based method with the current phase-field value. Then, it projects the elastic energy (`ψ_a`) and solves the phase-field equation.

This approach uses adaptive pseudo-time stepping (`dtau`) that increases when Newton convergence requires two or fewer iterations and decreases otherwise.

```python
from phasefieldx.Element.Phase_Field_Fracture.solver import solve

solve(
    Data, msh, final_time=1.0,
    V_u=V_u, V_Φ=V_Φ,
    bc_list_u=bc_u, bc_list_phi=bc_phi,
    f_list_u=None, T_list_u=tractions,
    update_boundary_conditions=update_bc,
    update_loading=update_load,
    dt=0.01, path="results/fracture"
)

```

*Source*: [[`solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver.py)](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Phase_Field_Fracture/solver/solver.py#L43-L61)

### Energy-Controlled Variational Solver (solver_ener_variational.py)

The **variational energy-controlled solver** ([`solver_ener_variational.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_ener_variational.py)) employs a monolithic (blocked) approach that simultaneously solves for displacement `u`, phase-field `Φ`, and a Lagrange multiplier `λ`. The multiplier enforces the constraint `c₁·γ + c₂·E_ext = τ`, where `γ` is the crack surface energy and `E_ext` is external work.

This solver uses `scifem.BlockedNewtonSolver` with an analytically assembled Jacobian for the three-field system `[u, Φ, λ]`. The variational formulation includes the derivative of the constraint with respect to `Φ`, ensuring rigorous energy control.

```python
from phasefieldx.Element.Phase_Field_Fracture.solver.solver_ener_variational import solve

solve(
    Data, msh, final_gamma=0.5,
    V_u=V_u, V_Φ=V_Φ,
    bc_list_u=bc_u, bc_list_phi=bc_phi,
    T_list_u=tractions,
    dtau=5e-4, path="results/ener_variational"
)

```

*Source*: [[`solver_ener_variational.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_ener_variational.py)](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Phase_Field_Fracture/solver/solver_ener_variational.py#L41-L58)

### Energy-Controlled Non-Variational Solver (solver_ener_non_variational.py)

The **non-variational energy-controlled solver** ([`solver_ener_non_variational.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_ener_non_variational.py)) implements the same energy constraint as the variational version but omits the derivative term `∂F_Φ/∂λ` from the Jacobian. This simplification yields a sparser block structure and faster assembly at the cost of weaker constraint enforcement.

Like its variational counterpart, it uses `scifem.BlockedNewtonSolver` but constructs a different `J_blocked` matrix that lacks the coupling block between the phase-field residual and the Lagrange multiplier.

```python
from phasefieldx.Element.Phase_Field_Fracture.solver.solver_ener_non_variational import solve

solve(
    Data, msh, final_gamma=0.5,
    V_u=V_u, V_Φ=V_Φ,
    bc_list_u=bc_u, bc_list_phi=bc_phi,
    T_list_u=tractions,
    dtau=5e-4, path="results/ener_nonvariational"
)

```

*Source*: [[`solver_ener_non_variational.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_ener_non_variational.py)](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Phase_Field_Fracture/solver/solver_ener_non_variational.py#L40-L58)

### Anisotropic Elasticity Solver (solver_anisotropic.py)

The **anisotropic solver** ([`solver_anisotropic.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_anisotropic.py)) handles elasticity problems requiring stiffness decomposition into tensile and compressive components (`ψ_a`, `ψ_b`, `σ_a`, `σ_b`). While primarily used as a sub-problem within fracture simulations to degrade only tensile energy, this solver functions independently for stand-alone anisotropic elasticity analyses.

It employs a standard `dolfinx.fem.petsc.NonlinearProblem` Newton solver without phase-field coupling.

```python
from phasefieldx.Element.Phase_Field_Fracture.solver.solver_anisotropic import solve

solve(
    Data, msh, final_time=1.0,
    V_u=V_u,
    bc_list_u=bc_u,
    T_list_u=tractions,
    dt=0.02, path="results/anisotropic"
)

```

*Source*: [[`solver_anisotropic.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_anisotropic.py)](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Phase_Field_Fracture/solver/solver_anisotropic.py#L27-L41)

## Stand-Alone Physics Solvers

Beyond fracture mechanics, phasefieldx provides specialized solvers for pure phase-field evolution, linear elasticity, and Allen-Cahn diffusion.

### Pure Phase-Field Solver (Element/Phase_Field/solver/solver.py)

Located in [`src/phasefieldx/Element/Phase_Field/solver/solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Phase_Field/solver/solver.py), this solver integrates the scalar phase-field equation without mechanical coupling. It uses a simple Newton iteration on a single scalar field, making it ideal for benchmarking phase-field evolution or testing new degradation functions in isolation.

### Linear Elasticity Solver (Element/Elasticity/solver/solver.py)

The isotropic elasticity solver ([`src/phasefieldx/Element/Elasticity/solver/solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Elasticity/solver/solver.py)) solves standard linear elasticity problems with a single displacement field `u`. Unlike the anisotropic variant, it does not split strain energy into tensile and compressive parts, using the standard isotropic constitutive law (`epsilon`). This solver serves for verification against analytical solutions and classical elasticity problems.

```python
from phasefieldx.Element.Elasticity.solver.solver import solve

solve(
    Data, msh, final_time=1.0,
    V_u=V_u,
    bc_list_u=bc_u,
    T_list_u=tractions,
    dt=0.01, path="results/elastic"
)

```

*Source*: [[`Elasticity/solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/Elasticity/solver.py)](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Elasticity/solver/solver.py#L27-L41)

### Allen-Cahn Solvers

The Allen-Cahn module provides two distinct approaches to phase-separation kinetics.

**Static Solver (solver_static.py):** Located in [`src/phasefieldx/Element/Allen_Cahn/solver/solver_static.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Allen_Cahn/solver/solver_static.py), this solver computes equilibrium configurations by solving the steady Allen-Cahn equation (`∂ψ/∂c = 0`) as a linear system. It is appropriate for computing equilibrium microstructures without temporal evolution.

```python
from phasefieldx.Element.Allen_Cahn.solver.solver_static import solve

solve(
    Data, msh, V_c=V_c,
    bc_list_c=bc_c,
    path="results/allen_cahn_static"
)

```

*Source*: [[`solver_static.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_static.py)](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Allen_Cahn/solver/solver_static.py#L27-L41)

**Dynamic Solver (solver.py):** The dynamic version ([`src/phasefieldx/Element/Allen_Cahn/solver/solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Allen_Cahn/solver/solver.py)) implements time-dependent Allen-Cahn kinetics using an implicit Euler scheme. It performs a linear solve at each time step to evolve microstructures temporally.

```python
from phasefieldx.Element.Allen_Cahn.solver.solver import solve

solve(
    Data, msh, V_c=V_c,
    bc_list_c=bc_c,
    dt=0.005, final_time=2.0,
    path="results/allen_cahn_dynamic"
)

```

*Source*: [[`Allen_Cahn/solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/Allen_Cahn/solver.py)](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Allen_Cahn/solver/solver.py#L27-L43)

## Key Differences in Coupling Strategies

Understanding the numerical architecture of these solvers requires examining three critical distinctions: coupling strategy, constraint enforcement, and time integration.

### Staggered vs. Monolithic Coupling

**Staggered solvers** (e.g., [`solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver.py)) solve the displacement field `u` and phase-field `Φ` sequentially within each load step. This approach decouples the nonlinear systems, making each solve smaller and more robust for large-scale problems, though it may require sub-stepping for strong coupling.

**Monolithic (blocked) solvers** (e.g., [`solver_ener_variational.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_ener_variational.py)) assemble the displacement, phase-field, and Lagrange multiplier into a single system matrix. This tight coupling enforces constraints exactly at the residual level but requires assembling the full Jacobian including cross-derivatives like `∂F_u/∂Φ`.

### Variational vs. Non-Variational Constraint Handling

The energy-controlled fracture solvers differ in how they enforce the crack surface constraint `c₁·γ + c₂·E_ext = τ`.

**Variational formulations** include the derivative of the constraint with respect to the phase-field `Φ` in the residual and Jacobian (`∂F_Φ/∂λ`). This guarantees that the energy constraint is satisfied at the variational level, maintaining consistency with the underlying physics.

**Non-variational formulations** omit the `∂F_Φ/∂λ` term from the Jacobian, simplifying the block structure. While this accelerates assembly, the constraint is only approximately satisfied, making it suitable for problems where exact energy control is less critical than computational speed.

### Time-Stepping and Adaptivity

Fracture solvers in `phasefieldx` implement **adaptive pseudo-time stepping** (`dtau`). The algorithm increases the load step when Newton convergence requires two or fewer iterations and decreases it when convergence fails, ensuring robust crack propagation tracking.

In contrast, the **Allen-Cahn dynamic solver** uses a user-specified fixed time step `dt` with implicit Euler integration, while the **elasticity solvers** typically perform single static solves where `final_time` serves merely as a pseudo-time parameter for loading ramps.

## Summary

- **Staggered fracture solvers** ([`solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver.py)) sequentially solve displacement and phase-field equations with adaptive time stepping, ideal for general fracture simulations.
- **Monolithic energy-controlled solvers** ([`solver_ener_variational.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_ener_variational.py) and [`solver_ener_non_variational.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_ener_non_variational.py)) use blocked Newton methods with Lagrange multipliers to enforce crack surface constraints, differing in whether they include the full variational derivative.
- **Anisotropic elasticity solver** ([`solver_anisotropic.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_anisotropic.py)) handles tensile/compressive energy splits independently or as a sub-problem within fracture models.
- **Stand-alone physics solvers** cover pure phase-field evolution ([`Element/Phase_Field/solver/solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/Element/Phase_Field/solver/solver.py)), isotropic linear elasticity ([`Element/Elasticity/solver/solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/Element/Elasticity/solver/solver.py)), and both static and dynamic Allen-Cahn diffusion (`Element/Allen_Cahn/solver/`).

## Frequently Asked Questions

### What is the difference between staggered and monolithic solvers in phasefieldx?

**Staggered solvers** solve the displacement and phase-field equations sequentially within each load step, making each individual solve smaller and more memory-efficient. **Monolithic solvers** assemble both fields (plus any Lagrange multipliers) into a single blocked system matrix, enforcing tighter coupling at the cost of larger Jacobian assemblies. Choose staggered for robustness in large problems; choose monolithic when exact constraint enforcement (like energy control) is required.

### When should I use the variational energy-controlled solver versus the non-variational version?

Use the **variational solver** ([`solver_ener_variational.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_ener_variational.py)) when you need rigorous enforcement of the crack surface energy constraint `c₁·γ + c₂·E_ext = τ` at the variational level, as it includes the full derivative `∂F_Φ/∂λ` in the Jacobian. Use the **non-variational solver** ([`solver_ener_non_variational.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_ener_non_variational.py)) when computational speed is prioritized over exact constraint satisfaction, as the simplified Jacobian (lacking `∂F_Φ/∂λ`) assembles faster but enforces the constraint only approximately.

### Can I use the anisotropic elasticity solver without phase-field fracture?

Yes. While [`solver_anisotropic.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_anisotropic.py) resides in the `Phase_Field_Fracture` module and implements the tensile/compressive energy split (`ψ_a`, `ψ_b`) required for fracture models, it functions as a **stand-alone anisotropic elasticity solver**. You can invoke it for pure mechanical problems requiring stiffness decomposition without activating any phase-field evolution, making it suitable for analyzing materials with distinct tensile and compressive responses.

### Which solver is best for pure phase-field diffusion without mechanics?

For pure phase-field evolution without displacement coupling, use the **Allen-Cahn dynamic solver** ([`src/phasefieldx/Element/Allen_Cahn/solver/solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Allen_Cahn/solver/solver.py)) for time-dependent microstructure evolution or the **Allen-Cahn static solver** ([`solver_static.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/solver_static.py)) for equilibrium configurations. If you specifically need the phase-field equation used in fracture models (e.g., for testing degradation functions), use the **pure phase-field solver** ([`src/phasefieldx/Element/Phase_Field/solver/solver.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/Element/Phase_Field/solver/solver.py)), which solves only the scalar phase-field equation without mechanical coupling.