A step-by-step, animated explainer

The Calculus of Variations

Ordinary calculus optimizes over numbers. This optimizes over entire functions. By the end you'll see it is the same machinery — gradient, zero-gradient, gradient descent — with the variable promoted from a point to a curve.

trial functions the extremal (least “cost”) a variation we test

The whole idea, in one sentence

A functional is just a function whose input is an entire function instead of a number; its functional derivative is the ordinary gradient with the index \(i\) replaced by a continuous coordinate \(x\) — so “set the gradient to zero” becomes “find the best-shaped curve.”

Everything below is built on one move, made early and never undone: stop thinking of a function \(y(x)\) as a rule, and start thinking of it as a point — a point with one coordinate for every value of \(x\). Once \(y\) is a point, a functional \(F[y]\) is an ordinary surface over that (infinite-dimensional) space, and all of calculus carries over. Seven steps, a live animation in every one, one promise: you'll watch functionals actually change, not just sliders move.

Step 1 — the move that makes everything work

A function is a point with infinitely many coordinates

Sample a function \(y(x)\) at \(N\) evenly spaced points on \([0,1]\). What you hold in your hand is now a vector \((y_1, y_2, \dots, y_N)\in\mathbb{R}^N\) — an ordinary multivariable input. A quantity built from it, \(F(y_1,\dots,y_N)\), is an ordinary multivariable function. Nothing exotic yet.

Now turn the dial: let \(N\to\infty\). The samples fuse into the continuous curve, the index \(i\) becomes the continuous coordinate \(x\), and \(F\) becomes a rule that eats the whole function and returns a number. We write it \(F[y]\) with square brackets and call it a functional. That is the entire conceptual leap — drag \(N\) and watch a vector become a function.

Figure 1From a point in \(\mathbb{R}^N\) to a curvelive · drag N

At small \(N\) the input is a point in \(\mathbb{R}^N\); push \(N\) up and the same data becomes a curve. A functional is what you get when you let this go all the way.

Step 2 — what a functional actually does

Feed it a whole curve, get back a single number

A functional collapses an entire function down to one number. Three workhorse examples, each an integral of the curve and possibly its slope:

\( \displaystyle A[y]=\int_0^1 y\,dx \qquad E[y]=\int_0^1 y^2\,dx \qquad L[y]=\int_0^1\!\sqrt{1+y'^2}\,dx \)

Press play and the curve below genuinely deforms — not a single knob turning, but the whole shape wandering through a family. Watch the output number breathe with it. Switch the functional and the same instant of curve produces a different number: that's the point — a functional is a rule from shapes to scalars.

Figure 2A whole curve collapses to a single numberlive · play & switch functional
A[y] = ∫ y dx
E[y] = ∫ y² dx
L[y] = ∫ √(1+y′²) dx

The shaded region is the integrand being accumulated. The bar is the output value on its natural range.

Step 3 — the gradient, reborn

The functional derivative: how sensitive the number is to a poke at \(x\)

In \(\mathbb{R}^N\), the gradient's \(i\)-th entry answers: if I nudge coordinate \(y_i\), how much does \(F\) change? Promote \(i\) to \(x\) and you get the functional derivative \(\frac{\delta F}{\delta y(x)}\): a whole function recording the sensitivity of \(F\) to a poke at each location \(x\).

We make this literal. Take \(F[y]=\int_0^1 y(x)^2\,dx\). Add a tiny localized bump \(\eta\) of strength \(\epsilon\) at a chosen spot, measure the change \(\Delta F\), and shrink \(\epsilon\). The linear part is the first variation \(\;\delta F=\int \frac{\delta F}{\delta y}\,\eta\,dx\). Then sweep the bump across the domain: the trail of sensitivities it leaves is the functional derivative — and for this \(F\) it traces exactly \(2y(x)\). Walk the four stages.

Figure 3Sweeping a bump traces out the functional derivativestaged · drag the bump

Drag the bump along the curve at any time. The dashed magenta curve is \(y+\epsilon\eta\); the gold trail is the swept sensitivity.

The formal version — first variation by direct expansion (article Eq. 9–11)

With \(F[y]=\int_0^1 y^2\,dx\) and variation \(\delta y=\epsilon\eta\):

\( F[y+\delta y]=\int_0^1 (y+\delta y)^2 dx = F[y] + \int_0^1 2y\,\delta y\,dx + \underbrace{\int_0^1 \delta y^2 dx}_{\to 0} \)

So \(\delta F = \int_0^1 2y(x)\,\delta y(x)\,dx\). Matching against \(\delta F=\int \frac{\delta F}{\delta y}\,\delta y\,dx\) reads off \(\frac{\delta F}{\delta y(x)} = 2y(x)\) — exactly the curve the sweep traces.

Step 4 — the shortcut

Euler–Lagrange: compute that sensitivity curve without sweeping

Sweeping bumps is honest but slow. Suppose the functional is the integral of one fixed recipe — call it the integrand \(L(x,y,y')\) (physicists call it the Lagrangian), an ordinary function of the position \(x\), the height \(y\), and the slope \(y'\): \(F[y]=\int_a^b L(x,y,y')\,dx\). (The square brackets in \(F[y]\) still mean “whole function in, one number out”; the parentheses in \(L(x,y,y')\) are an ordinary function of three numbers.) For any functional of this form, the Euler–Lagrange equation hands you the functional derivative in closed form:

\( \displaystyle \frac{\delta F}{\delta y(x)} \;=\; \frac{\partial L}{\partial y} \;-\; \frac{d}{dx}\frac{\partial L}{\partial y'} \)

That formula isn't handed down from above — it is forced by a single demand: at the optimal curve, nudging it in any admissible direction must not change \(F\) to first order. That is the exact function-space echo of \(\nabla F=0\). Here is the whole thing from scratch, in five moves; the engine is the same perturb-and-expand trick you used by hand in Step 3.

From first principles

Deriving Euler–Lagrange

Goal: find the function \(y(x)\) on \([a,b]\), with the endpoints pinned at \(y(a)=y_a\) and \(y(b)=y_b\), that makes \(F[y]=\int_a^b L(x,y,y')\,dx\) stationary.

Move 1 Perturb the winner

Suppose \(y(x)\) is the optimum. Jiggle it in some arbitrary direction \(\eta(x)\), scaled by a single knob \(\varepsilon\). Both the curve and its slope move:

\[ y(x)\;\longrightarrow\; y(x)+\varepsilon\,\eta(x), \qquad y'(x)\;\longrightarrow\; y'(x)+\varepsilon\,\eta'(x). \]

Every competitor must still hit the same two endpoints, so the variation is not allowed to move them: \(\eta(a)=\eta(b)=0\). Apart from that, \(\eta\) is completely free — remember this, it does all the work in Move 5.

Move 2 Collapse to one variable

Feed the jiggled curve into the functional. Since \(\varepsilon\) is now the only thing that varies, the functional becomes an ordinary single-variable function:

\[ \Phi(\varepsilon)\;=\;F[\,y+\varepsilon\eta\,]\;=\;\int_a^b L\big(x,\;y+\varepsilon\eta,\;y'+\varepsilon\eta'\big)\,dx. \]

\(y\) is stationary exactly when \(\Phi\) is stationary at \(\varepsilon=0\) — and it must be, for every direction \(\eta\). The infinite-dimensional problem just became first-year calculus: \(\Phi'(0)=0\).

Move 3 Differentiate under the integral

We need \(\Phi'(\varepsilon)=\dfrac{d}{d\varepsilon}\displaystyle\int_a^b L\big(x,\,y+\varepsilon\eta,\,y'+\varepsilon\eta'\big)\,dx\), and two things deserve care: moving the \(\tfrac{d}{d\varepsilon}\) inside, and the chain rule itself.

Move the derivative inside. The limits \(a,b\) are fixed, and with \(L\) continuously differentiable and \(y,\eta\in C^1[a,b]\) the integrand and its \(\varepsilon\)-partial are continuous on the compact strip \([a,b]\times(-\delta,\delta)\). Leibniz's rule then licenses the swap, with no boundary terms since the domain itself doesn't move:

\[ \Phi'(\varepsilon)=\int_a^b \frac{\partial}{\partial\varepsilon}\,L\big(x,\,y+\varepsilon\eta,\,y'+\varepsilon\eta'\big)\,dx. \]

Chain rule, slot by slot. Read \(L=L(x,u,v)\) as a function of three independent slots; along the perturbation \(u=y+\varepsilon\eta\) and \(v=y'+\varepsilon\eta'\), so only \(u,v\) carry \(\varepsilon\), with \(\partial u/\partial\varepsilon=\eta\) and \(\partial v/\partial\varepsilon=\eta'\). Crucially, \(\partial L/\partial y\) and \(\partial L/\partial y'\) mean the partials of \(L\) in its second and third slots (each holding the other two fixed) — not derivatives with respect to the function \(y\). So

\[ \frac{\partial}{\partial\varepsilon}\,L(x,u,v)=\frac{\partial L}{\partial y}\bigg|_{(x,u,v)}\!\eta+\frac{\partial L}{\partial y'}\bigg|_{(x,u,v)}\!\eta'. \]

Evaluate at \(\varepsilon=0\). There \(u=y(x)\) and \(v=y'(x)\), so the partials land on the unperturbed curve:

\[ \Phi'(0)=\int_a^b\!\left(\frac{\partial L}{\partial y}\bigg|_{(x,y,y')}\!\eta+\frac{\partial L}{\partial y'}\bigg|_{(x,y,y')}\!\eta'\right)dx. \]

This equality is an identity — true for every admissible \(\eta\). It defines the first variation \(\delta F(y;\eta)=\frac{d}{d\varepsilon}\big|_{0}\,F[y+\varepsilon\eta]\), the directional (Gâteaux) derivative of \(F\) at \(y\) along \(\eta\). Stationarity from Move 2 is the separate demand that this identity equal zero for every \(\eta\): \(\;\Phi'(0)=0\).

Move 4 Integrate by parts to free \(\eta\)

The term in \(\eta'\) is the nuisance: it ties the direction to its own slope. Integration by parts trades \(\eta'\) for \(\eta\):

\[ \int_a^b \frac{\partial L}{\partial y'}\,\eta'\,dx =\underbrace{\left[\frac{\partial L}{\partial y'}\,\eta\right]_a^b}_{=\,0} -\int_a^b \frac{d}{dx}\!\left(\frac{\partial L}{\partial y'}\right)\eta\,dx. \]

The boundary term vanishes because \(\eta(a)=\eta(b)=0\) (Move 1 earns its keep). Substitute back, and every appearance of the direction is now a single clean factor of \(\eta\):

\[ \Phi'(0)=\int_a^b\left(\frac{\partial L}{\partial y}-\frac{d}{dx}\frac{\partial L}{\partial y'}\right)\eta(x)\,dx=0. \]

Move 5 Kill \(\eta\) — the fundamental lemma

Read that last line as an inner product: \(\langle g,\eta\rangle=0\), where \(g(x)=\frac{\partial L}{\partial y}-\frac{d}{dx}\frac{\partial L}{\partial y'}\). It holds for every admissible direction \(\eta\). A fixed function that is orthogonal to all directions can only be the zero function: if \(g\) were, say, positive on any little sliver, just choose \(\eta\) to be a small positive bump sitting exactly on that sliver and zero elsewhere — the integral would come out strictly positive, a contradiction. Hence \(g\equiv 0\). It is the function-space twin of a plain fact in \(\mathbb{R}^n\): a vector perpendicular to every vector is \(\mathbf{0}\).

The bracket that survives is the functional derivative — the same per-point gradient density you swept out bump-by-bump in Step 3:

\[ \frac{\delta F}{\delta y(x)}=\frac{\partial L}{\partial y}-\frac{d}{dx}\frac{\partial L}{\partial y'} \]

Setting it to zero is nothing more exotic than \(\nabla F=0\), enforced at one continuous coordinate \(x\) at a time.

Take the bending energy \(E[y]=\int_0^1 \tfrac12\,y'^2\,dx\). Here \(\partial L/\partial y=0\) and \(\partial L/\partial y'=y'\), so \(\;\frac{\delta E}{\delta y}=-\,y''(x)\): the sensitivity at each point is minus the curvature there. Reshape the curve by dragging its control points and watch the bottom panel — the functional derivative — respond live. Where the curve is straight, it reads zero.

Figure 4Euler–Lagrange live: a curve and its sensitivity \(-y''\)live · drag the curve

Top: the trial curve \(y(x)\) (drag the gold control points). Bottom: its functional derivative \(\delta E/\delta y=-y''(x)\). Zero-crossings — the inflection points — are marked.

The formal version — Euler–Lagrange worked examples (article Eq. 14–15)

If \(L\) has no \(y'\): \(F[y]=\int_0^1 x^3 e^{-y}\,dx \Rightarrow \frac{\delta F}{\delta y}=\frac{\partial L}{\partial y}=-x^3 e^{-y}\).

With both terms: \(F[y]=\int_0^1 x^2 y^3 y'^4\,dx \Rightarrow \frac{\delta F}{\delta y}=3x^2y^2y'^4-\frac{d}{dx}\!\left(4x^2y^3y'^3\right)=3x^2y^2y'^4-8xy^3y'^3.\)

The recipe never changes: differentiate \(L\) treating \(y\) and \(y'\) as independent, then subtract the total-\(x\) derivative of the \(y'\)-part.

Step 5 — the payoff

Set it to zero: gradient descent that optimizes a shape

A point is stationary when its gradient is zero. A function is stationary when its functional derivative is zero everywhere. And just as we slide downhill against \(-\nabla F\) to reach a minimum, we can flow a whole curve against \(-\frac{\delta F}{\delta y}\) to reach the extremal:

\( \displaystyle \frac{\partial y}{\partial t} \;=\; -\frac{\delta E}{\delta y} \;=\; y''(x) \)

That is gradient descent in function space — and for bending energy it is literally the heat equation. Pin the endpoints, start from any wiggle, press play: the curve diffuses, the energy falls monotonically, the functional derivative vanishes, and what's left is the straight line. You just watched “the shortest path is a line” happen as an optimization. Drag the endpoints and re-run.

Figure 5Gradient flow relaxes any wiggle to the straight linelive sim · drag the endpoints

Teal: the evolving curve. Faint trail: where it started. Gold dashed: the straight-line extremal it's converging to. Drag either endpoint.

The formal version — the article's classic: arc length ⇒ straight line (Eq. 16–19)

For true arc length \(L[y]=\int_a^b\sqrt{1+y'^2}\,dx\): \(\partial L/\partial y=0\) and \(\partial L/\partial y'=\frac{y'}{\sqrt{1+y'^2}}\). Euler–Lagrange gives

\( \frac{d}{dx}\frac{y'}{\sqrt{1+y'^2}}=0 \;\Rightarrow\; \frac{y'}{\sqrt{1+y'^2}}=C \;\Rightarrow\; y'=\text{const}=A \;\Rightarrow\; y=Ax+B. \)

A straight line — exactly the extremal the bending-energy flow above relaxes to. Two different functionals, same stationary curve; we animate bending energy because its flow (the heat equation) is numerically gentle and visibly diffusive.

Step 6 — the cost decides the shape

Change the functional again: the most even route across a landscape

The straight line was optimal because the cost was bending. Change the cost and the optimal shape changes with it — that is the whole point. So make it physical: a hiking route \(y(x)\) across a landscape, from a fixed start to a fixed end, that wants to be as even as possible — gaining and losing as little height as it can. The evenest imaginable climb is a constant grade: a straight line in altitude from the start's height \(z_0\) to the end's height \(z_1\). So penalize the squared deviation of the trail's altitude from that steady grade \(\bar z(x)=z_0+(z_1-z_0)x\):

\( \displaystyle F[y]=\int_0^1 \Big[\big(h(x,y(x))-\bar z(x)\big)^2+\tfrac{\lambda}{2}\,y'^2\Big]\,dx \)

Euler–Lagrange gives \(2\,(h-\bar z)\,\partial_y h=\lambda\,y''\), and flowing against \(-\,\delta F/\delta y\) settles any starting trail onto the terrain that best tracks the grade. When the two ends sit at the same height, \(\bar z\) is flat and the route becomes a level contour — climbing nothing, weaving between the hills and the sinks. When the ends sit at different heights, \(\bar z\) tilts and the same principle picks the steady descent down a valley floor — all drop, no wasted up-and-down. Swap the surface under the figure to see both; drag to view either from any angle.

Figure 6The most even route, on a surface you can swap3d · drag to rotate
surface Level traverse

Use the ‹ surface › box to swap the landscape calculus of variations runs on. Level traverse: a winding zero-level corridor between two ends at the same height — the even route is a flat contour. Valley descent: a winding valley whose floor drops from a high start to a low end — now the even route is the steady, constant-grade descent along the floor. Drag to orbit; press Most even trail to watch the route relax. Bottom: altitude \(h\) along the trail (white) against the steady grade (gold) it is trying to match.

The formal version — Euler–Lagrange for the grade-deviation functional, and the \(\lambda\to0\) iso-grade limit

With integrand \(L=(h-\bar z)^2+\tfrac{\lambda}{2}y'^2\) (and \(\bar z(x)=z_0+(z_1-z_0)x\) a fixed line, since the endpoints are pinned), the pieces are \(\partial L/\partial y=2(h-\bar z)h_y\) and \(\partial L/\partial y'=\lambda y'\). Euler–Lagrange, \(\;\frac{\delta F}{\delta y}=2(h-\bar z)h_y-\lambda y''=0\), i.e. \(\lambda y''=2(h-\bar z)h_y\). The flow \(y_t=-2(h-\bar z)h_y+\lambda y''\) is what the figure integrates: the reaction term pushes each point toward terrain whose altitude equals the grade (it vanishes where \(h=\bar z\)); \(\lambda y''\) lightly smooths.

Let \(\lambda\to0\) and the condition is purely \(h(x,y)=\bar z(x)\): the trail rides the terrain's iso-grade set. If \(z_0=z_1\) that is a contour line — constant altitude, zero gain or loss. If \(z_0\neq z_1\) it is the locus descending at exactly the average grade, so the trail's total climb-plus-descent equals the unavoidable net drop \(|z_1-z_0|\) and not a metre more — the steady valley descent. Same machinery as Step 5's bending energy; a different integrand; a different optimal shape.

Step 7 — where it pays off in probability

Constraints & maximum entropy: where the Gaussian comes from

Constrained optimization over functions works exactly like Lagrange multipliers over numbers: add \(-\lambda\,(\text{constraint})\) to the integrand and apply Euler–Lagrange. Maximize the entropy \(H[f]=-\int f\log f\,dx\) of a density subject to what you know about it. Euler–Lagrange on the Lagrangian immediately yields the stationary form

\( \displaystyle f(x)=e^{-1-\lambda_0-\lambda_1 x-\lambda_2 x^2} \)

— an exponential family. The multipliers are then tuned so the constraints hold. Know only the support \([a,b]\)? Only \(\lambda_0\) is active and \(f\) flattens to the uniform distribution. Know the mean and variance on the whole line? All three multipliers engage and \(f\) becomes the Gaussian. Press a button and watch the density climb in entropy while the constraint meters lock onto their targets.

Figure 7Constraints select the maximum-entropy densitylive sim · choose the constraints

Teal fill: the current density \(f(x)\), always of the exponential form above. Gold dashed: the analytic maximum-entropy target. Meters show the active constraints converging.

Putting it together

Every step was the same calculus you already know, with one substitution made on the first line — the variable is a function, the index is continuous. Line them up:

Ordinary calculus on \(\mathbb{R}^N\)Variational calculus on functions
point \((y_1,\dots,y_N)\)function \(y(x)\)
function \(F(y_1,\dots,y_N)\)functional \(F[y]\)
gradient \(\partial F/\partial y_i\)functional derivative \(\delta F/\delta y(x)\)
\(dF=\sum_i \frac{\partial F}{\partial y_i}\,dy_i\)\(\delta F=\int \frac{\delta F}{\delta y}\,\delta y\,dx\)
set \(\nabla F=0\) for the best pointEuler–Lagrange \(=0\) for the best function
Lagrange multipliersmultipliers / multiplier functions
gradient descent \(\dot y=-\nabla F\)gradient flow \(\partial_t y=-\delta F/\delta y\)

The kernel holds up: a functional is a function of a function, its derivative is the gradient with a continuous index, and zeroing that derivative finds the best-shaped curve. The dressing — \(\delta\) instead of \(\partial\), \(\int\) instead of \(\sum\), “extremal” instead of “minimizer” — is just formal clothes over machinery you've used since first-year calculus.

Variational inference / VAEs

Maximizing the ELBO is optimizing a functional over a density \(q\). The “variational” in variational Bayes is this exact calculus.

The exponential family

Maximum entropy under moment constraints yields \(e^{\sum\lambda_k T_k(x)}\) — Gaussian, exponential, Boltzmann — straight from Euler–Lagrange.

Least action in physics

Trajectories extremize \(\int L\,dt\). Euler–Lagrange on the action is Newton's laws; geodesics minimize path length.

Diffusion as gradient flow

The heat equation you saw in Step 5 is gradient flow of an energy — the same lens behind score-based diffusion and many smoothing methods.