Last modified: 2025-03-31 21:01:35
< 2025-03-30 2025-04-01 >Two issues here:
MeshNode has some weird leakageVoxelNode clips the bounding boxFor some reason the MeshNode leakage only appears when evaluating it for VoxelNode or
Mesher. I think this suggests it only happens at exact 0 coordinates or something.
I wonder if this is because Claude keeps trying to use P.lt() and that doesn't exist
so I always rewrite it to P.lte() on the basis that exact equality doesn't matter.
Actually it turns out my latest on this is not pushed and I am working on a different computer today. So, look at something else.
So the next thing to look at is making Peptide do automatic differentiation.
So each node type will give a function that returns a node that evaluates to its first derivative.
One important question is how we handle taking the first derivative of a function that yields a vector. I had kind of imagined that we would just return a vector showing the derivative with respect to p_x, p_y, and p_z, but Claude is telling me that you should use a Jacobian matrix?
OK, the issue is that the derivative of a scalar function is a vec3, but the derivative of
a vector function is a mat3.
This is annoying because it means Peptide will need more robust matrix support. Also you have the problem that you need to return a different type of derivative depending on whether the variable you're differentiating with respect to is a vector or a scalar.
Ideas are:
I think go for the second.
And how do I take the derivative of a matrix constant? Maybe instead of returning a vec3 Peptide node, I instead return 3 Peptide expressions? It will obviously (?) be less efficient to evaluate on GPU, but means I don't need everything to have an extra dimension.
One thing is that the derivative code will probably have a lot of common subtrees with the original code. Maybe that doesn't matter, because we'll calculate them in different places anyway.
There are some functions that I can't really work out how to differentiate, like step() and mod().
For now we're going with first derivative of step() is always 0, and first derivative of mod(a,b)
is equal to first derivative of a. Not sure how legit this is.
I'm having trouble getting tests of second derivative of vlength() to pass.
The issue is that the first derivative of vlength() involves dividing the vector by its length,
and if its length is 0 then that is no good, so Claude has clamped the length to be above some minimum value.
The issue is the first derivative is not defined for vec3(0,0,0). So does that mean the test is bad? No,
Claude doesn't seem to think so.
I'm going to try out Cursor agent mode to see if it can solve it.
No, it's no better. It seems to be the same as normal Cursor except it automatically searches for useful context. Maybe I don't find agent mode useful because I already know how to give Cursor the context it needs?
I could do with a way to dump out Peptide expressions so that I can read them.
It was because it was using P.max() to make the length division not divide by zero, and the derivative
of max() is discontinuous. Switching to a continuous version has fixed it!
However this kind of issue makes me think that automatic differentiation is actually a really hard problem and that this isn't going to work for me.
For the expression
max(vecX(p), vecY(p))
I'm getting first derivative is:
x = add(mul(step(sub(vecX(p), vecY(p)), 0), vecX(Vec3(1, 0, 0))), mul(sub(1, step(sub(vecX(p), vecY(p)), 0)), vecY(Vec3(1, 0, 0))))
y = add(mul(step(sub(vecX(p), vecY(p)), 0), vecX(Vec3(0, 1, 0))), mul(sub(1, step(sub(vecX(p), vecY(p)), 0)), vecY(Vec3(0, 1, 0))))
z = add(mul(step(sub(vecX(p), vecY(p)), 0), vecX(Vec3(0, 0, 1))), mul(sub(1, step(sub(vecX(p), vecY(p)), 0)), vecY(Vec3(0, 0, 1))))
Which simplifies down to:
x = step(p.x-p.y, 0)
y = 1 - step(p.x-p.y, 0)
z = 0
Is that right? I need to get my head around how step() is meant to work.
step(edge, x) = 1 if x > edge, 0 otherwise.
0 > a-b if b > a.
So step(a-b, 0) == step(a, b).
So I think the max() derivative may be using step() backwards? Yep, inverting the arguments fixes it.
So, for sake of argument, let's say we have working derivatives now. How do we use them for calculating normals?
It might be good to have a button to toggle between central differences and automatic differentiation, to see the difference. And forget edge detection for now because it's broken anyway.
Whoa! This is actually working!
I can't tell whether it is faster or slower than using central differences. Automatic differentiation seems to be ever so slightly faster, but less than 10% better. And they're not pixel-identical, but close enough to be fine. I expect if anything the central differences version is slightly wrong.
Weirdly, I think the shader compilation time is worse with central differences? Even though there is more code in
the automatic differentiation one? Is this because every call site has to have a copy-paste of the code, and central differences
has 6 extra call sites for map() but automatic differentiation only has 1 call site for map2()?
Yes it does compile faster with automatic differentiation! So that is good.
But there is some glitch if you enable radius on a box, the derivative becomes broken. So hopefully I can fix this,
otherwise can revert to central differences.
It's not to do with radius, it's just the existence of a box that breaks it. The normals on a box are broken.
Also: how is this working? I'm returning the first derivative... but that's not the same as a normal! Is this only working by accident? Is my lighting code all wrong? I think the lighting might only work by accident, my old central differences method is calculating the gradient rather than the normal. Good enough, don't worry about it for now.
So why is the normal of a box wrong?
I was trying to debug this but getting enormous expressions with lots of constant operations. Obviously the GLSL compiler can do constant folding perfectly well, but it's easy enough to add to Peptide so I may as well.
Maybe the issue is that because it's a gradient rather than a normal, it's 0 everywhere on the surface of a box, instead of pointing outwards? No, that's not it. The gradient has 3 values, you'd expect to have 0 in 2 of them and 1 in the 3rd.
Oh, maybe a gradient of a distance actually is the same as the normal of a surface.
Claude reckons the issue is that BoxNode uses vabs() which is
discontinuous at the surface, so doesn't give you good normals. Suggests
changing from abs(x) to sqrt(x^2 + eps), which is the same as
abs(x) when x is much larger than eps, but is smooth near 0.
But it doesn't make any difference.
I made the "show coordinates on hover" thing show the computed surface normal, and it... looks correct? Is it maybe right in JavaScript but wrong in GLSL?
In GLSL the start of the derivative code is:
v2 = p;
v0 = abs(v2);
v3 = u_Box_4_size;
f0 = 2.0000000000000000;
v1 = v3 / f0;
v4 = v0 - v1;
v1 = vec3(0.0000000000000000, 0.0000000000000000, 0.0000000000000000);
v0 = max(v4, v1);
f0 = length(v0);
f0 = f0 * f0;
f1 = 0.0000000000000000;
f0 = f0 + f1;
f0 = sqrt(f0);
v0 = v0 / f0;
And when we get to the line v0 = v0 / f0, for points on or inside
the cube we get vec3(0,0,0) / 0.
The Peptide code is:
let d = P.vsub(P.vabs(p), halfSize);
expr = P.add(P.vlength(P.vmax(d, P.vconst(new Vec3(0.0)))),
P.min(P.max(P.vecX(d), P.max(P.vecY(d), P.vecZ(d))), P.zero()));
and v4 is obviously d.
The JavaScript code looks equivalent so should also be
doing 0/0. But I made Vec3.div() log if the denominator is 0,
and it doesn't log anything.
OK!! Figured it out.
The JavaScript code halts ray-marching when it gets near the surface, and the GLSL code halts ray-marching when it gets "just below" the surface (because I wanted to make sure 0-sized surfaces disappear). If I make the GLSL code stop just before the surface then I get good normals on the box.
So the answer is that my normals are wrong on the inside of the box.
I think I consider that a bug to be honest. The distance field still has a gradient on the inside of the box, so the first derivative should still work.
So what is actually going on with that division? What things is it trying to divide?
Oh! It's the thing where vlength() uses sqrt(x^2 + eps) as
a smooth abs(). eps is chosen as 1e-20, but we
write everything with 16 decimal places in GLSL. So it comes out as
precisely 0. I increased that to 1e-10 and now is working, great
success.