M2-C2: Central-Difference Normals

Taught: 2026-04-30 Milestone: Phase 1 M2 — Pedestal SDF Result: PASS (after teacher format was rewritten mid-session; the conversational walkthrough version is below)

Why this concept matters here

M1 had a single sphere with an analytical normal: normalize(p). Free, exact, primitive-specific. M2 introduces composition (min(d_sphere, d_box)), which breaks analytical normals — there is no longer a single primitive's formula to reach for. This concept establishes the central-difference gradient as the universal normal computation: works for any SDF, including compositions, smin blends, displacements, anything.

The walkthrough

Step 1: Where you are right now — analytical normals for primitives in isolation

For a sphere centered at the origin, the surface is the set of points where length(p) - r == 0. If you grab any point on that surface and ask "which way is out?", the answer is just normalize(p) — the position vector itself, normalized. This is a closed-form, zero-cost answer that drops out of the parametric form of the sphere. For a box, you do a small dance with step() against the absolute-value position to figure out which face you're closest to, and you return that face's axis-aligned normal. Both of these are exact, both cost essentially nothing, and both are what you'd reach for in Houdini if you were authoring a primitive VOP. The catch — and it's a big one — is that each formula is married to its specific primitive's parametric form. A sphere normal formula has no idea what to do at a point on a box.

Step 2: The dispatch problem when you compose two SDFs with min

When you write d = min(sdSphere(p), sdBox(p)), the surface you're rendering is the union — a sphere with a box poking out (or vice versa). The march loop hits some surface point p_hit and says "give me a normal." Your first instinct is fine: figure out which primitive p_hit belongs to (whichever returned the smaller distance) and call that primitive's analytical normal function. This works at the cost of an extra dispatch — you have to re-evaluate both SDFs to know which one is "active," then branch into the right normal formula. Annoying, but tractable for two primitives. For N primitives it becomes O(N) work per normal, every shaded pixel. Already this is starting to smell wrong, but it's not yet broken.

Step 3: The seam problem — where union surfaces become ambiguous

Now consider the seam. Picture the sphere and box overlapping, and look at the curve where their surfaces meet. At a point exactly on that seam, both sdSphere(p) and sdBox(p) return zero — you're on both surfaces simultaneously. Which analytical formula do you use? The sphere's, which says "point radially outward"? The box's, which says "point along the nearest face axis"? They give different answers, and neither one is the actual normal of the union surface, because the union surface at the seam is a crease — its true normal is a discontinuous jump, but as you approach the seam from either side the analytical formulas give you stable but wrong answers because they don't know about the other primitive. You'll see a visible lighting tear along the seam. This is fixable in the union case if you're careful, but it's a warning shot.

Step 4: The smin destruction — where analytical normals die for good

smin doesn't just take the minimum, it blends the two distance fields in a small region near the seam, producing a smooth peanut shape. Read carefully: the surface that smin produces is not on the sphere, and not on the box. It's a brand-new surface that exists nowhere in either parent SDF's parametric definition. There is no "sphere formula" or "box formula" that describes where that blended surface lives, and so there is no analytical normal formula for it either. You can't dispatch your way out of this — neither primitive's normal function describes the actual shape you're shading. If you tried, the blend zone would light up wrong: smooth geometry, but with one of the two parent normals being faked across it, producing a visible discontinuity exactly where the eye is most sensitive. This is the wall. Analytical is over.

Step 5: The conceptual leap — SDF is a scalar field, and its gradient IS the normal

Forget primitives for a second. An SDF is a function f(p) that assigns a number to every point in 3D space — positive outside the surface, negative inside, zero on it. That's a scalar field, exactly like a density field or a temperature field in Houdini. Now ask the Houdini question: what's the gradient of a scalar field? It's a vector at every point, pointing in the direction the scalar increases fastest, with magnitude equal to the rate of increase.

Apply that here. The SDF increases as you move away from the surface (distance grows). So the gradient at any point in space points away from the nearest surface. At a point ON the surface, the gradient points directly outward — which is exactly the definition of the surface normal. This is true regardless of what the SDF is made of. Sphere, box, union, smin, fractal Mandelbulb, displaced terrain — if you can evaluate f(p), you can take its gradient, and the gradient is the normal. You stop asking "what's the formula for this shape's normal" and start asking "can I evaluate the SDF at p ± ε." That's the shift.

Step 6: How to actually compute the gradient when you don't have an analytical derivative

In closed form you'd write ∇f = (∂f/∂x, ∂f/∂y, ∂f/∂z). But once f is a tree of min/smin/transforms, taking that derivative symbolically is a nightmare and not what you want to do in a shader. Instead, you approximate each partial derivative with a central difference — sample f at p + ε and p - ε along the axis, take the difference, divide by 2ε. The 2ε divisor cancels out when you normalize, so in WGSL the typical form is:

fn calc_normal(p: vec3<f32>) -> vec3<f32> {
    let e = vec2<f32>(NORMAL_EPS, 0.0);
    return normalize(vec3<f32>(
        scene_sdf(p + e.xyy) - scene_sdf(p - e.xyy),
        scene_sdf(p + e.yxy) - scene_sdf(p - e.yxy),
        scene_sdf(p + e.yyx) - scene_sdf(p - e.yyx),
    ));
}

What this code is doing: six SDF evaluations, two per axis, taking the difference along each axis to estimate the gradient, then normalizing the resulting vector. Notice it makes zero assumption about what scene_sdf is made of — it just calls it. This is the universality. The same six lines work for one sphere, two primitives unioned, ten primitives smin'd together, or a noise-displaced surface.

Step 7: A worked numerical example so the arithmetic is concrete

Take the simplest possible case so you can verify against the analytical answer: unit sphere at origin, f(p) = length(p) - 1.0. Camera ray hits the surface at p_hit = (1, 0, 0) (the +X pole). Use ε = 0.001.

X-axis partial:

• f(1.001, 0, 0) = sqrt(1.001² + 0 + 0) - 1 = 1.001 - 1 = 0.001
• f(0.999, 0, 0) = sqrt(0.999² + 0 + 0) - 1 = 0.999 - 1 = -0.001
• difference: 0.001 - (-0.001) = 0.002

Y-axis partial:

• f(1, 0.001, 0) = sqrt(1 + 0.000001) - 1 ≈ 0.0000005
• f(1, -0.001, 0) — same answer by symmetry, ≈ 0.0000005
• difference: ≈ 0

Z-axis partial: identical to Y by symmetry, difference ≈ 0.

Raw gradient vector: (0.002, ~0, ~0). Normalize: (1, 0, 0).

That matches the analytical answer normalize(p) = (1, 0, 0) exactly. The reassurance here is that central differences don't trade accuracy for generality on cases analytical also handles — you get the same answer, just at the cost of 6 SDF evals instead of one normalize. That cost buys you the ability to handle smin and any future composition without rewriting your normal code.

Step 8: Cost model — when does this actually run, and what does it scale with

The normal calculation runs once per shaded pixel (at the hit point), not once per march step. For a 1280×720 surface with 70% coverage, that's roughly 645k normal computations per frame, each doing 6 SDF evals. The cost scales with the complexity of scene_sdf — if your scene SDF is a single sphere, six evals of length(p) - r is nothing. If your scene SDF is a smin tree of fifty primitives, six evals is fifty × six = 300 primitive distance evaluations per pixel, which starts to matter. The good news: it scales linearly with scene complexity, not with primitive count squared, and the analytical alternative for smin is literally impossible, so even an expensive numerical normal is the only game in town.

Step 9: Epsilon choice — why NORMAL_EPS = 0.001 and not something else

Epsilon is the only knob, and it's a tradeoff. Too small, and f(p+ε) - f(p-ε) becomes dominated by floating-point noise — you're subtracting two nearly-equal numbers, catastrophic cancellation, and the resulting "normal" jitters frame to frame. Too large, and you're sampling the SDF too far from the actual hit point, smoothing over fine geometric detail (a surface bump narrower than 2ε will be invisible to the gradient). For a unit-scale scene (camera ~3 units back, primitives ~1 unit in size), 0.001 is the standard choice — three orders of magnitude smaller than scene features, safely above f32 noise. M2 will hard-code const NORMAL_EPS = 0.001 in the shader. If you ever scale the scene way up or way down, this constant needs to scale with it.

The mental model in one sentence

A normal is the gradient of the distance field, and the gradient of any distance field — no matter how composed — can be measured by sampling it at six nearby points; analytical primitive normals are an optimization that stops working the instant you blend.

Explain-back question

Concrete scenario: sphere at origin (radius 1) and box at (1.5, 0, 0) with half-extents (0.5, 0.5, 0.5), composed with smin (smoothing radius say 0.2). The composed surface is a peanut. A march ray hits in the blend zone between the two primitives — the smooth waist where neither parent surface actually exists.

If you tried to use analytical normals (pick whichever primitive's SDF returned the smaller value, plug into that primitive's analytical normal formula), what would the rendered image look like, and why is central-difference the only thing that fixes it? Tie your answer to what the gradient is actually measuring at that blend-zone point.

User's answer

If you used analytical solution, firstly it would be impossible to correctly determine which is the right primitive to pick for normal calculation. Even if we were to pick one, the normals would be discrete and completely different depending on the primitive picked. The solution is to use gradient since SDF are continuous functions and we can use gradient computation for normal computation to solve this issue. So for each hit shading, we will be doing 6 computation per primitive!

Judgment

PASS. Got:

• Dispatch impossibility — picking which primitive's analytical formula to use is the first wall
• Discontinuous normals — even if you pick one, you get the wrong answer (and implicitly that this would show as visible lighting tears)
• Gradient because SDF is continuous — the conceptual leap from "primitive formulas" to "field gradient" is intact

Minor framing note: the cost is "6 calls to scene_sdf, each internally touching all N primitives → total = 6×N primitive evals," not "6 work per primitive." Same total math, cleaner mental model. Not load-bearing.