# How to Calculate Error Norms and Check Convergence in PhaseFieldX

> Learn how to calculate L², H¹, and Lᵖ error norms and monitor convergence in PhaseFieldX. Utilize dedicated modules for accurate error analysis and solver feedback.

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

---

**PhaseFieldX provides dedicated modules in [`src/phasefieldx/norms.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/norms.py) and [`src/phasefieldx/errors_functions.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/errors_functions.py) to compute L², H¹, and Lᵖ error norms, while convergence monitoring is handled by the Newton solver in [`src/phasefieldx/solvers/newton.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/solvers/newton.py) using PETSc SNES residual norms.**

PhaseFieldX is an open-source finite element framework built on DOLFINx for simulating phase-field fracture and multi-physics problems. When validating numerical solutions against analytical references or debugging solver behavior, you need rigorous tools to calculate error norms and check convergence of nonlinear iterations.

## Computing Error Norms in PhaseFieldX

### Core Norm Implementations in norms.py

The foundation for all error calculations resides in [`src/phasefieldx/norms.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/norms.py). These functions implement standard functional norms using DOLFINx forms with MPI-aware reduction, returning a single global `float` across all ranks.

- `norm_L2(field, msh, dx=ufl.dx)` computes $\sqrt{\int \|field\|^2 \, dx}$
- `norm_Lp(field, msh, p=2, dx=ufl.dx)` computes $(\int \|field\|^p \, dx)^{1/p}$  
- `norm_semiH1(field, msh, dx=ufl.dx)` computes $\sqrt{\int \|\nabla field\|^2 \, dx}$
- `norm_H1(field, msh, dx=ufl.dx)` computes $\sqrt{\int \|field\|^2 \, dx + \int \|\nabla field\|^2 \, dx}$

These are pure wrappers around DOLFINx scalar assembly, ensuring accurate global norms across distributed meshes.

### High-Level Error Functions in errors_functions.py

For practical verification, [`src/phasefieldx/errors_functions.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/errors_functions.py) provides convenience wrappers that accept two `dolfinx.Function` objects—your numerical solution `uh` and an exact/reference solution `u_ex`—and delegate to the norm utilities:

- `compute_error_L2`, `compute_error_H1`, `compute_error_semiH1` call the corresponding norm on the difference `field_a - field_b`.
- `eval_error_LP` supports arbitrary $L^p$ norms.
- `error_L2_higher_order_space` builds a higher-order function space, interpolates both fields there, and integrates the squared difference—essential when the exact solution lives in a richer space than the FE approximation.
- `eval_error_L2_normalized` returns the relative L² error $\|u_h - u_{ex}\|_{L^2} / \|u_{ex}\|_{L^2}$ with safeguards against division-by-zero.

## Monitoring Solver Convergence

### Newton Solver Configuration in newton.py

Nonlinear solver convergence is managed by the `NewtonSolver` class in [`src/phasefieldx/solvers/newton.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/solvers/newton.py), which wraps PETSc's SNES (Scalable Nonlinear Equations Solvers). The wrapper exposes three critical attributes for controlling convergence criteria:

- `solver.rtol` – relative tolerance for the residual norm.
- `solver.atol` – absolute tolerance for the residual norm.
- `solver.convergence_criterion` – defaults to `"residual"`, meaning SNES reports convergence when the residual norm falls below the specified tolerances.

Setting `solver.report = True` enables verbose logging of the residual norm at each iteration via `logger.info` calls.

### Residual Norm Logging and Inspection

During each SNES iteration, the Newton solver retrieves the residual norm via `snes.getFunctionNorm()`. This value is printed to the log when reporting is enabled, allowing you to track how quickly the nonlinear residual decreases. For phase-field fracture simulations, the residual history is often written to `.conv` files, which can be post-processed using scripts like [`examples/PhaseFieldFracture/plot_1718.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/examples/PhaseFieldFracture/plot_1718.py) to generate convergence plots.

## Complete Workflow Example

The following example demonstrates how to calculate error norms and check convergence in a complete PhaseFieldX simulation workflow:

```python
import dolfinx.mesh, dolfinx.fem
import ufl, numpy as np
from mpi4py import MPI
from phasefieldx.errors_functions import compute_error_L2, compute_error_H1
from phasefieldx.norms import norm_L2, norm_H1
from phasefieldx.solvers.newton import NewtonSolver

# ----------------------------------------------------------------------

# 1️⃣  Set up mesh, function spaces and exact solution

# ----------------------------------------------------------------------

mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 32, 32)
V = dolfinx.fem.FunctionSpace(mesh, ("Lagrange", 1))

# Numerical solution (already computed by a solver)

uh = dolfinx.fem.Function(V)

# ... fill `uh` via your simulation ...

# Exact solution (example: u(x,y) = sin(pi*x)*sin(pi*y))

exact_expr = ufl.sin(np.pi * ufl.x) * ufl.sin(np.pi * ufl.y)
u_ex = dolfinx.fem.Function(V)
u_ex.interpolate(dolfinx.fem.Expression(exact_expr, V.element.interpolation_points))

# ----------------------------------------------------------------------

# 2️⃣  Compute error norms

# ----------------------------------------------------------------------

err_L2 = compute_error_L2(uh, u_ex, mesh)          # → float

err_H1 = compute_error_H1(uh, u_ex, mesh)          # → float

rel_err = err_L2 / norm_L2(u_ex, mesh)             # manual relative error

print(f"L2 error = {err_L2:.3e}, H1 error = {err_H1:.3e}, Relative L2 = {rel_err:.3e}")

# ----------------------------------------------------------------------

# 3️⃣  Run Newton solver and monitor convergence

# ----------------------------------------------------------------------

# Assuming `problem` supplies residual and Jacobian forms

solver = NewtonSolver(problem)
solver.solver.rtol = 1e-8                           # tighter relative tolerance

solver.solver.atol = 1e-9                           # absolute tolerance

solver.solver.report = True                         # prints residual each iteration

solver.solve()                                      # logs "Residual norm u: …"

```

## Summary

- **Error norms** in PhaseFieldX are computed via [`src/phasefieldx/norms.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/norms.py), which provides MPI-aware implementations of L², Lᵖ, H¹-semi, and full H¹ norms.
- **Error calculation** between numerical and exact solutions is streamlined through [`src/phasefieldx/errors_functions.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/errors_functions.py), supporting absolute, relative, and higher-order space projections.
- **Convergence monitoring** is built into the `NewtonSolver` class in [`src/phasefieldx/solvers/newton.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/solvers/newton.py), leveraging PETSc SNES residual norms with configurable relative and absolute tolerances.

## Frequently Asked Questions

### What norm types does PhaseFieldX support?

PhaseFieldX supports L², Lᵖ (for arbitrary p), H¹-semi, and full H¹ norms through the functions in [`src/phasefieldx/norms.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/norms.py). These norms are computed using DOLFINx forms with MPI-aware reduction, ensuring accurate global values across distributed meshes.

### How do I compute relative error norms?

Use the `eval_error_L2_normalized` function from [`src/phasefieldx/errors_functions.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/src/phasefieldx/errors_functions.py), which returns $\|u_h - u_{ex}\|_{L^2} / \|u_{ex}\|_{L^2}$ with built-in safeguards against division-by-zero. Alternatively, manually divide the absolute error computed by `compute_error_L2` by the norm of the exact solution using `norm_L2(u_ex, mesh)`.

### Can I monitor convergence without using the Newton solver?

While the built-in convergence monitoring is integrated into `NewtonSolver` via PETSc SNES, you can implement custom convergence checks by evaluating the residual norm directly within your own solver loop. The `snes.getFunctionNorm()` method (or equivalent DOLFINx/PETSc calls) provides the residual magnitude needed to assess convergence criteria manually.

### Where are convergence histories stored in PhaseFieldX?

For phase-field fracture simulations, convergence histories are typically written to `.conv` files during the simulation. These files can be post-processed using scripts such as [`examples/PhaseFieldFracture/plot_1718.py`](https://github.com/castillonmiguel/phasefieldx/blob/main/examples/PhaseFieldFracture/plot_1718.py) to generate visualizations of residual norms versus iteration count or time steps.