To understand this question, please first understand the question and answer here: Create a formula that creates a curve between two points
We are essentially transcending a 2d problem into a 3d problem.
I have 4 points arranged as the corners of a 1 unit square on the x-y-plane. Each point also has an arbitrary Z value, a x-gradient, a y-gradient.
To be clear, by "x-gradient" I mean that you draw a tangent line on the X-Z-Plane from the point and get (the change in z) divided by (the change in X).
First describe a surface that exists at those four points and has the correct gradients at those points.
Secondly, create a function that returns the Z value for any two values of X and Y between 0 and 1.
Am I right if I say that you have 3 constraints at each "corner": its $z$ (height) its $x$ gradient and its $y$ gradient ? i.e., 12 informations.
Why don't you try surfaces depending on (at least) 12 parameters, for example
$$z=ax^4+by^4+cx^3y+dxy^3+ex^2y^2+fx^3+gy^3+hx^2y+ixy^2+jx+ky+l$$
Imposing the different constraints will give you a system of 12 linear equations in 12 unknowns.
It can be also a surface of the form $z=f(x)g(y)$, $f$ and $g$ depending resp. on 6 and 7 parameters resp. (like above) (the 7th parameter is necessary because you can simultaneously $f$ by a constant $k$ and multiply $g$ by $k$). The equations are no longer linear, but the constraints on derivatives are more easily expressed.
You can also turn to what is known under the name of "Bezier patches" (or "Hermite patches").