Surface Intersection of procedural shapes to find vertex normals?

If you have a vertex sitting on a procedural shape, like a perlin mountain/caves or a Sin Cos jellyfish ball, and you wish to sample some points around that vertex's position, not on the mesh itself, but on the function that made the mesh...

What is the Best method to do it??? I think I need 3 points adjacent to the vertex to triangulate the normal?

For example, surface formula "P" perlin(Xq,Yr)*perlin(Zs,t)*sqrt(xx*yy*zz)= 0.400

for the moment i can see these solutions:

A/generate some Lines adjacent to the vertex, and sample the intersection of the lines with formula P, (they have to be equidistant from the vertex?) and take 3 points in that fashion.

B/ make a larger of formula P, = to 0.401 (the surface will be slightly larger/above previous one) and the take the closest point on the second formula, to the vertex on the first one? seems ok but how do i calculate the closest point on the larger formula compared to the vertex?

C/make a sphere formula around the vertex, say 10 pixels diameter, and take three points from the intersection of the sphere and the formula P.

I can see these in theory but i think the maths is difficult! i got an E in advanced maths! is there any code/ theory on about normals from intersections of points with surface formulas? what would be the fastest?