Is there a way to decompose a polynomial $xy+f(x,y)$ into the product of two convergent power series $x+g(x,y)$ and $y+h(x,y)$? Suppose we have a polynomial $xy+f(x,y)$, where $f(x,y)$ is a polynomial in $\mathbb C[x,y]$ whose the lowest degree term has degree at least 3. My question is, are we always able to decompose $xy+f(x,y)$ into the product of two convergent power series $x+g(x,y)$ and $y+h(x,y)$ in a neighborhood of $(0,0)$, where terms in $g$ and $h$ have order higher than 1?
I have no idea about the convergence of power series with two or more variables. Any solution or reference will be appreciated! Please also note that it's not enough to simply decompose $xy+f(x,y)$ formally into two $x+g(x,y)$ and $y+h(x,y)$. We also need the convergence of $g$ and $h$! 
 A: Note: In Rick Miranda's book "Algebraic Curves and Riemann Surfaces", Lemma 2.3 on page 68 gives an argument that this factorization is possible in terms of formal power series, and then says that "it is an easy exercise" to show that these actually converge. His factorization is not in any way uniquely specified, and the formal factorization is not unique and can certainly be made to diverge, so it is not quite clear what he had in mind. Here is an alternative proof using the implicit function theorem:
We first claim that there exists a local analytic function $y=x^2 z(x)$ whose graph is contained in the zero set of $F(x,y) = xy + f(x,y)=0$. In order to do so, we change variables to $(x,z)$ by $y=x^2 z$ and divide the equation by $x^3$ to get
$$
P(x,z) = \frac{F(x,x^2 z)}{x^3} = z + \frac{f(x,x^2 z)}{x^3} = 0.
$$
Since every term in $f$ has degree $\ge 3$, we know that $P(x,z)$ is a polynomial in $x$ and $z$, and that $P(0,z) = z-a$ where $-a$ is the coefficient of $x^3$ in $f(x,y)$. By the implicit function theorem there exists a local analytic function $z(x)$ such that $P(x,z(x))=0$, with $z(0) = a$. Then $h(x) = x^2 z(x)$ is analytic with $F(x,h(x)) = 0$. This implies that for every fixed $x$, the polynomial $y \mapsto F(x,y)$ has a zero at $y=h(x)$, so it has a linear factor of $y-h(x)$, implying that the function
$$
G(x,y) = \frac{F(x,y)}{y-h(x)} = \frac{xy + f(x,y)}{y-x^2z(x)} = x+g(x,y)
$$
is actually analytic near $(0,0)$, with $g(x,y)$ containing only terms of degree $\ge 2$. This shows
$$
xy + f(x,y) = (y-h(x))(x+g(x,y))
$$
with both $g$ and $h$ containing only higher-order terms. Note that this type of factorization is unique (by the uniqueness part of the implicit function theorem.)
Note also that this factorization can be made more symmetric by observing that $F(x,y)$ also has a factor of the form $x-g(y)$, by the exact same argument as above with variables switched, so we can factor
$$
F(x,y) = xy + f(x,y) = (x-g(y))(y-h(x))u(x,y)
$$
with $g(0)=g'(0)=h(0)=h'(0)=0$ and $u(0,0)=1$ analytic near $(0,0)$. This shows that the zero set of $F$ is locally the union of the graphs given by $x=g(y)$ and $y=g(x)$ which intersect transversally.
ADDENDUM: The fact that $G(x,y)$ is analytic in both $x$ and $y$ follows from the Cauchy Integral Formula (and this must be a standard result, but I am not an expert in several complex variables.) Since $h(x) = O(x^2)$, we can find $r>0$ such that $|h(x)|<r$ for $|x|=r$. Then Cauchy's Integral Formula applied to the polynomial function $y \mapsto G(x,y)$ gives
$$
G(x,y) = \frac{1}{2\pi i} \int_{|\zeta|=r} \frac{G(x,\zeta)}{\zeta - y} \, d\zeta
$$
for all $|x|,|y|<r$. The function $x \mapsto G(x,\zeta) = \frac{F(x,\zeta)}{\zeta - h(x)}$ is analytic for all $|\zeta|=r$ since the denominator never vanishes on this circle, so $x \mapsto G(x,y)$ is analytic in $|x|<r$.
