I'm writing this mostly for my own benefit, since I'm currently working on the same things.
First, we show that if $h(z)$ is harmonic in a simply connected domain $D$ (a domain being a connected open set, simply connected meaning 'no holes'), then there is $F(z)$ analytic (synonymous with holomorphic; use 'real analytic' for a real-valued funcion) with $\Re F(z) = h(z)$.
As the previous poster pointed out, the first observation is that for harmonic $h(z)$, the function $\psi = \frac{\partial h}{\partial x} - i\frac{\partial h}{\partial y}$ is also harmonic. Again, as mentioned, this follows from the C.R. equations, but one must also note that we have continuity of the partials (there are many examples where the C.R. equations hold, yet the function is not analytic).
Now, since $D$ is simply connected and the function $\psi$ is analytic, it has a primitive in $D$ (simple connectivity is crucial here, as this is where Cauchy's theorem holds). Call this primitive $G(z)$, so
$$G'(z) = \frac{\partial h}{\partial x} - i\frac{\partial h}{\partial y}.\tag{1}$$
But we know that $G$ is analytic, so $G = u+iv$. Thus,
\begin{gather*}
G' = u_x + iv_x = v_y-iu_y
\end{gather*}
Comparing this last equation to (1) (using $u_x+iv_x$ for the first comparison and $v_y-iu_y$ for the second), we find
\begin{gather*}
\Re \frac{\partial G}{\partial x} = \frac{\partial h}{\partial x} = \frac{\partial \Re G}{\partial x} \cr
\Re \frac{\partial G}{\partial y} = \frac{\partial h}{\partial y} = \frac{\partial \Re G}{\partial y}
\end{gather*}
Thus,
\begin{gather*}
\frac{\partial}{\partial x}(\Re G(z) - h(z)) = 0 \cr
\frac{\partial}{\partial y}(\Re G(z) - h(z)) = 0
\end{gather*}
Since $D$ is a connected domain, we must have $h(z) = \Re G(z) + c$, so setting $F(z) = G(z) + c$ gives us $h(z) = \Re F(z)$, as desired.
We now have, given a harmonic function $h$, an analytic function whose real part is $h$. This required some work and wasn't totally trivial; in particular, we needed simple connectivity of the domain. What's important for us is the total converse: given $h(z) = \Re F(z)$ for $F$ analytic, $h$ is harmonic. You can prove this using the Cauchy-Riemann equations
Now, we want to show that the composition of a harmonic function with an analytic one is harmonic. Suppose then $\phi(z)$ is harmonic in $V$, and let $f(z)$ be analytic in $U$ with $f(U) \subset V$. Then $u(f(z))$ is harmonic for $z \in U$.
To prove this, keep in mind that harmonicity is a local property. By the previous result, we have (locally) some function $\varphi(z)$, analytic, such that $\phi(z) = \Re \varphi(z)$. Let then $z_0 \in U$ with $f(z_0) = w_0 \in V$. Then $\varphi(f(z)) = \varphi(w)$ will be analytic for $w$ near $w_0$, so
$$u(f(z)) = \Re \varphi(f(z)). $$
By the total converse, $u(f(z))$ is harmonic for $z$ close enough to $z_0$.
Best Answer
Possible idea. Let $f=f(x,y)$. If $f$ and $1/f$ are harmonic then $$ \frac{\partial^2 f}{\partial x^2} + \frac{\partial^2 f}{\partial y^2} = 0 $$ $$ \frac{\partial^2 1/f}{\partial x^2} + \frac{\partial^2 1/f}{\partial y^2} = 0 $$ where $$ \frac{\partial 1/f}{\partial x} = -\frac{\partial f}{\partial x} \frac{1}{f^2} \qquad \frac{\partial^2 1/f}{\partial x^2} = -\frac{\partial^2 f}{\partial x^2} \frac{1}{f^2} +2\left(\frac{\partial f}{\partial x}\right)^2\frac{1}{f^3} $$
$$ \frac{\partial 1/f}{\partial y} = -\frac{\partial f}{\partial y} \frac{1}{f^2} \qquad \frac{\partial^2 1/f}{\partial y^2} = -\frac{\partial^2 f}{\partial y^2} \frac{1}{f^2} +2\left(\frac{\partial f}{\partial y}\right)^2\frac{1}{f^3} $$ So $$ -\frac{\partial^2 f}{\partial x^2} \frac{1}{f^2} +2\left(\frac{\partial f}{\partial x}\right)^2\frac{1}{f^3}-\frac{\partial^2 f}{\partial y^2} \frac{1}{f^2} +2\left(\frac{\partial f}{\partial y}\right)^2\frac{1}{f^3} = \frac{2}{f^3}\left(\left(\frac{\partial f}{\partial x}\right)^2+\left(\frac{\partial f} {\partial y}\right)^2\right) $$
$$ \frac{2}{f^3}\left(\left(\frac{\partial f}{\partial x}\right)^2+\left(\frac{\partial f} {\partial y}\right)^2\right)=0 $$ Since $f \not \equiv 0$, we have $$ 0=\left(\frac{\partial f}{\partial x}\right)^2+\left(\frac{\partial f} {\partial y}\right)^2 =\left(\frac{\partial f}{\partial x}+\frac{\partial f} {\partial y}i\right)\left(\frac{\partial f}{\partial x}-\frac{\partial f} {\partial y}i\right) $$