Calculus of variations, Euler–Lagrange and physics-informed losses

Before this lesson: lesson 02 (the directional derivative), lesson 03 (critical points and second-order conditions), lesson 10 (integration by parts). For section 3 it helps to know positive definite matrices (linear algebra lesson 07) and condition numbers (linear algebra lesson 13). Time: about 80 minutes, plus the notebook.

---

Why this matters #

Many laws of physics say that some integral is as small as possible, or at least stationary. An elastic body at rest minimizes its potential energy. The Poisson equation is the condition for the minimum of the Dirichlet energy. Newton's law and the wave equation make an "action" stationary (Hamilton's principle). The calculus of variations gives the derivative of a number $J[u]$ with respect to a whole function $u$, and the rule that turns "$J$ is stationary" into a differential equation.

For ML for science this gives two kinds of physics losses for a network $u_\theta$:

• an energy (Ritz) loss minimizes the integral itself (the Deep Ritz method);
• a residual (PINN) loss penalizes the squared residual of the differential equation (physics-informed neural networks).

They need different derivatives, treat boundary conditions differently and fail in different ways. The same "gradient of a functional" also gives maximum-entropy distributions, the mean-field update of variational inference, gradient flows, and the adjoint equations of weather data assimilation (4D-Var) and seismic full-waveform inversion.

1. Predict first #

2. The picture #

A functional maps a function to a number. Two examples on $[a, b]$:

• the length of the graph of $u$: $\;J[u] = \int_a^b \sqrt{1 + u'(x)^2}\,dx$;
• the Dirichlet energy with a load $f$: $\;J[u] = \int_a^b \big(\tfrac12 u'(x)^2 - f(x)\,u(x)\big)\,dx$. For an elastic string, $\tfrac12 u'^2$ is the stored energy of the stretch and $-fu$ is the work of the load.

The general form is $J[u] = \int_a^b L(x, u, u')\,dx$. The function $L(x, u, p)$ is the Lagrangian; $p$ stands for the slope $u'$.

A functional on a grid is a function of many variables #

Put $u$ on a grid to see $J$ as an ordinary function. Take $x_i = a + ih$ for $i = 0, \dots, N$, with $h = (b - a)/N$, and values $U_i \approx u(x_i)$. For the Dirichlet energy a simple discretization is

with the boundary values $U_0$ and $U_N$ fixed. A minimizer has $\partial J_h/\partial U_i = 0$ at every interior point. Let us compute this partial derivative step by step.

1. $U_i$ appears in two slopes: $p_{i-1/2} = (U_i - U_{i-1})/h$ (the cell on its left) and $p_{i+1/2} = (U_{i+1} - U_i)/h$ (the cell on its right).
2. The left term $h \cdot \tfrac12 p_{i-1/2}^2$ has derivative $h \cdot p_{i-1/2} \cdot \frac1h = p_{i-1/2}$.
3. The right term $h \cdot \tfrac12 p_{i+1/2}^2$ has derivative $h \cdot p_{i+1/2} \cdot \big(-\frac1h\big) = -p_{i+1/2}$.
4. The load term gives $-h f(x_i)$. Together:

The last step uses the second difference (lesson 07): $(U_{i+1} - 2U_i + U_{i-1})/h^2 \approx u''(x_i)$. Three things to see here.

• The factor $h$. $U_i$ controls a piece of length $h$ of the integral, so its partial derivative is $h$ times a density. Divided by $h$, the gradient becomes a function of $x$: the functional derivative $\delta J/\delta u = -u'' - f$.
• Where $d/dx$ comes from. $U_i$ sits in two neighbouring slopes, so the gradient contains the difference of $\partial L/\partial p = p$ across the cell. That difference, divided by $h$, is $-\frac{d}{dx}\partial L/\partial p$. This is integration by parts on a grid.
• The equation. Setting the gradient to zero gives $-u'' = f$, the Poisson equation, with the fixed boundary values. The minimum of a functional is the solution of a differential equation.

Interactive figure: loading…

Variations are directions in function space #

In the continuous version, we perturb $u$ in a direction $\eta$. The direction is itself a function, and it must be zero where $u$ is fixed (here at both ends), so that $u + \epsilon\eta$ still has the right boundary values. Then $\varphi(\epsilon) = J[u + \epsilon\eta]$ is an ordinary function of one number $\epsilon$. If $u$ is a minimizer, $\varphi$ has its minimum at $\epsilon = 0$, so $\varphi'(0) = 0$. This is the directional derivative $\nabla f(x)^\top v$ of lesson 02, with a function $\eta$ in place of the vector $v$.

Interactive figure: loading…

3. Definitions and main results #

Setup. $L(x, u, p)$ is a twice continuously differentiable ($C^2$) function of three real variables. The functional is $J[u] = \int_a^b L(x, u(x), u'(x))\,dx$. The admissible functions are $u \in C^2[a, b]$ with $u(a) = \alpha$ and $u(b) = \beta$. We write $L_u = \partial L/\partial u$, $L_p = \partial L/\partial p$, $L_{pp} = \partial^2 L/\partial p^2$, and so on, always evaluated along $(x, u(x), u'(x))$.

First variation. For $\eta \in C^2[a, b]$ with $\eta(a) = \eta(b) = 0$, the function $u + \epsilon\eta$ is admissible for every $\epsilon$. The first variation is

To get the integral, differentiate under the integral sign and use the chain rule: $\frac{d}{d\epsilon} L(x, u + \epsilon\eta, u' + \epsilon\eta') = L_u\,\eta + L_p\,\eta'$.

Functional derivative. If a function $\delta J/\delta u$ satisfies $\delta J[u; \eta] = \int_a^b \frac{\delta J}{\delta u}(x)\, \eta(x)\,dx$ for every such $\eta$, it is the functional derivative of $J$ at $u$. It is the gradient of $J$ for the $L^2$ inner product $\langle g, \eta \rangle = \int g\,\eta\,dx$, in the same way as $\nabla f$ satisfies $Df(x)\,v = \nabla f(x)^\top v$. On a grid of spacing $h$: $\partial J_h/\partial U_i \approx h\, \frac{\delta J}{\delta u}(x_i)$.

The fundamental lemma #

Lemma (fundamental lemma of the calculus of variations). Let $g$ be continuous on $[a, b]$, and let $\int_a^b g\,\eta\,dx = 0$ for every smooth $\eta$ that is zero outside some closed interval inside $(a, b)$. Then $g(x) = 0$ for every $x \in [a, b]$.

The idea: if $g$ were positive somewhere, a bump $\eta$ placed exactly there would make $\int g\eta > 0$.

Proof of the fundamental lemma

1. Suppose $g(x_0) \ne 0$ for some $x_0 \in (a, b)$. Replacing $g$ by $-g$ if needed, $g(x_0) > 0$.
2. By continuity, $g > 0$ on a closed interval $[x_0 - \delta, x_0 + \delta] \subset (a, b)$ for some $\delta > 0$.
3. Take the smooth bump $\eta(x) = \exp\big(-1/(\delta^2 - (x - x_0)^2)\big)$ inside this interval and $\eta = 0$ outside. It is smooth, positive inside and zero outside.
4. Then $g\eta > 0$ inside the interval and $g\eta = 0$ outside, so $\int g\eta\,dx > 0$. This contradicts the assumption.
5. So $g = 0$ on $(a, b)$, and by continuity also at $a$ and $b$. $\blacksquare$

The Euler–Lagrange equation #

Theorem (Euler–Lagrange equation). Let $L \in C^2$, and let the admissible $u \in C^2$ be a local minimizer of $J$: $J[u] \le J[v]$ for every admissible $v$ with $\max|v - u| + \max|v' - u'|$ small enough. Then

The proof uses integration by parts (lesson 10) in this form:

Proof.

1. Fix $\eta \in C^2$ with $\eta(a) = \eta(b) = 0$. For small $|\epsilon|$, $u + \epsilon\eta$ is admissible and close to $u$, so $\varphi(\epsilon) = J[u + \epsilon\eta]$ has a local minimum at $\epsilon = 0$. Hence $\varphi'(0) = 0$.
2. $\varphi'(0) = \int_a^b (L_u\,\eta + L_p\,\eta')\,dx$. (Differentiation under the integral is allowed: $L$ is $C^2$ and $[a, b]$ is a closed bounded interval.)
3. Integrate the second term by parts (the formula above). The boundary term $\big[L_p\,\eta\big]_a^b$ is $0$ because $\eta(a) = \eta(b) = 0$. ($\frac{d}{dx}L_p$ is continuous because $L \in C^2$ and $u \in C^2$.)
4. So $\int_a^b \big(L_u - \tfrac{d}{dx}L_p\big)\,\eta\,dx = 0$ for every such $\eta$.
5. $g = L_u - \frac{d}{dx}L_p$ is continuous, so the fundamental lemma gives $g = 0$. $\blacksquare$

Step 4 also shows that

The derivative $\frac{d}{dx}$ is a total derivative along the curve: $L_p$ depends on $x$ directly and through $u(x)$ and $u'(x)$. By the chain rule,

So the Euler–Lagrange equation $L_u - L_{px} - L_{pu}u' - L_{pp}u'' = 0$ is a second-order ODE when $L_{pp} \ne 0$. It comes with two boundary conditions. Example: for $L = \tfrac12 x\,u'^2$ we have $L_p = xu'$, so $\frac{d}{dx}L_p = u' + xu''$, not only $xu''$.

Natural boundary conditions #

Now fix only $u(a) = \alpha$ and leave $u(b)$ free. Then the directions $\eta$ must have $\eta(a) = 0$, but $\eta(b)$ can be anything. Step 3 of the proof now keeps one boundary term, at $x = b$:

1. First take all $\eta$ with $\eta(b) = 0$. As before, the Euler–Lagrange equation holds, so the integral is $0$ for every $\eta$.
2. Then take an $\eta$ with $\eta(b) = 1$. The equation above becomes $L_p\big(b, u(b), u'(b)\big) = 0$.

This is the natural boundary condition. You do not impose it; the minimization produces it. For $L = \tfrac12 u'^2 + \dots$ it is $u'(b) = 0$, a Neumann condition. Prescribed values such as $u(a) = \alpha$ are called essential (Dirichlet) conditions.

When is a solution a minimum? #

The Euler–Lagrange equation is a necessary condition: maxima and saddle points of $J$ satisfy it too. A sufficient condition:

Theorem. If $L(x, u, p)$ is convex in $(u, p)$ for every $x$, then every admissible solution of the Euler–Lagrange equation is a global minimizer.

Proof. Let $v$ be admissible and $\eta = v - u$, so $\eta(a) = \eta(b) = 0$. Convexity says that $L$ lies above its tangent plane: $L(x, v, v') \ge L(x, u, u') + L_u\,\eta + L_p\,\eta'$ at every $x$. Integrate:

Without convexity, the second variation $\varphi''(0)$ decides, like the Hessian in lesson 03. A necessary condition for a minimum is $L_{pp} \ge 0$ along $u$ (the Legendre condition).

Several variables #

For $u: \Omega \subset \R^d \to \R$ and $J[u] = \int_\Omega L(x, u, \nabla u)\,dx$, write $p = \nabla u \in \R^d$ and let $\nabla_p L$ be the vector of partial derivatives of $L$ in $p$. The same steps work with the divergence theorem in place of integration by parts:

with $n$ the outward unit normal. The Euler–Lagrange equation is

and on a part of the boundary $\partial\Omega$ where $u$ is free, the natural condition is $\nabla_p L \cdot n = 0$.

Classic cases #

1. Shortest path. Arc length: $L = \sqrt{1 + p^2}$. $L_u = 0$, so the equation says $\frac{d}{dx}L_p = 0$, that is, $L_p = u'/\sqrt{1 + u'^2} = C$, a constant. The map $p \mapsto p/\sqrt{1 + p^2}$ is strictly increasing, so $u'$ is constant: a straight line. $L$ is convex in $p$ ($L_{pp} = (1 + p^2)^{-3/2} > 0$) and does not depend on $u$, so the line is the global minimizer.
2. Dirichlet energy and the Poisson equation. $E[u] = \int_\Omega \big(\tfrac12|\nabla u|^2 - f u\big)\,dx$ with $u = 0$ on $\partial\Omega$. Here $L_u = -f$ and $\nabla_p L = p = \nabla u$, so the equation is $-f - \nabla \cdot \nabla u = 0$:

$L$ is convex in $(u, p)$, so the solution is the minimizer (Dirichlet's principle). Where no boundary value is imposed, the natural condition is $\partial u/\partial n = 0$. In 1-D: $-u'' = f$.

Hamilton's principle: stationary, not minimal

In mechanics the integral runs over time. The action $S[q] = \int \big(\tfrac12 m\dot q^2 - V(q)\big)\,dt$ has $L_q = -V'(q)$ and $L_{\dot q} = m\dot q$, so the Euler–Lagrange equation is Newton's law $m\ddot q = -V'(q)$. For a 1-D elastic medium with density $\rho$ and stiffness $\mu$, $S[u] = \iint \big(\tfrac12\rho u_t^2 - \tfrac12\mu u_x^2\big)\,dx\,dt$ gives the wave equation $\rho u_{tt} = \mu u_{xx}$ (the elastic wave equation of seismology is the 3-D version). This $L$ is not convex in $(u_t, u_x)$: the motion makes the action stationary, not minimal.

Energy (Ritz) losses and residual (PINN) losses #

Take $-\Delta u = f$ in $\Omega$, $u = 0$ on $\partial\Omega$, and a neural network $u_\theta$.

Ritz (energy) loss. Minimize $E[u_\theta]$ over $\theta$. The Deep Ritz method estimates the integral by Monte Carlo, with new uniform points $x_j \in \Omega$ at every step:

The Dirichlet condition is either a penalty $\beta\int_{\partial\Omega} u_\theta^2\,dS$, or it is built into the network, $u_\theta(x) = d(x)\,N_\theta(x)$ with a function $d$ that is $0$ on $\partial\Omega$. The classical Ritz method uses $u_c = \sum_k c_k \varphi_k$ with fixed basis functions $\varphi_k$; the finite element method is this method with piecewise polynomials.

Upper bound. If the trial functions are admissible (essential conditions exact) and the integral is exact, they form a subset of the admissible set. So $\min_\theta E[u_\theta] \ge \min_u E[u]$: the Ritz value approaches the true minimum from above. For the Dirichlet energy the gap is exactly an error: $E[u_\theta] - E[u^*] = \tfrac12\int_\Omega |\nabla(u_\theta - u^*)|^2\,dx$, because $\delta E[u^*; \cdot] = 0$ and $E$ is quadratic.

PINN (strong-form residual) loss. With collocation points $x_j \in \Omega$ and boundary points $y_k \in \partial\Omega$:

	Energy loss (Ritz, Deep Ritz)	Residual loss (PINN)
Applies to	PDEs that are Euler–Lagrange equations of a functional (symmetric problems: Poisson, elasticity)	any PDE, also non-variational ones (advection, Navier–Stokes)
Derivatives of $u_\theta$	half the order: first derivatives for Poisson	full order: second derivatives for Poisson
Loss at the solution	the minimum energy, not $0$	$0$; the residual is an error indicator
Boundary conditions	essential ones must be imposed; natural ones come free	every condition needs a loss term (or must be built in)
Regularity	weak solutions; kinks and corner singularities are fine	needs a smooth $u_\theta$ (strong solution)
Conditioning	like $A$: condition number $\kappa$	like $A^\top A$: condition number $\kappa^2$

Why the residual loss is worse conditioned. Discretize $-u''$ on $N = 100$ interior points to get the matrix $A = \mathrm{tridiag}(-1, 2, -1)/h^2$. It is symmetric positive definite.

1. The energy $\tfrac12 U^\top A U - F^\top U$ has gradient $AU - F$ and Hessian $A$, with $\kappa(A) \approx 4.1 \cdot 10^3$.
2. The residual loss $\tfrac12\|AU - F\|^2$ has gradient $A^\top(AU - F)$ and Hessian $A^\top A$. Its eigenvalues are the squares of those of $A$, so $\kappa(A^\top A) = \kappa(A)^2 \approx 1.7 \cdot 10^7$.
3. Plain gradient descent needs a number of steps proportional to $\kappa$ (linear algebra lesson 14 shows how CG reduces this to $\sqrt\kappa$).

Both losses have the same minimizer $U^* = A^{-1}F$. Only the speed differs. For neural networks the picture is less simple, but Krishnapriyan et al. (see "Watch and read") show that the differential operator in the loss makes PINN training hard.

Interactive figure: loading…

Two more facts. A linear problem $\mathcal A u = f$ is the Euler–Lagrange equation of $\tfrac12\langle u, \mathcal A u\rangle - \langle f, u\rangle$ only if $\mathcal A$ is symmetric (self-adjoint) in $L^2$; an advection term $b \cdot \nabla u$ is not, so such problems have no energy of this form. And a ReLU network is piecewise linear, so $\Delta u_\theta = 0$ almost everywhere: a PINN for a second-order PDE needs smooth activations (tanh, SiLU, sine). Deep Ritz needs only first derivatives, but E and Yu still use a smoother activation, $\max(x^3, 0)$, for accuracy.

The functional derivative and gradient descent #

Since $\partial J_h/\partial U_i \approx h\,\delta J/\delta u(x_i)$, a gradient-descent step $U \leftarrow U - \tau\,\nabla J_h$ is $U_i \leftarrow U_i - \tau h\,\frac{\delta J}{\delta u}(x_i)$. This is an explicit Euler step of length $\tau h$ for the gradient flow

For the Dirichlet energy $\delta J/\delta u = -u'' - f$, so this flow is the heat equation $u_t = u'' + f$. Gradient descent smooths $u$ the way heat spreads. Explicit Euler for the heat equation is stable only for time steps up to $h^2/2$, so $\tau h \le h^2/2$, that is $\tau \le h/2$. The number of steps to reach a fixed time grows like $1/h^2$. This is the conditioning problem again: plain gradient descent on a fine grid is slow.

Interactive figure: loading…

Answers to section 1. (1) $L = p^2$ gives $L_u = 0$ and $L_p = 2u'$, so the equation is $2u'' = 0$. With the boundary values, $u = x$, and the minimum is $\int_0^1 1\,dx = 1$. No other admissible $u$ does better: by the Cauchy–Schwarz inequality, $1 = \big(\int_0^1 u'\,dx\big)^2 \le \int_0^1 u'^2\,dx$. (2) Yes: the Dirichlet energy $\int\big(\tfrac12 u'^2 - f u\big)\,dx$. Its Euler–Lagrange equation is $-u'' = f$, and it needs only $u'$. This is the Ritz loss.

4. Worked example #

Minimize $J[u] = \int_0^1 \big(\tfrac12 u'^2 - u\big)\,dx$ with $u(0) = 0$ and $u(1)$ free.

Step 1: Euler–Lagrange. $L = \tfrac12 p^2 - u$, so $L_u = -1$ and $L_p = p$. The equation $L_u - \frac{d}{dx}L_p = -1 - u'' = 0$ gives $u'' = -1$. This is the Poisson equation $-u'' = f$ with $f = 1$.

Step 2: boundary conditions. Essential: $u(0) = 0$. Natural, at the free end: $L_p = u'(1) = 0$.

Step 3: solve. Integrate twice: $u = -x^2/2 + Ax + B$. From $u(0) = 0$: $B = 0$. From $u'(1) = -1 + A = 0$: $A = 1$. So $u^*(x) = x - x^2/2$.

Step 4: value. $u^{*\prime} = 1 - x$. The two parts are $\int_0^1 \tfrac12(1 - x)^2\,dx = \tfrac16$ and $\int_0^1 (x - x^2/2)\,dx = \tfrac12 - \tfrac16 = \tfrac13$. So $J[u^*] = \tfrac16 - \tfrac13 = -\tfrac16$. $L$ is convex in $(u, p)$, so this is the global minimum.

Step 5: Ritz with one trial function, $u = cx$. It satisfies $u(0) = 0$, but not $u'(1) = 0$. $J(c) = \int_0^1 (\tfrac12 c^2 - cx)\,dx = c^2/2 - c/2$. The minimum is at $c = \tfrac12$, with $J = \tfrac18 - \tfrac14 = -\tfrac18 = -0.125$. This is above $-\tfrac16 \approx -0.1667$, as the upper bound says. The slope at $x = 1$ is $\tfrac12$, not $0$: the natural condition holds only approximately.

Step 6: Ritz with two trial functions, $u = c_1 x + c_2 x^2$. Then $u' = c_1 + 2c_2 x$ and

Setting both partial derivatives to zero gives $c_1 + c_2 = \tfrac12$ and $c_1 + \tfrac43 c_2 = \tfrac13$. Subtract: $\tfrac13 c_2 = -\tfrac16$, so $c_2 = -\tfrac12$ and $c_1 = 1$. This is exactly $u^*$, with $J = -\tfrac16$, and $u'(1) = 0$ holds although nobody imposed it. The trial space contains the solution, so Ritz finds it.

Step 7: the PINN loss for the same problem is

The condition $u'(1) = 0$ needs its own term and its own weight $\lambda_1$. The Ritz loss contains it already; only $u(0) = 0$ must be imposed, for example with $u_\theta(x) = x\,N_\theta(x)$.

5. In code #

NumPy and SciPy check the worked example: an Euler–Lagrange residual by finite differences (with the total derivative), the Ritz coefficients, the natural boundary condition on a grid, and the two condition numbers.

import numpy as np
from scipy import integrate

# (a) Euler-Lagrange residual L_u - d/dx L_p, with the total derivative d/dx by central differences
def el_residual(Lu, Lp, u, du, x, h=1e-4):
    dLp = (Lp(x + h, u(x + h), du(x + h)) - Lp(x - h, u(x - h), du(x - h))) / (2 * h)
    return Lu(x, u(x), du(x)) - dLp

Lu = lambda x, u, p: 0.0 * x            # L = x p^2 / 2: L_u = 0, L_p = x p
Lp = lambda x, u, p: x * p
xs = np.array([1.5, 2.0, 2.5])
print(el_residual(Lu, Lp, np.log, lambda x: 1 / x, xs).round(6))           # [0. 0. 0.]: u = log x solves it
print(el_residual(Lu, Lp, lambda x: x, lambda x: 1 + 0 * x, xs).round(6))  # [-1. -1. -1.]: u = x does not

# (b) worked example: J[u] = int_0^1 (u'^2/2 - u) dx, u(0) = 0, u(1) free
J = lambda u, du: integrate.quad(lambda x: 0.5 * du(x) ** 2 - u(x), 0, 1)[0]
print(round(J(lambda x: x - x**2 / 2, lambda x: 1 - x), 6))   # -0.166667

# Ritz with phi_1 = x, phi_2 = x^2: J(c) = c^T K c / 2 - b^T c, K_jk = int phi_j' phi_k', b_j = int phi_j
phi = [lambda x: x, lambda x: x**2]
dphi = [lambda x: 1.0, lambda x: 2 * x]
K = np.array([[integrate.quad(lambda x: dphi[j](x) * dphi[k](x), 0, 1)[0] for k in range(2)] for j in range(2)])
b = np.array([integrate.quad(phi[j], 0, 1)[0] for j in range(2)])
c1 = b[0] / K[0, 0]
print(c1, round(-0.5 * c1 * b[0], 6))     # 0.5 -0.125   (one trial function)
c = np.linalg.solve(K, b)
print(c, round(-0.5 * c @ b, 6))         # [ 1.  -0.5] -0.166667   (two: exactly u*)

# (c) discrete energy with a free right end: unknowns U_1..U_N, U_0 = 0, trapezoid weights for f = 1
N = 100; h = 1 / N
A = (np.diag(np.full(N, 2.0)) - np.eye(N, k=1) - np.eye(N, k=-1)) / h
A[-1, -1] = 1 / h                        # the last unknown sits in one cell only
F = np.full(N, h); F[-1] = h / 2
U = np.linalg.solve(A, F)                # gradient of J_h = 0
x = np.arange(1, N + 1) * h
print(np.abs(U - (x - x**2 / 2)).max().round(6), round((U[-1] - U[-2]) / h, 4))   # 0.0 0.005 (= h/2 -> 0)

# (d) conditioning: energy loss (Hessian A) vs residual loss (Hessian A^T A), N = 100
n = 100; hh = 1 / (n + 1)
A1 = (2 * np.eye(n) - np.eye(n, k=1) - np.eye(n, k=-1)) / hh**2
print(f"{np.linalg.cond(A1):.0f} {np.linalg.cond(A1.T @ A1):.3g}")   # 4134 1.71e+07

In (a), the term $-u'$ in $-\frac{d}{dx}(xu') = -u' - xu''$ is the $L_{px}$ part of the total derivative. A code that uses only $L_{pp}u''$ would report a wrong residual. In (c) nobody imposed $u'(1) = 0$; the discrete slope at the free end is $h/2$, which goes to $0$ with $h$.

6. Where it shows up #

• Physics-informed ML. PINNs, Deep Ritz, energy-based losses for elasticity, and weak-form ("variational") PINNs. Weak-form PINNs test the residual against a finite set of functions $\eta$ after an integration by parts; for a problem with an energy, this is the first variation set to zero. Finite elements use the same weak form with piecewise polynomials. See the ML for science track.
• Seismic and weather models. The elastic wave equation follows from Hamilton's principle. Full-waveform inversion and 4D-Var minimize a misfit functional over an Earth model or an initial state. Their adjoint equations come from the first variation of a Lagrangian with a multiplier.
• Regularization of inverse problems and images. Tikhonov smoothing $\tfrac12\int|\nabla m|^2$ adds $-\Delta m$ to the Euler–Lagrange equation and smooths everything. Total variation $\int|\nabla m|$ adds $-\nabla\cdot(\nabla m/|\nabla m|)$ and keeps sharp edges (faults, fronts).
• Probability and information theory. Maximum-entropy distributions (the Gaussian for a given mean and variance) come from setting the functional derivative of the entropy $-\int p\log p$, plus Lagrange multipliers for the constraints, to zero. The mean-field variational inference update $q_j \propto \exp\big(\mathbb E_{-j}[\log p(x, z)]\big)$ comes from the functional derivative of the ELBO in $q_j$.
• Gradient flows. Besides the heat equation (section 3), the Fokker–Planck equation of Langevin dynamics is the gradient flow of the KL divergence in the Wasserstein metric (not derived here).

7. Common mistakes #

• A partial instead of a total derivative. $\frac{d}{dx}L_p$ includes $L_{px}$, $L_{pu}u'$ and $L_{pp}u''$. For $L = \tfrac12 x u'^2$: $\frac{d}{dx}(xu') = u' + xu''$, not $xu''$.
• Wrong variations. $\eta$ must be zero exactly where $u$ is prescribed. Imposing both $u(b) = \beta$ and the natural condition $L_p(b) = 0$ asks for too much: the problem usually has no solution.
• "Euler–Lagrange means minimum." It is necessary only. Check convexity, or the second variation. Actions in mechanics are stationary, not minimal.
• Reading the Ritz loss as an error. Its optimum is the minimum energy ($-\tfrac16$ in the worked example), not $0$. Measure the error against a reference solution, or look at the residual.
• Boundary conditions in Ritz. Essential conditions must be imposed. A penalty weight that is too small leaves them violated; one that is too large makes the loss ill-conditioned. Natural conditions need no term.
• Fixed quadrature points. The network can make the integrand small at those points only. Resample the points at every step (as in Deep Ritz) or use enough of them.
• ReLU networks in a PINN for a second-order PDE. $\Delta u_\theta = 0$ almost everywhere.
• Forgetting the $1/h$. $\partial J_h/\partial U_i \approx h\,\delta J/\delta u(x_i)$. Comparisons with the functional derivative, step sizes and convergence tests must use $(1/h)\,\partial J_h/\partial U$.

8. Check yourself #

1. Give the Euler–Lagrange equation for $L = \tfrac12 u'^2 + \tfrac12 k^2 u^2 - f(x)\,u$.
   Answer

   $L_u = k^2 u - f$ and $L_p = u'$, so $k^2 u - f - u'' = 0$, that is, $-u'' + k^2 u = f$ (a screened Poisson equation).
2. Minimize $J[u] = \int_0^1 (u'^2 + u^2)\,dx$ with $u(0) = 0$ and $u(1) = 1$. Give the minimizer and $J$, and compare with $u = x$.
   Answer

   $L_u = 2u$ and $L_p = 2u'$, so $2u - 2u'' = 0$, that is, $u'' = u$. With the boundary values, $u = \sinh x/\sinh 1$. Integration by parts with $u'' = u$ gives $J$ $= \int(u'^2 + uu'')\,dx$ $= [uu']_0^1$ $= u(1)\,u'(1)$ $= \cosh 1/\sinh 1$ $= \coth 1$ $\approx 1.3130$. For $u = x$: $J = 1 + \tfrac13 \approx 1.3333$, which is larger, as it must be ($L$ is convex).
3. $L = \tfrac12 x u'^2$ on $[1, e]$ with $u(1) = 0$ and $u(e) = 1$. Give the Euler–Lagrange equation and its solution.
   Answer

   $L_u = 0$ and $L_p = xu'$, so $-\frac{d}{dx}(xu') = 0$, that is, $-u' - xu'' = 0$. Then $xu' = C$, $u = C\log x + D$, and the boundary values give $u = \log x$.
4. $J[u] = \int_0^1 \big(\tfrac12 u'^2 - f u\big)\,dx$ with no boundary conditions at all. Which conditions does a minimizer satisfy, and what must $f$ satisfy for a minimizer to exist?
   Answer

   $-u'' = f$ with the natural conditions $u'(0) = u'(1) = 0$ at both ends. Integrating the equation, $\int_0^1 f\,dx = -[u']_0^1 = 0$. If $\int f \ne 0$, the constants $u = c$ give $J = -c\int f$, which goes to $-\infty$ for a suitable sign of $c$, so there is no minimizer. If $\int f = 0$, the minimizer is unique up to an added constant.
5. The 1-D Laplacian matrix on a grid has $\kappa(A) \approx 4.1 \cdot 10^3$. Plain gradient descent on the energy loss needs about $k$ steps. About how many does it need on the residual loss $\tfrac12\|AU - F\|^2$, and why?
   Answer

   About $\kappa(A)$ times more, so roughly $4 \cdot 10^3\,k$. The residual loss has Hessian $A^\top A$ with $\kappa(A^\top A) = \kappa(A)^2 \approx 1.7 \cdot 10^7$, and the number of gradient-descent steps grows like the condition number.

Recall cardsAnswer in your head first, then reveal. Add the cards to your review deck to see them again just before you would forget.

1. What is a functional? Give one example.

   Show answer

2. Define the first variation $\delta J[u; \eta]$.

   Show answer

3. State the fundamental lemma of the calculus of variations.

   Show answer

4. Write the Euler–Lagrange equation for $J[u] = \int_a^b L(x, u, u')\,dx$.

   Show answer

5. What condition does a minimizer satisfy at a free end $x = b$?

   Show answer

6. Which PDE is the Euler–Lagrange equation of the Dirichlet energy $\int_\Omega (\tfrac12|\nabla u|^2 - fu)\,dx$?

   Show answer

7. On a grid with spacing $h$, how is $\partial J_h/\partial U_i$ related to $\delta J/\delta u$?

   Show answer

8. For $-u'' = f$, which derivatives of the network does a Ritz loss need, and which does a PINN loss need?

   Show answer

9. Why is the residual loss $\tfrac12\|AU - F\|^2$ worse conditioned than the energy $\tfrac12U^\top AU - F^\top U$?

   Show answer

10. Gradient descent on the discretized Dirichlet energy is a time step of which equation?

    Show answer

Add these cards to my reviewOpen the review →

9. Interview angle #

• "Derive the Euler–Lagrange equation." Take the first variation, integrate by parts with $\eta(a) = \eta(b) = 0$, and apply the fundamental lemma. Add the natural boundary condition at a free end.
• "Prove that the shortest curve between two points is a straight line." $L_p = u'/\sqrt{1 + u'^2}$ is constant, so $u'$ is constant; convexity of $L$ makes the line the global minimum.
• "You must solve a Poisson equation with a neural network. Which loss do you use?" Ritz: first derivatives only, Neumann conditions for free, but it needs an energy and imposed Dirichlet conditions. PINN: any PDE, but second derivatives, a term for every condition, and worse conditioning.
• "Why do PINNs often fail to train?" The differential operator in the loss makes it ill-conditioned, the loss terms need weights, and networks fit low frequencies first. Remedies include curriculum training, causal (time-ordered) weighting, and boundary conditions built into the network.
• "What is the functional derivative of $\int p\log p\,dx$? Find the maximum-entropy density with a given mean and variance." $\log p + 1$. With Lagrange multipliers for $\int p = 1$, the mean and the variance, $\log p$ is a quadratic in $x$: the density is a Gaussian.

10. Practice #

The notebook for this lesson runs in Colab. It has four exercises, each with a check cell and a folded solution:

1. Write euler_lagrange_residual with a finite-difference total derivative and use it to test candidate solutions, including $L = \tfrac12 x u'^2$.
2. Write discrete_functional_gradient for the Dirichlet energy on a grid and compare $(1/h)\,\partial J_h/\partial U$ with $\delta J/\delta u$ for several $N$.
3. Solve the 1-D Poisson problem by minimizing the discrete energy and the discrete residual, with gradient descent and with conjugate gradients, and count the iterations.
4. Ritz with polynomial trial functions $x, x^2, \dots, x^K$ for a problem with a free end: watch the natural condition $u'(1) = 0$ appear without being imposed.

11. Watch and read #

Introduction to Variational Calculus - Deriving the Euler-Lagrange EquationGood Vibrations with Freeball · 25 min

Functionals, the variation $u + \epsilon\eta$, and a step-by-step derivation of the Euler–Lagrange equation with integration by parts.

Calculus of Variations ft. Flammable Mathsvcubingx · 21 min

A visual introduction: functionals, variations as small bumps, and the shortest path between two points.

Approximate Solutions - The Ritz MethodGood Vibrations with Freeball · 27 min

The Ritz method on a beam: pick trial functions, write the potential energy, and minimize over the coefficients.

Physics Informed Neural Networks (PINNs) [Physics Informed Machine Learning]Steve Brunton · 35 min

How a PINN puts the PDE residual and the boundary conditions into the loss, and where it struggles.

• Peter Olver, The Calculus of Variations, free lecture notes: §3 (the Euler–Lagrange equation; the fundamental lemma is Lemmas 3.3–3.4), §4 (boundary conditions, including natural boundary conditions), §5 (the second variation), §6 (several variables; Example 6.1 is the Dirichlet integral).
• E and Yu, The Deep Ritz Method: §2 (trial functions and Monte Carlo quadrature inside SGD) and §3.1 (the Poisson equation with a boundary penalty).
• Raissi, Perdikaris and Karniadakis, Physics Informed Deep Learning (Part I): the residual (PINN) loss with collocation points.
• Krishnapriyan, Gholami, Zhe, Kirby and Mahoney, Characterizing possible failure modes in physics-informed neural networks (NeurIPS 2021): §3 (failures for convection and reaction-diffusion), §4 (the PDE term makes the loss landscape hard to optimize), §5 (curriculum regularization).