So I have a vector function f(x,y,z) -> (P.x,P.y,P.z).

I know I have to calculate the normal analytically as I’m displacing the vertex in the shader. I just have no idea how to go about it. I’ve tried getting the Jacobian
in Maple and Matlab, but aren’t getting the right results. Could someone with a little Calculus knowledge guide me in the right direction…?

IOW, change cos() to sin(), multiply by frequency * directions[i].x (or .y), and negate (although that doesn’t actually matter because (-A)×(-B)=A×B).

Obviously, there’s a lot of potential to remove common subexpressions from that.

Also, for large numbers of waves, it’s common to use an inverse FFT rather than explicitly evaluating and summing sine waves. Sampling m waves at n points requires O(mn) time, while a FFT can do n waves at n points in O(nlog(n)).

BTW, I just noticed this, which may be a source of misunderstanding:

No you don’t; you have a function f(u,v) -> (P.x,P.y,P.z).

The difference doesn’t really matter much here; it just means that the Jacobian won’t be square (it will be 3x2). Although it may have caused some confusion if you were expecting the parameter (pos) to be 3D.