@Rouben Rostamian Thank you for your replies.
The function u(x,y) is the solution to a PDE (elliptical) which represents a 2D surface over its domain. The point (xo,yo) lies on the bounday of the domain. I wish to calculate the Neumann boundary condition (the normal derivative) at (xo,yo), i.e. du/dn = n . del(u) . I may have confused the issue by showing u in 'bold'. Actually, u is a scalar function not a vector function.
The boundary does not need to be a level curve (although it can be); for example, the solution to a well posed Laplace equation with a Dirichlet boundary whose value varies over its length. This could represent (at equilibrium) the Heat equation with different temperatures imposed along the boundary.