Laplace equation in transformed coordinate system I am writing a solver for the 2D Laplace equation in the transformed coordinate system. I have derived the final form of the Laplace equation, however I am not sure if it is right because the FEM solution (in $x$, $y$) and FDM solution (in $u$, $v$) are slightly different.
The coordinate transformation is given with the following expressions:
$$u=x$$
$$v=\frac{y-b(x)}{t(x)-b(x)}$$
where $u$,$v$ are the new coordinates, $b(x)$ is the function describing bottom of the domain, $t(x)$ is the function describing top of the domain. The situation is displayed in the figure.

The final form of the Laplace equation I've got is:
$$\Delta f_{u,v}=\frac{\partial^2f}{\partial u^2}+2\frac{\partial v}{\partial x}\frac{\partial^2f}{\partial u\partial v}+\left[\left(\frac{\partial v}{\partial x}\right)^2+\left(\frac{\partial v}{\partial y}\right)^2\right]\frac{\partial^2f}{\partial v^2}=0$$
$\frac{\partial u}{\partial x}=1$; $\frac{\partial u}{\partial y}=0$;
$\frac{\partial v}{\partial x}=\frac{-b'(x)\left(t(x)-b(x)\right)-(y-b(x))(t'(x)-b'(x))}{\left(t(x)-b(x)\right)^2}$; $\frac{\partial v}{\partial y}=\frac{1}{t(x)-b(x)}$
Additioanlly I am not sure what to do with the $x$ and $y$ coordinates that remain in $\frac{\partial v}{\partial x}$ and $\frac{\partial v}{\partial y}$. (I've just switched them to $u$ and $v$ based on the inverse transformations $x=u$ and $y=v(t(x)-b(x))+b(x)$):
$\frac{\partial v}{\partial x}_{x \rightarrow u, y\rightarrow v}=\frac{-b'(u)\left(t(u)-b(u)\right)-([v(t(u)-b(u))+b(u)]-b(u))(t'(u)-b'(u))}{\left(t(u)-b(u)\right)^2}$; $\frac{\partial v}{\partial y}_{x \rightarrow u, y\rightarrow v}=\frac{1}{t(u)-b(u)}$
I that ok?
EDIT:
I have used the laplacian formulation in the form $\Delta f=div(grad f)$ in the tensor notation which according to this pdf gives:
$$\Delta f=\frac{1}{J}\frac{\partial}{\partial \zeta_j}\left[J\left(\frac{\partial f}{\partial \zeta_k}\frac{\partial \zeta_k}{\partial x_i}\frac{\partial \zeta_j}{\partial x_i}\right)\right]$$ where $x_1=x$, $x_2=y$ and $\zeta_1=u$, $\zeta_2=v$.
After expansion of the indexes I've got:
$$\Delta f=\frac{1}{J}\left\{\frac{\partial}{\partial u}\left[J\left(\frac{\partial f}{\partial u} \frac{\partial u}{\partial x}\frac{\partial u}{\partial x}+\frac{\partial f}{\partial u}\frac{\partial u}{\partial y}\frac{\partial u}{\partial y}+\frac{\partial f}{\partial v}\frac{\partial v}{\partial x}\frac{\partial u}{\partial x}+\frac{\partial f}{\partial v}\frac{\partial v}{\partial y}\frac{\partial u}{\partial y}\right) \right]+ \frac{\partial }{\partial v}\left[J\left( \frac{\partial f}{\partial u} \frac{\partial u}{\partial x}\frac{\partial v}{\partial x}+\frac{\partial f}{\partial u}\frac{\partial u}{\partial y}\frac{\partial v}{\partial y}+\frac{\partial f}{\partial v}\frac{\partial v}{\partial x}\frac{\partial v}{\partial x}+\frac{\partial f}{\partial v}\frac{\partial v}{\partial y}\frac{\partial v}{\partial y}\right) \right] \right\}$$
And since $\frac{\partial u}{\partial x}=1$ and $\frac{\partial u}{\partial y}=0$ the laplacian takes the following form
$$\Delta f=\frac{1}{J}\left\{\frac{\partial}{\partial u}\left[J\left(\frac{\partial f}{\partial u} +\frac{\partial f}{\partial v}\frac{\partial v}{\partial x}\right) \right]+ \frac{\partial }{\partial v}\left[J\left( \frac{\partial f}{\partial u}\frac{\partial v}{\partial x}+\frac{\partial f}{\partial v}\left(\frac{\partial v}{\partial x}\right)^2+\frac{\partial f}{\partial v}\left( \frac{\partial v}{\partial y}\right)^2\right) \right] \right\}$$
Here some more of my doubts arise: if the Jacobian is  $J=t(x)-b(x)$ should it be treated as constant when $\frac{\partial}{\partial u}$ or $\frac{\partial}{\partial v}$ is applied? I have the same concern about applying $\frac{\partial}{\partial u}$ or $\frac{\partial}{\partial v}$ to $\frac{\partial v}{\partial x}$ or $\frac{\partial v}{\partial y}$. I have assumed that since there is no $u$ or $v$ in those exprerssions they are not differentiated. Therefore no first order derivatives of $f$ appeared.
According to MathWonk's answer the result is: $$\Delta f_{u,v}=\frac{\partial^2f}{\partial u^2}\left(\left(\frac{\partial u}{\partial x}\right)^2+ \left(\frac{\partial u}{\partial y}\right)^2\right )+2\frac{\partial^2f}{\partial u\partial v}\left(\frac{\partial u}{\partial x}\frac{\partial v}{\partial x}
+\frac{\partial u}{\partial y}\frac{\partial v}{\partial y}\right )+
\frac{\partial^2f}{\partial v^2}\left(\left(\frac{\partial v}{\partial x}\right)^2+ \left(\frac{\partial v}{\partial y}\right)^2\right )+\frac{\partial f}{\partial u}\left(\frac{\partial^2u}{\partial x^2}+ \frac{\partial^2u}{\partial y^2} \right )+\frac{\partial f}{\partial v}\left(\frac{\partial^2v}{\partial x^2}+ \frac{\partial^2v}{\partial y^2} \right )$$
which after substituting derivatives of $u$ and $v$ yields:
$$\Delta f_{u,v}=\frac{\partial^2f}{\partial u^2}
+2\frac{\partial^2f}{\partial u\partial v}\frac{\partial v}{\partial x}+
\frac{\partial^2f}{\partial v^2}\left(\left(\frac{\partial v}{\partial x}\right)^2+ \left(\frac{\partial v}{\partial y}\right)^2\right )+
\frac{\partial f}{\partial v}\left(\frac{\partial^2v}{\partial x^2} \right )$$
which is a completely different expression than mine. Where did I make a mistake?
 A: I am editing my original answer so that it addresses the issues you raised in your edited version of your question.
I think you are confused by two rival methods for computing the Laplacian in curvilinear coordinates. One method (method B) works best when x and y are known functions of the curvilinear coordinates u and v, and the other method (C) works  best when u and v are known functions of x and y. In your problem you can easily invert either set of transformation equations, so you are free to choose whichever method you prefer.
A. Note that your original proposed answer is clearly incorrect because it omits some first-order derivatives.
B. In your edited posting, you introduce a completely different  method for computing the Laplacian,  involving the expression J. That method works best when you are given x and y as explicit functions of u and v.  (The opposite of what is described in your original posting.)
Thus, to use that J-method you should first invert your set of coordinate transformations. If you do that, the formulas you quoted involving J will work fine.  P.S. In those formulas you should treat J   as a function of u and v. The calculations will be tedious but eventually you will get the same answer as the next method.
C. In your original posting, u and v are known functions of x and y.
In this situation my original answer works well.  Here is a recap of that original answer.
The strategy is to use vector identities to break up the calculation of the laplacian into smaller steps.
(1) Definition of Laplacian:  If f= f(u,v) then $\triangle f= div (grad f)$
For the innermost gradients you can use
(2) Chain rule for gradients:
grad f= $f_u grad \ u + f_v grad \ v = f_u \nabla u + f_v \nabla v$
and for the outer-most divergence you can use
(3) Product Rule for divergence:
$ div ( h \ grad u) = grad h  \cdot grad u + h \ div \ grad u =\nabla h \cdot \nabla g  + h \triangle u$
and  a similar expansion  for  $div ( h \ grad v)$.
The final result is $$ \triangle f=f_{uu} \nabla u \cdot \nabla u + 2 f_{u,v} \nabla u \cdot \nabla v+ f_{vv} \nabla v \cdot \nabla v + f_u \triangle u + f_v \triangle v$$
Here the functions $u$ and $v$ are originally expressed in terms of x and y. The calculation of all the expressions $\nabla u, \nabla v, \triangle u, \triangle v$   can be done using $x$ and $y$ derivatives.
Finally these expressions should  ultimately be re-expressed in terms of $u$ and $v$ at the end.
Finally, note that you can check your answer by testing a few simple cases by choosing various expressions for $f$ and computing $\triangle f$. For example, does your Laplace formula correctly compute that  $\triangle x =0$ and $\triangle y=0$?
