Solve a first order non-linear PDE, or just find a residue of the solution In the process of solving a combinatorial problem in my research, I came up with recursive formula, and assembled a two variable generating function so that the rescursion is encoded by the PDE:
$$
x\frac{\partial}{\partial x}(f(x,y)^2) = \frac{\partial}{\partial y} f(x,y)
$$
with initial condition $f(x,0) = x + x^{-1}$.
Actually, I don't need the whole solution $f$, I just need the residue around $x=0$, that is, the function of $y$ that is the coefficient of $x^{-1}$ in the Laurent expansion.
I don't know much about PDEs, so I'm not even sure where to start, or if this is a good way to solve it.
I am hoping to prove that the residue is something like
$$
\sum_n \frac{4^n}{(n+1)!n!} y^{2n}
$$
 A: This is not a complete answer, but may lead you or someone else to find a complete solution. I just thought it's too long for a comment. 
Let's see if we can make some progress on finding a solution to
$$ F = 2 xp f - q = 0$$
where $p = f_x$ and $q = f_y$. This quasilinear PDE must be solved together with $f(x,0) = x + 1/x$. If we parametrize this as $x = x_o(s) = s$, $y = y_o(s) = 0$ and $f = f_o(s) = s+1/s$, then the solution to the Lagrange-Charpit relations,
\begin{align*}
\mathrm{d}x/\mathrm{d}t & =  2 x f \\
\mathrm{d}y/\mathrm{d}t & = -1 \\
\mathrm{d}f/\mathrm{d}t & = 0,
\end{align*} is given by
$$\log(x/s) = 2  \, ( s + 1/s) \, t, \quad y = -t, \quad f = s + 1/s $$
Now, to express $f = f(x,y)$, one may eliminate $s$ and $t$ from the expressions for $x(s,t)$ and $y(s,t)$ but, as you can see, this is quite difficult since it involves the solution of some ugly-looking trascendental equation.

One could also work with 
$$ \frac{\mathrm{d}x}{2xf} = \frac{\mathrm{d}y}{-1} = \frac{\mathrm{d}f}{0} $$ 
and try to solve it via characteristics. But, I'm dim-witted when it comes to applying the initial condition. Indeed, from the last relations we get $f = c_1$ and $\log{x} + 2 c_1 y = c_2$. Put now $c_2$ as a function of $c_1$ and we obtain the solution in the implicit form:
$$ \log{x} + 2  y f = G(f)  $$
where $G$ is an arbitrary function of its argument. By imposing the initial condition one gets:
$$ \log{x} + 0  = G(x+1/x)$$
define the dummy variable $\sigma = x+1/x \implies x = \frac{1}{2} \left( \sigma \pm \sqrt{\sigma^2 - 4} \right)$ for which one concludes that the solution is implicitly given by

$$ \log{x} + 2y f = \log\left[\frac{1}{2} \left( f \pm \sqrt{f^2 - 4} \right)\right] $$

or, alternatively

$$  2 x \exp(2y f) = f \pm \sqrt{f^2-4}  $$

This gives you the value(s) of $f$ for each $x$ and $y$, but I don't know if you can extract anything practical from it. Edit: I just checked in Mathematica that the solution above is indeed correct. 

Hope someone can take it from here or is of any help for you.
