|
| 1 | +# The Adjoint Method |
| 2 | + |
| 3 | +Many simulation tasks ask a question that runs in the opposite direction to a forward solve: which material value, load, actuator, or design parameter makes |
| 4 | +the simulated result agree with an observation? |
| 5 | +The quantity being optimized is a loss on a state obtained only after solving an equilibrium problem. |
| 6 | +This is sensitivity analysis through a simulation, and it is what the adjoint layer in SymX is for. |
| 7 | + |
| 8 | +The layer deliberately keeps the forward model unchanged. |
| 9 | +You declare and solve a `NewtonsMethod` exactly as in the [global-optimization layer](newtons_method.md). |
| 10 | +Then you declare a separate loss and give both objects to `AdjointNewton`. |
| 11 | +The result is an efficient gradient with respect to all loss parameters, independent of how many of them there are. |
| 12 | + |
| 13 | +SymX is perfect for this type of problems as it has access to the symbolic representation of the entire problem and can handle with the annoying task of calculating cumbersome cross derivatives. |
| 14 | + |
| 15 | +## Equilibrium sensitivities |
| 16 | + |
| 17 | +Let $u$ be the simulation DoFs, $p$ a set of parameters, and $P(u,p)$ the physical potential. |
| 18 | +The converged forward state $u^*(p)$ satisfies |
| 19 | + |
| 20 | +$$ |
| 21 | +P_u(u^*(p), p) = 0. |
| 22 | +$$ |
| 23 | + |
| 24 | +We want the total derivative of a loss $L(u^*(p),p)$. |
| 25 | +Differentiating the equilibrium condition would require a linear solve for every parameter: |
| 26 | + |
| 27 | +$$ |
| 28 | +P_{uu}\frac{d u^*}{d p} + P_{up} = 0. |
| 29 | +$$ |
| 30 | + |
| 31 | +Instead, SymX solves one adjoint system, |
| 32 | + |
| 33 | +$$ |
| 34 | +P_{uu}^{\mathsf T}\lambda = -L_u^{\mathsf T}, |
| 35 | +$$ |
| 36 | + |
| 37 | +and evaluates the gradient as |
| 38 | + |
| 39 | +$$ |
| 40 | +\frac{dL}{dp} = L_p + \lambda^{\mathsf T} P_{up}. |
| 41 | +$$ |
| 42 | + |
| 43 | +That one additional solve is why the method is especially useful when the number of parameters is large. |
| 44 | +The forward Hessian is reused as the adjoint system matrix; choose the adjoint linear solver just as deliberately as the one used by Newton. |
| 45 | + |
| 46 | +## Declare the forward problem and loss |
| 47 | + |
| 48 | +The forward `GlobalPotential` owns the state DoFs. |
| 49 | +The loss is another `GlobalPotential`, whose DoFs are the quantities to optimize. |
| 50 | +Both potentials use the same mapped-stencil declarations, so a loss can naturally read the converged state and any parameter arrays. |
| 51 | + |
| 52 | +```cpp |
| 53 | +spGlobalPotential physical = GlobalPotential::create(); |
| 54 | +// physical->add_potential(...); // defines P(u, p) |
| 55 | +physical->add_dof(state); // u |
| 56 | + |
| 57 | +auto newton = NewtonsMethod::create(physical, context); |
| 58 | + |
| 59 | +spGlobalPotential loss = GlobalPotential::create(); |
| 60 | +// loss->add_potential(...); // defines L(u, p) |
| 61 | +loss->add_dof(parameters, "parameters"); // p |
| 62 | + |
| 63 | +LinearSolveSettings adjoint_linear; |
| 64 | +adjoint_linear.linear_solver = LinearSolver::DirectLLT; |
| 65 | +auto adjoint = AdjointNewton::create(newton, loss, adjoint_linear); |
| 66 | +``` |
| 67 | +
|
| 68 | +The physical potential may depend on `parameters` without registering them as simulation DoFs. |
| 69 | +Conversely, the loss does not register the state as a loss parameter: it reads it from the forward model. |
| 70 | +Register every independent parameter set that should receive a gradient with `loss->add_dof`. |
| 71 | +
|
| 72 | +## Evaluate a value and gradient |
| 73 | +
|
| 74 | +`run_forward_and_evaluate_L()` solves the forward equilibrium and evaluates the loss. |
| 75 | +Once that state is current, `evaluate_dL_dp_no_solve()` computes its adjoint gradient without solving the forward problem again. |
| 76 | +
|
| 77 | +```cpp |
| 78 | +AdjointValueResult value = adjoint->run_forward_and_evaluate_L(); |
| 79 | +if (!value.success) { |
| 80 | + // value.status explains whether the forward or loss evaluation failed. |
| 81 | +} |
| 82 | +
|
| 83 | +AdjointGradientResult gradient = adjoint->evaluate_dL_dp_no_solve(); |
| 84 | +if (!gradient.success) { |
| 85 | + // gradient.status explains the failure. |
| 86 | +} |
| 87 | +
|
| 88 | +double L = value.value; |
| 89 | +const Eigen::VectorXd& dL_dp = gradient.gradient; |
| 90 | +``` |
| 91 | + |
| 92 | +The `run_forward_and_evaluate_dL_dp()` convenience method performs both operations. |
| 93 | +The `*_no_solve` methods are intentional: they let an outer optimizer reuse an equilibrium state, but require the caller to keep that state valid. |
| 94 | +A failed value is infinity and a failed gradient is empty, so always inspect `success` before consuming either result. |
| 95 | + |
| 96 | +## Verify an implementation |
| 97 | + |
| 98 | +Finite differences remain the best compact check of a new inverse model. |
| 99 | +The adjoint helper evaluates centered finite differences, restoring both the parameters and simulation state even if an evaluation fails. |
| 100 | + |
| 101 | +```cpp |
| 102 | +auto analytic = adjoint->run_forward_and_evaluate_dL_dp(); |
| 103 | +auto finite_difference = adjoint->evaluate_dL_dp_finite_differences(1e-3); |
| 104 | + |
| 105 | +if (analytic.success && finite_difference.success) { |
| 106 | + double relative_error = |
| 107 | + (analytic.gradient - finite_difference.gradient).norm() / |
| 108 | + finite_difference.gradient.norm(); |
| 109 | +} |
| 110 | +``` |
| 111 | +
|
| 112 | +For temporarily fixed scalar components, use `force_gradient_to_zero(indices)`. |
| 113 | +Indices address the flattened vector of all loss DoF sets and are validated by the API. |
| 114 | +
|
| 115 | +The complete workflow, including a target deformation, finite-difference check, and L-BFGS solve, is developed in [Inverse Simulation](inverse_simulation.md). |
0 commit comments