Source-linked AI summary
A penalty method for PDE-constrained optimization in inverse problems
Tristan van Leeuwen, Felix J. Herrmann
TL;DR
PDE-constrained inverse problems must infer parameters from partial measurements while balancing the storage demands of all-at-once optimization against the nonlinear, PDE-solving burden of reduced methods. This paper develops a quadratic-penalty method that eliminates states through an augmented system, achieving comparable computational complexity while retaining a larger search space. Numerical results indicate reduced nonlinearity and lower sensitivity to noise and the initial model.
Problem
PDE-constrained inverse problems require scalable optimization from partial measurements, but all-at-once methods require storing all experiment states and reduced methods enforce constraints at every iteration.
Method
The method uses a quadratic penalty formulation and eliminates the state by solving an overdetermined system that combines discretized PDE equations with measurements.
Results
The penalty method may be less nonlinear, less sensitive to noise and the initial model, and sometimes produce better parameter reconstructions than the conventional reduced approach.
Takeaways & Limitations
The penalty formulation is a viable alternative to the conventional reduced approach for PDE-constrained inverse problems.
Takeaways & Limitations
Augmented Lagrangian alternatives require storing multipliers, which is infeasible for the large-scale problems considered.
Abstract
from arXiv · showhide
Many inverse and parameter estimation problems can be written as PDE-constrained optimization problems. The goal, then, is to infer the parameters, typically coefficients of the PDE, from partial measurements of the solutions of the PDE for several right-hand-sides. Such PDE-constrained problems can be solved by finding a stationary point of the Lagrangian, which entails simultaneously updating the paramaters and the (adjoint) state variables. For large-scale problems, such an all-at-once approach is not feasible as it requires storing all the state variables. In this case one usually resorts to a reduced approach where the constraints are explicitly eliminated (at each iteration) by solving the PDEs. These two approaches, and variations thereof, are the main workhorses for solving PDE-constrained optimization problems arising from inverse problems. In this paper, we present an alternative method that aims to combine the advantages of both approaches. Our method is based on a quadratic penalty formulation of the constrained optimization problem. By eliminating the state variable, we develop an efficient algorithm that has roughly the same computational complexity as the conventional reduced approach while exploiting a larger search space. Numerical results show that this method indeed reduces some of the non-linearity of the problem and is less sensitive the initial iterate.
1. Introduction
PDE-constrained inverse problems infer physical parameters from partial measurements, but standard all-at-once and reduced formulations face scalability or nonlinearity challenges. The paper develops a quadratic-penalty alternative that eliminates states through an augmented physics-and-data system while retaining a larger search space.
- Motivation: Inverse problems infer PDE parameters from partial measurements of solutions across multiple right-hand-sides.The parameters typically appear as PDE coefficients, with applications including geophysics, medical imaging, and non-destructive testing.
- Problem setting: Multiple experiments produce block vectors of states, sources, and observations, with block-diagonal PDE and sampling matrices.Typical seismic problems involve large dimensions for parameters, states, experiments, and measurements.
- Existing approaches: All-at-once methods jointly optimize parameters, states, and adjoint variables but may be infeasible because storing all experiment states is too costly.Reduced methods avoid this storage but require PDE solves and restrict the search space.
- Proposed approach: The proposed method derives an efficient algorithm from a quadratic penalty formulation and eliminates states by solving an overdetermined system combining PDE equations with measurements.This approach is motivated by recovering states from incomplete measurements and extends the equation-error idea beyond invertible measurement operators.
- Contribution: A finite penalty parameter can approximate the constrained solution to arbitrary finite precision, with experiments suggesting that one fixed value is often sufficient.The method need only enforce dominant constraints when unresolved null-space components do not require high precision.
2. All-at-once and reduced methods
The paper contrasts all-at-once optimization, which delays constraint satisfaction but requires storing all variables, with reduced optimization, which satisfies PDE constraints at every iteration but can be computationally expensive and highly nonlinear.
- All-at-once approach: All-at-once methods seek a stationary point of the Lagrangian by jointly updating parameters, states, and adjoint variables.Newton-like methods can avoid explicit PDE solves and satisfy constraints only approximately until convergence.
- All-at-once approach: All-at-once optimization is often infeasible for large-scale applications because it requires simultaneously storing all variables.Its advantage is avoiding explicit PDE solves, while constraints are only approximately satisfied during optimization.
- Reduced approach: Reduced formulations eliminate the state with u_red(m) = A(m)^−1q, yielding a smaller unconstrained problem whose iterates satisfy the PDE constraints.They can be solved with black-box nonlinear optimization methods.
- Reduced approach: Reduced optimization requires PDE solves at each update and strictly enforces the constraint, potentially making the parameter dependence more nonlinear.These properties make the approach computationally expensive and may hinder optimization.
3. Penalty and augmented Lagrangian methods
Penalty methods replace constrained optimization with unconstrained objectives that penalize constraint violations, while augmented Lagrangian methods additionally update multipliers. The paper motivates a differentiable quadratic penalty formulation whose state subproblem has an efficiently computable closed-form solution.
- Penalty methods: Penalty formulations add a nonnegative function of constraint violation, with sufficiently large λ making the unconstrained and constrained solutions coincide.Quadratic penalties yield differentiable objectives, but their minimizers generally approach the constrained solution as λ increases.
- Penalty methods: Quadratic penalty algorithms commonly solve a sequence of unconstrained problems with increasing λ, although large λ can cause ill-conditioning.For PDE-constrained inverse problems, invertibility of A(m) supplies enough constraints to avoid the cited limiting issue.
- Penalty methods: Nonsmooth penalties can reach a constrained solution at a finite λ, but they sacrifice differentiability and still typically use continuation.This trades arbitrarily large penalty parameters for a nonsmooth optimization problem.
- Augmented Lagrangian methods: Augmented Lagrangian methods add a quadratic penalty to the Lagrangian and update multipliers and λ, but storing multipliers is infeasible for the target problems.ADMM is a basic example of this approach.
- Paper’s formulation: The paper uses a quadratic penalty because it produces a differentiable unconstrained problem whose optimization over u has a closed-form solution that can be computed efficiently.This structure makes the formulation suitable for large-scale PDE-constrained inverse problems.
4. A reduced penalty method
The reduced penalty method eliminates the state by solving a data- and physics-informed overdetermined system, yielding a reduced objective whose derivatives are computable. Its finite-penalty iterates approximate stationary points of the original constrained problem while retaining a broader, less nonlinear search space.
- Formulation: The quadratic penalty formulation eliminates u and defines a reduced objective φλ(m) by solving an overdetermined system that combines the PDE and measurements.Unlike the conventional reduced approach, the state estimate is not obtained from the PDE alone.
- Formulation: The modified system A^T A + λ^-1P^T P is a low-rank modification of the original PDE and incorporates data into the state solve.This makes the state estimate depend on both the physics, current model, and observations.
- Algorithm: The reduced penalty algorithm computes the gradient and Hessian of φλ and obtains the adjoint state without an additional PDE solve.Forward and adjoint solves are performed simultaneously through the normal equations.
- Analysis: For fixed m, the penalty states uλ and vλ are related to the conventional reduced states ured and vred, with uλ = ured + O(λ^-1).The analysis also gives a scaling interpretation for large and small penalty parameters.
- Analysis: At each iteration, the constraint residual is O(λ^-1), and successful termination yields an O(λ^-1) Lagrangian-gradient discrepancy from a stationary point.A finite λ can therefore provide finite accuracy when the penalty contribution is sufficiently small relative to the model gradient.
- Scope and accuracy: In ill-posed problems, acceptable reconstructions may have relatively large parameter or gradient errors because only components associated with dominant Hessian eigenvectors are recoverable.Further reducing the gradient mainly refines increasingly weak eigenmodes and may barely change reconstruction or data fit.
5. Algorithm
The algorithm solves the penalty method’s augmented systems using direct or iterative linear algebra, with spectral and numerical evidence suggesting that these systems can remain tractable. Its storage and leading computational costs match the reduced approach, although continuation in λ may add iterations.
- Implementation: Independent experiments produce a block-diagonal system, so the corresponding linear systems can be solved independently and in parallel.The blocks correspond to the K experiments and each has size N × N.
- Linear solvers: Sparse augmented systems can be solved by QR or Cholesky factorization, while industry-scale applications typically require preconditioned iterative methods such as LSQR, LSMR, or CGLS.When P^T P is nearly diagonal, the augmented system has a similar sparsity pattern to the original.
- Iterative methods: The augmented system is a rank-L modification of A^T A, with eigenvalue shifts whose coefficients sum to L.At worst, one eigenvalue is shifted by Lλ^-1.
- Iterative methods: For the tested parabolic and Helmholtz examples, the augmented system has a lower condition number than the original system, while eigenvalue distributions change little.The authors therefore expect iterative methods to perform similarly, although effective preconditioning remains outside the paper’s scope.
- Complexity: The penalty, reduced, and all-at-once approaches have summarized per-iteration costs, with penalty and reduced methods requiring storage and computation of the same order in K, N, and M.Their PDE solves can also be performed independently and in parallel.
- Complexity: A continuation strategy in λ adds cost because the penalty and reduced methods may require different numbers of iterations, even if augmented systems are solved as efficiently as the original PDE.The paper explicitly leaves their iteration-count comparison unresolved before numerical testing.
6. Case studies
The case studies evaluate the penalty method on 1D and 2D inverse problems, including tests of nonlinearity, computational behavior, noise sensitivity, and sensitivity to initialization. Across these experiments, the method follows the predicted penalty behavior and often matches or improves on the reduced approach.
- Computational comparison: The penalty method can achieve practically the same or better results than the reduced method at lower computational cost.In the reported QN comparison, it converges in fewer iterations and uses fewer PDE solves per iteration.
- Convergence: The penalty method converges to the predicted O(λ−1) stationarity behavior, and continuation in λ reaches the desired tolerance.This behavior is reported for both the 1D experiment and the 2D Gauss–Newton experiment.
- Non-linearity: The reduced misfit has strong nonlinearity and spurious stationary points, whereas the penalty misfit is better behaved for small λ.For larger λ, the penalty misfit increasingly resembles the reduced misfit.
- Noise sensitivity: With 10% and 20% Gaussian noise, the penalty approach is not overly sensitive and produces very similar or slightly better results than the reduced approach.The noise experiments use the QN method.
- Sensitivity to the initial guess: For a difficult initial guess, small λ lets the penalty method reach a data-fitting solution similar to that obtained from a better initial guess, unlike the reduced method.The reduced method stalls at relatively high data misfit, while the penalty method avoids the reported cycle-skipping behavior in this experiment.
7. Discussion
The discussion identifies practical uncertainty in selecting and continuing the penalty parameter, while also describing reduced initialization sensitivity and exploitable Hessian structure. Further work is needed to characterize robustness and develop a reliable continuation strategy.
- Penalty-parameter selection: Choosing a suitable penalty parameter λ a priori remains unresolved, although a small fraction of ∥PA−1∥2 at the initial iterate gives good results.The authors propose continuation when the method stalls, but do not establish a robust strategy.
- Continuation: A continuation strategy is needed to guarantee finding a Lagrangian stationary point at a preset tolerance, and its robust implementation requires further study.The suggested strategy detects stalling and then reduces λ.
- Robustness: Relaxing the PDE and data constraints permits errors in physics and data while reducing sensitivity of the final reconstruction to the initial guess.The authors state that the robustness of this effect still requires characterization.
- Hessian structure: The penalty Gauss–Newton Hessian is full rank and admits a natural sparse approximation Hλ ≈ λG^TG, unlike the reduced Gauss–Newton Hessian.The reduced Hessian has rank at most ML and lacks the same natural sparse approximation.
8. Conclusions
The paper presents a quadratic-penalty method that eliminates states through an augmented PDE-and-measurement system, retaining reduced-like computational complexity while exploring a larger search space. Theory and numerical examples show finite-λ stationarity, reduced nonlinearity, and lower sensitivity to noise and initialization.
- Method: The method uses a quadratic penalty formulation and efficiently eliminates state variables by solving an overdetermined linear system.The formulation is unconstrained in parameters and states before state elimination.
- Method: The penalty approach has similar computational complexity to the reduced method while exploring a larger search space because it does not enforce PDE constraints exactly.This is the principal algorithmic trade-off identified in the conclusions.
- Theory: For sufficiently large λ, the method can theoretically find a constrained-problem Lagrangian stationary point within a prescribed tolerance at finite precision.Although the theory considers λ increasing, finite λ suffices for finite-precision accuracy.
- Algorithmic properties: The augmented system combines the discretized PDE and measurements, and its structure and eigenvalues are not dramatically altered relative to the original PDE.This supports efficient elimination of the state variables.
- Numerical behavior: Numerical examples obtain very good results with a single relatively small λ, while continuation can gradually increase λ to reach the desired tolerance.The continuation result is described as viable rather than fully established.
- Numerical conclusions: The penalty method can reduce nonlinearity and is less sensitive to noise and the initial model than the conventional reduced approach.Some experiments also produce better parameter reconstructions than the reduced method.