The common method to calculate vertex normals for a mesh is to use the neighbour vertices to determine the orientation of the surface with subtraction of the positions and the cross product to get the normal vector.
Calculating the vertex normals for a mesh in a shader is usually impossible, because there is no access to the neighbours. So, if you have a mesh with no normals you have to calculate the normals on the CPU side.
But all this changes to a quite easy task, if you are working with dynamic meshes and you have a formula to describe the mesh in the following form: f(u, v) = (x, y, z). With that formula you are able to reach every point on the mesh surface, hence its easy to do all in the shader. Assuming the input mesh is a grid in the xy-plane, we first scale and offset the grid to get the input space we want, then we calculate the new position of the vertex by using the x and y coordinates of the current vertex as u and v in the formula. After that we get two neighbours by calculating two other positions with a small offset in u and v direction. A pseudo vertex shader code could look like this:
u = pos.x * scale + offset ;
v = pos.y * scale + offset ;
pos = f(u, v);
neighbour1 = f(u+smallvalue, v);
neighbour2 = f(u, v+smallvalue);
tangent = neighbour1 - pos;
bitangent = neighbour2 - pos;
normal = cross(tanget, bitanget);
nice, eh? now go looking for some mesh functions ;)
If you haven't enough yet, you can improve the shading smoothnes by calculating the neighbours and normals in the pixel shader. Thats especially good, if the resolution of the input grid is not so high. For that, after calculating the position, pass the u and v values to the pixel shader and do all remainding stuff there. Read the forums discussion from where the idea is and the download here.
Finally a picture of the supersmooth pixel shader implementation:

