Perturbation Theory Through the Lens of Calculus
How small corrections connect deep fractal zooms, pendulums, quantum physics, and better engineering decisions
By Michael Stebel
In my earlier articles on perturbation theory, the Taylor series, and the Infinity Engine, I explored how a carefully computed reference orbit can make deep Mandelbrot rendering practical.
That raises a broader question: why does the same mathematical idea also help us understand an atom in a magnetic field, a swinging pendulum, or a spacecraft’s trajectory?
The connection is calculus. When a problem changes slightly, derivatives can tell us how its solution changes. Perturbation theory turns that observation into a problem-solving method: start with something understood, introduce a controlled change, and calculate the corrections.
The interesting part is that we can often calculate those corrections even when we cannot solve the full modified problem in a convenient closed form.
Start with the derivative
Suppose we know a function at a point x. For a small displacement h, its first derivative gives a local prediction:
f(x + h) ≈ f(x) + h f′(x)
The derivative is the rate of change at the starting point. Multiply it by the displacement to get our first correction. The second derivative adds information about curvature:
f(x + h) ≈ f(x) + h f′(x) + ½h² f″(x)
For example, near x = 0, the exponential function becomes:
eʰ ≈ 1 + h + ½h²
At h = 0.1, the exact value is approximately 1.105171. The linear prediction is 1.1; adding the quadratic correction gives 1.105. One extra term reduces the absolute error by about a factor of 30.
This is the familiar calculus picture: a tangent line captures slope, and the next term begins to capture bending. A finite Taylor approximation also has a remainder. Its accuracy depends on the displacement and the derivatives in the region we use.
Figure 1. The approximations agree near the expansion point. Their errors grow as we move away; the quadratic term captures curvature that the tangent line misses.
Expand the solution itself
In a physics or engineering problem, we may not have an explicit function to expand. We may have an equation whose solution is unknown.
Introduce a dimensionless parameter ε, pronounced “epsilon,” that controls a departure from a simpler problem. For a regular perturbation expansion, write the unknown solution as:
y(ε) = y₀ + εy₁ + ε²y₂ + …
Here y₀ solves the original problem. The next terms describe how that solution responds as the modification is introduced. If y depends analytically on ε near zero, these coefficients are Taylor coefficients:
y₁ = y′(0), y₂ = ½y″(0)
We now take the derivative with respect to the parameter controlling the problem. When the solution is a trajectory, the coefficients can themselves be functions of time. MIT’s classical mechanics notes develop this approach for dynamical systems.
A small algebraic example makes the method concrete:
y + εy² = 1
At ε = 0, the solution is y₀ = 1. Substitute the proposed series into the equation and equate coefficients of like powers of ε. The constant terms give y₀ = 1. The terms proportional to ε give y₁ + 1 = 0. The terms proportional to ε² give y₂ + 2y₁ = 0.
So:
y(ε) ≈ 1 − ε + 2ε²
At ε = 0.05, this predicts 0.955. The exact root connected continuously to y = 1 is approximately 0.954451.
We could solve this particular quadratic directly. For example, the same procedure extends to differential equations where an exact answer is much harder to obtain.
Implicit differentiation gives another view. Differentiating the original equation with respect to ε produces:
y′ + y² + 2εyy′ = 0
At ε = 0 and y = 1, we recover y′(0) = −1. The perturbation coefficient measures the solution’s sensitivity to the added term.
The bridge back to Mandelbrot rendering
For the Mandelbrot iteration, choose a reference parameter C and its reference orbit Zₙ:
Zₙ₊₁ = Zₙ² + C
A nearby pixel uses c = C + δc, and its orbit is zₙ = Zₙ + δzₙ. Substitute these expressions and cancel the reference equation. What remains is:
δzₙ₊₁ = 2Zₙδzₙ + (δzₙ)² + δc
The derivative of z² is 2z, which explains the linear term. The second derivative supplies the quadratic term. Because the map is quadratic, this one-step expansion has no higher powers of δz.
This sharpens an important point from my earlier Taylor-series explanation: the full recurrence above is an exact algebraic identity. Dropping (δzₙ)² produces a first-order approximation. Keeping it preserves the original quadratic dynamics, subject to numerical arithmetic errors.
Get Michael Stebel’s stories in your inbox
Join Medium for free to get updates from this writer.
Mandelbrot Metal’s Infinity Engine retains that quadratic term. Its published architecture shares an accurately computed reference orbit and tracks pixel-relative deviations. This can reduce the cost of repeated high-precision work, but it still requires precision management and recovery when a numerical path becomes unreliable. I cover the engineering details in my Infinity Engine article.
Figure 2. Keeping the quadratic term makes the rearrangement algebraically exact. It does not eliminate floating-point error or make every nearby orbit easy to compute.
A pendulum reveals the physics
The ideal pendulum follows a nonlinear differential equation:
θ̈ + (g/L) sin θ = 0
Here θ is the angle, θ̈ is its second derivative with respect to time, g is gravitational acceleration, and L is the pendulum length.
For small swings, sin θ ≈ θ. That gives simple harmonic motion and the familiar period:
T₀ = 2π√(L/g)
The model predicts a period independent of amplitude. A larger swing takes longer, however, because the sine function departs from its tangent line.
Perturbation theory captures the leading correction:
T ≈ T₀(1 + θ₀²/16)
The amplitude θ₀ must be in radians. At a 60-degree amplitude, the correction predicts a period about 6.85% longer than T₀; the exact ideal-pendulum result is about 7.32% longer. The approximation is much closer than the amplitude-independent model, while its remaining error is visible.
The first amplitude correction is quadratic. Reversing the direction of the initial displacement leaves the period unchanged, so a term linear in the signed amplitude would violate that symmetry. Nelson and Olsson’s pendulum paper derives the period expansion.
Figure 3. The correction explains why larger swings take longer. The curves describe an ideal, undamped pendulum; the horizontal axis uses degrees, while the calculation uses radians.
Atoms and spacecraft
In quantum mechanics, the known problem might be an atom without an external field. Adding a weak interaction changes the Hamiltonian, the operator describing the system’s energy:
Ĥ(ε) = H₀ + εV
For a nondegenerate energy level, the first-order energy correction is:
Eₙ(ε) ≈ Eₙ⁽⁰⁾ + ε⟨ψₙ⁽⁰⁾|V|ψₙ⁽⁰⁾⟩
The bracketed expression is the expectation value of the added interaction in the original state. It gives the first-order energy shift using that state. If symmetry makes this contribution zero, a higher-order term may supply the leading shift.
This approach helps describe weak-field energy shifts and departures from ideal harmonic molecular vibrations. Its validity depends on the interaction and relevant energy-level separations; degenerate or nearly degenerate levels need special care. MIT’s quantum mechanics notes explain the method and its limits.
For a spacecraft, an idealized two-body orbit offers another useful starting point. Real trajectories also respond to other gravitating bodies, solar radiation pressure, and, where relevant, atmospheric drag. NASA’s trajectory-design work accounts for these additional effects.
Analytical perturbation methods can explain departures from an ideal orbit, while mission software also uses numerical integration. The shared idea is to understand the baseline and then account for what changes it. A small acceleration can produce a substantial trajectory difference when it acts for long enough.
Beyond physics
The same calculus supports practical sensitivity analysis. Consider a model whose output depends on several inputs. A small input change gives the first-order prediction:
Δoutput ≈ J Δinput
J is the Jacobian: a matrix collecting the relevant partial derivatives. It tells us which combinations of inputs have the greatest local effect.
For a thermal or structural model, this helps estimate how changes in material properties or geometry affect computed temperature or displacement. For uncertainty analysis, it can translate small input errors into approximate output errors. These are local estimates, and their usefulness depends on the model staying in the same operating regime.
Optimization adds another layer: the output may be the best design a solver finds. How should that design change when demand, cost, or a constraint changes slightly? Differentiating the solution map addresses that question. Agrawal and Boyd demonstrate this for optimization problems used in queueing-system design and machine learning. Their paper directly links parameter perturbations to solution sensitivity.
How small is small enough
The size of ε alone cannot answer that. A small parameter multiplied by a large sensitivity may produce a large correction. A weak effect can also accumulate over time.
Check whether successive retained terms decrease in the range being used, then compare against a more accurate calculation or measurement. A smaller next term is useful evidence, but it does not guarantee a universal error bound.
Some perturbation series diverge but remain useful asymptotic expansions: the first few terms can improve the answer, while adding too many eventually makes it worse. Other problems are singular, so an ordinary expansion fails in part of the domain. MIT’s mechanics notes discuss growing time-dependent corrections, and its quantum notes show that useful perturbation series need not converge.
What keeps drawing me back to this mathematics is how much it reveals about a system. The correction tells us which feature of the original model matters, how strongly the system responds, and where our current explanation begins to fail.
That is the connection between a derivative in calculus and the mathematics behind a deep fractal zoom: understanding change well enough to calculate it, and checking how far that understanding remains reliable.
Resources that I recommend
- Asymptotics and perturbation methods — Prof. Steven Strogatz — https://youtube.com/playlist?list=PL5EH0ZJ7V0jV7kMYvPcZ7F9oaf_YAlfbI&si=aAHk30pWnaulgsTa
- MIT Paper on perturbation theory in quantum mechanics — https://ocw.mit.edu/courses/8-06-quantum-physics-iii-spring-2018/a0889c5ca8a479c3e56c544d646fb770_MIT8_06S18ch1.pdf
- Agrawal and Boyd paper — https://web.stanford.edu/~boyd/papers/pdf/diff_llcvx.pdf
- NASA trajectory and design work — https://www.nasa.gov/ames-engineering/spaceflight-division/flight-dynamics/trajectory-design/
Michael









