Calculus for Simulations: Cheat Sheet
Every physics simulation ultimately reduces to differential equations — and solving them requires a toolbox of calculus operators, coordinate systems, numerical schemes, and rotation math. This page collects the essentials in one place.
1. Vector Calculus Operators
Four fundamental differential operators on scalar fields f(x,y,z) and vector fields F = (F_x, F_y, F_z):
∇f = (∂f/∂x, ∂f/∂y, ∂f/∂z)
Divergence (vector → scalar): net outflow per unit volume
∇·F = ∂F_x/∂x + ∂F_y/∂y + ∂F_z/∂z
Curl (vector → vector): local rotation
∇×F = (∂F_z/∂y − ∂F_y/∂z, ∂F_x/∂z − ∂F_z/∂x, ∂F_y/∂x − ∂F_x/∂y)
Laplacian (scalar → scalar): curvature / diffusion
∇²f = ∂²f/∂x² + ∂²f/∂y² + ∂²f/∂z²
Vector Laplacian: ∇²F = (∇²F_x, ∇²F_y, ∇²F_z) (Cartesian only!)
2. Coordinate Systems
Cylindrical (r, φ, z)
∇f = (∂f/∂r, (1/r)∂f/∂φ, ∂f/∂z)
∇²f = (1/r)∂/∂r(r ∂f/∂r) + (1/r²)∂²f/∂φ² + ∂²f/∂z²
dV = r dr dφ dz
Spherical (r, θ, φ)
∇f = (∂f/∂r, (1/r)∂f/∂θ, (1/r sinθ)∂f/∂φ)
∇²f = (1/r²)∂/∂r(r²∂f/∂r) + (1/r²sinθ)∂/∂θ(sinθ ∂f/∂θ) + (1/r²sin²θ)∂²f/∂φ²
dV = r² sinθ dr dθ dφ
3. Integral Theorems
∫∫∫_V ∇·F dV = ∮∮_S F · n̂ dA
Stokes' Theorem: surface integral of curl = line integral around boundary
∫∫_S (∇×F) · n̂ dA = ∮_C F · dl
Green's Theorems:
First: ∫∫∫ (φ∇²ψ + ∇φ·∇ψ) dV = ∮∮ φ(∇ψ·n̂) dA
Second: ∫∫∫ (φ∇²ψ − ψ∇²φ) dV = ∮∮ (φ∇ψ − ψ∇φ)·n̂ dA
These theorems underpin conservation laws in simulations: mass conservation (divergence of velocity), electromagnetic boundary conditions (Stokes for Faraday's law), and finite-volume methods (Gauss's theorem converts PDEs to surface fluxes).
4. PDE Classification
Second-order linear PDEs in 2D: Au_xx + 2Bu_xy + Cu_yy + (lower order) = 0. Classification depends on the discriminant B² − AC:
| Type | Discriminant | Example PDE | Physics | Numerical method |
|---|---|---|---|---|
| Elliptic | B² − AC < 0 | ∇²u = f (Poisson) | Steady-state heat, electrostatics | Jacobi, Gauss-Seidel, multigrid |
| Parabolic | B² − AC = 0 | ∂u/∂t = α∇²u (diffusion) | Heat conduction, diffusion | FTCS, Crank-Nicolson, ADI |
| Hyperbolic | B² − AC > 0 | ∂²u/∂t² = c²∇²u (wave) | Waves, acoustics, advection | Leapfrog, upwind, Lax-Wendroff |
5. Finite Difference Schemes
Approximate derivatives on a discrete grid with spacing Δx (space) and Δt (time):
∂u/∂x ≈ (u_{i+1} − u_{i−1}) / (2Δx) + O(Δx²)
Second derivative (central, 2nd order):
∂²u/∂x² ≈ (u_{i+1} − 2u_i + u_{i−1}) / Δx² + O(Δx²)
2D Laplacian (5-point stencil):
∇²u ≈ (u_{i+1,j} + u_{i−1,j} + u_{i,j+1} + u_{i,j−1} − 4u_{i,j}) / Δx²
Time Integration Schemes for ∂u/∂t = α∇²u
| Scheme | Formula | Stability | Accuracy |
|---|---|---|---|
| FTCS (Forward Time, Central Space) | u_i^{n+1} = u_i^n + αΔt/Δx² · (u_{i+1}^n − 2u_i^n + u_{i−1}^n) | Conditional: r = αΔt/Δx² ≤ 0.5 | O(Δt, Δx²) |
| Crank-Nicolson | Average of FTCS at n and n+1 (implicit) | Unconditionally stable | O(Δt², Δx²) |
| ADI (Alternating Direction Implicit) | Split 2D into two 1D implicit steps (x then y) | Unconditionally stable | O(Δt², Δx²) |
c·Δt/Δx ≤ 1 (Courant number C ≤ 1)
Von Neumann stability (parabolic):
r = αΔt/Δx² ≤ 0.5 for FTCS
Leapfrog (2nd order for wave eq):
u_i^{n+1} = 2u_i^n − u_i^{n−1} + C²(u_{i+1}^n − 2u_i^n + u_{i−1}^n)
where C = cΔt/Δx
6. Complex Numbers and Phasors
Modulus: |z| = √(a² + b²) = r
Argument: arg(z) = atan2(b, a) = θ
Conjugate: z* = a − bi
Inverse: 1/z = z*/|z|²
Multiplication: z₁·z₂ = r₁r₂ · e^{i(θ₁+θ₂)} (multiply moduli, add angles)
Division: z₁/z₂ = (r₁/r₂) · e^{i(θ₁−θ₂)}
Power: z^n = r^n · e^{inθ} (De Moivre's theorem)
Euler's identity: e^{iπ} + 1 = 0
Phasors represent sinusoidal steady-state as complex amplitudes: v(t) = Re[V̂ · e^{jωt}] where V̂ = V₀e^{jφ}. This converts differential equations in time to algebraic equations in ω — the basis of AC circuit analysis, transfer functions, and frequency-domain signal processing.
7. Quaternion Math
Quaternions are the standard rotation representation in 3D simulations — they avoid gimbal lock, interpolate smoothly (slerp), and compose rotations efficiently.
where i² = j² = k² = ijk = −1
Conjugate: q* = w − xi − yj − zk
Norm: |q| = √(w² + x² + y² + z²)
Inverse: q⁻¹ = q*/|q|² (for unit quaternion: q⁻¹ = q*)
Multiplication (Hamilton product):
(q₁q₂).w = w₁w₂ − x₁x₂ − y₁y₂ − z₁z₂
(q₁q₂).x = w₁x₂ + x₁w₂ + y₁z₂ − z₁y₂
(q₁q₂).y = w₁y₂ − x₁z₂ + y₁w₂ + z₁x₂
(q₁q₂).z = w₁z₂ + x₁y₂ − y₁x₂ + z₁w₂
NOT commutative: q₁q₂ ≠ q₂q₁ (just like rotations)
Axis-Angle ↔ Quaternion
q = cos(θ/2) + sin(θ/2)(n_x·i + n_y·j + n_z·k)
Inverse: θ = 2·acos(w), n̂ = (x,y,z)/sin(θ/2)
Rotate vector v by quaternion q:
v' = q · (0 + v) · q* (sandwich product)
Slerp (Spherical Linear Interpolation)
Explicit formula (avoiding pow):
Ω = acos(q₀ · q₁) (dot product of 4-vectors)
slerp(t) = (sin((1−t)Ω)/sinΩ) · q₀ + (sin(tΩ)/sinΩ) · q₁
If q₀ · q₁ < 0: negate one quaternion to take the shorter arc
If Ω ≈ 0: fall back to nlerp (normalised linear interpolation)
Exponential and Logarithm Map
exp(q) = (cos|v|, v̂·sin|v|)
log(q) for unit quaternion q = (w, v):
log(q) = (0, v̂·acos(w))
Usage: exp/log allow quaternion "velocity" (angular velocity ω → q(t) = exp(½ωt)) and smooth interpolation in the tangent space.
8. Key Identities
∇×(∇f) = 0 (curl of gradient is zero)
∇·(∇×F) = 0 (divergence of curl is zero)
∇×(∇×F) = ∇(∇·F) − ∇²F (vector triple product)
∇·(fF) = f(∇·F) + F·(∇f) (product rule for divergence)
∇×(fF) = f(∇×F) + (∇f)×F (product rule for curl)
Integration by parts (1D):
∫_a^b u dv = [uv]_a^b − ∫_a^b v du
Taylor expansion:
f(x+h) = f(x) + hf'(x) + (h²/2!)f''(x) + (h³/3!)f'''(x) + O(h⁴)
Useful limits:
lim_{x→0} sin(x)/x = 1
lim_{n→∞} (1 + 1/n)^n = e
Derivative of rotation matrix:
dR/dt = ω× · R (ω× = skew-symmetric matrix of angular velocity)
Reynolds transport theorem (moving control volume):
d/dt ∫∫∫_V(t) f dV = ∫∫∫_V ∂f/∂t dV + ∮∮_S f (v·n̂) dA