Solution (source code)

= Solution

Let $D=(1-x)^2+y^2$. Direct <partial differentiation>, away from $(1,0)$, gives
$$
u_x=v_y=\frac{2((1-x)^2-y^2)}{D^2},\qquad
u_y=-v_x=-\frac{4(1-x)y}{D^2}.
$$
These are the <Cauchy-Riemann equations>. Equivalently, with $z=x+iy$, the analytic mapping is $u+iv=(1+z)/(1-z)$, the inverse <Cayley transform between the half-plane and disk>.

The supplied harmonicity result can also be checked by the <multivariable chain rule>. The <Cauchy-Riemann equations> imply $\Delta u=\Delta v=0$, $\nabla u\cdot\nabla v=0$, and $|\nabla u|^2=|\nabla v|^2$. Hence
$$
\Delta\Psi=(u_x^2+u_y^2)(\Phi_{uu}+\Phi_{vv})=0.
$$
This is <conformal invariance of harmonicity> for the given mapping.