I agree with John Bentin, that the question has little sense in the formal category.
Let me restrict myself to the simpler question of factorization into linear factors
(so everyting will be over complex numbers).
A formal infinite product $(1+a_1x)(1+a_2x)(1+a_3x)...$ makes no sense, because we have to
sum infinite series to obtain the coefficient at $x$.
For this reason I disagree with solution of Will Sawin. The product $(1+a_1x)(1+a_2x^2)...$
does make sense, formally, but after you factor each multile into linear factors, it does not,
even formally!
So let us restrict to convergent series and products. Then the next question, of course, is
"where are they supposed to converge?"
If we want convergence in the whole plane, a satisfactory answer can be given. This is not exactly
Weierstrass, but Weierstrass combined with Hadamard:
If
$$\rho:=\limsup\frac{n\log n}{-\log|a_n|}\; \; \; < 1,$$
then
$$1+a_1x+a_2x^2+\ldots=(1+c_1x)(1+c_2x)\ldots,$$
and both sides are convergent in the ordinary sense in the whole plane. Notice that
some convergence of
$$\sum c_n$$
is necessary just to make sense of the RHS. If this series of $c_n$ is ABSOLUTELY convergent,
then then the right hand side is convergent in the whole plane, so it defines the left hand side
(the power series, and for this power series $\rho\leq 1$ must hold.
So we see that the condition $\rho<1$ is almost best possible. The case $\rho=1$ can be also
completely studied, but this is somewhat more complicated, and for this I refer to B. Levin, Distribution of zeros of entire functions,
AMS, 1964.
Let me repeat the main point. For the product $(1+c_1x)(1+c_2x)\ldots$ to make sense, even FORMALLY,
the series of $c_n$ must converge. Now, of course one can discuss various meanings of this convergence. But if it converges in the most usual sense, that is absolutely, the product converges in the whole plane, and we are within the subject of entire functions.
Exact condition on the coefficients of the power series can be written, and for this I refer
to Levin's book. So a necessary and sufficient condition on the series which guarantees that
it has a product expansion with linear factors can be given.
Passing from complex to real product is routine, of course.
I will try to by more helpful than I was in the comments. First a general observation.
Your infinite product contains socalled
Euler functions,
$\phi(q)=\prod_{n=1}^{\infty}(1-q^n)$
which have the
Lambert series expansion
$\ln\phi(q)=-\sum_{n=1}^{\infty}\frac{1}{n}\frac{q^{n}}{1-q^{n}}$.
Part of your infinite product can be expanded in this way,
$P(x,y)\equiv\prod_{n=1}^{\infty}\frac{1}{(1-x^{n/2}y^{n/2})(1-x^{n/2}y^{-n/2})}=
\frac{1}{\phi(\sqrt{xy})\phi(\sqrt{x/y})}=\exp\left[\sum_{n=1}^{\infty}\frac{1}{n}g(x^n,y^n)\right]$.
The function $g(x,y)$ is given by
$g(x,y)=\frac{(xy)^{1/2}}{1-(xy)^{1/2}}+\frac{(x/y)^{1/2}}{1-(x/y)^{1/2}}=\sum_{n=1}^{\infty}x^{n/2}(y^{n/2}+1/y^{n/2}).$
This is not quite what you have written. I have difficulty verifying your expression. For example, the limit $x\rightarrow 0$ of the left-hand side of your expression is $1+\sqrt{x}(\sqrt{y}+1/\sqrt{y})$, but the same limit of the right-hand side is $1+\sqrt{x}(y+1/y)$.
UPDATE 1: Your corrected expression is still problematic; the left-hand side contains the infinite product
$\prod_{n=1}^{\infty}(1+x^n)=\exp\left[-\sum_{n=1}^{\infty}(-1)^{n}\frac{1}{n}G(x^{n})\right]$, with $G(x)=x/(1-x)$.
The factor $(-1)^n$ in the sum over $n$ seems inconsistent with the right-hand side of your expression --- but in fact it is consistent (see update 2)
.
UPDATE 2: One more identity is needed, in addition to the Euler function identities, to complete the identification:
$\prod_{n=1}^{\infty}(1+x^n)=\exp\left[\sum_{n=1}^{\infty}\frac{1}{n}\frac{x^n}{1-x^{2n}}\right]$.
We then have, quite generally,
$\prod_{n=1}^{\infty}\frac{(1+x^n)}{(1-q_1^n)(1-q_2^n)}=\exp\left(\sum_{n=1}^{\infty}\frac{1}{n}[F(x^n)+g(q_1^n)+g(q_2^n)]\right)$
with the functions $F(x)=x/(1-x^2)$, $g(q)=q/(1-q)$.
Your equation (1) corresponds to $q_1=\sqrt{x/y}$, $q_2=\sqrt{xy}$.
I know, your function $f(x,y)$ looks much more lengthy, but it is really just $F(x)+g(\sqrt{xy})+g(\sqrt{x/y})$.
Best Answer
One nice relation can be written as follow: $$\sum\limits_{n = - \infty }^ \infty z^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^n z^{n^3}= \sum\limits_{n = - \infty }^ \infty (-1)^n z^{2n^3} \phi(z^{6n}) $$
$\quad z \in \mathbb{C},\; |z|=1 $
Where $$\phi(q)=\sum\limits_{m = - \infty }^ \infty (-1)^m q^{m^2}$$
More generally if we write $$\sum\limits_{n = - \infty }^ \infty z^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^nz^n q^{n^2} h^{n^3}= \sum\limits_{n = - \infty }^ \infty \sum\limits_{m = - \infty }^ \infty (-1)^{n+m} z^{2n} q^{2(n^2+m^2)} h^{2n(n^2+3m^2)} \tag{1} $$
The proof of relation (1):
$$F(z)=\sum\limits_{n = - \infty }^ \infty z^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^nz^n q^{n^2} h^{n^3}$$
$$F(-z)=\sum\limits_{n = - \infty }^ \infty (-z)^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^n (-z)^n q^{n^2} h^{n^3}$$
$$F(-z)=\sum\limits_{n = - \infty }^ \infty (-1)^nz^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^n (-1)^n z^n q^{n^2} h^{n^3}$$ $$F(-z)=\sum\limits_{n = - \infty }^ \infty (-1)^nz^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty z^n q^{n^2} h^{n^3}$$
Because of $F(-z)=F(z)$ , we can write that:
$$F(z)=\sum\limits_{n = - \infty }^ \infty (-1)^nz^{2n} q^{2n^2} h^{2n^3} A_n(q,h)$$
We need to find $A_n(q,h)$
$$\sum\limits_{n = - \infty }^ \infty (-1)^nz^{2n} q^{2n^2} h^{2n^3} A_n(q,h)=\sum\limits_{n = - \infty }^ \infty z^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^nz^n q^{n^2} h^{n^3}$$
Let's use the transformation
$z=ZQ^{2}h^{3}$
$q=Qh^{3}$
$$\sum\limits_{n = - \infty }^ \infty (-1)^nZ^{2n} Q^{2n^2+4n} h^{2n^3+6n^2+6n} A_n(Qh^3,h)=\sum\limits_{n = - \infty }^ \infty Z^n Q^{n^2+2n} h^{n^3+3n^2+3n} .\sum\limits_{n = - \infty }^ \infty (-1)^nZ^n Q^{n^2+2n} h^{n^3+3n^2+3n}$$
If we multiply both side by $Z^2Q^2h^2$
$$\sum\limits_{n = - \infty }^ \infty (-1)^nZ^{2n+2} Q^{2n^2+4n+2} h^{2n^3+6n^2+6n+2} A_n(Qh^3,h)=\sum\limits_{n = - \infty }^ \infty Z^{n+1} Q^{n^2+2n+1} h^{n^3+3n^2+3n+1} .\sum\limits_{n = - \infty }^ \infty (-1)^nZ^{n+1} Q^{n^2+2n+1} h^{n^3+3n^2+3n+1}$$
$$\sum\limits_{n = - \infty }^ \infty (-1)^nZ^{2(n+1)} Q^{2(n+1)^2} h^{2(n+1)^3} A_n(Qh^3,h)=\sum\limits_{n = - \infty }^ \infty Z^{n+1} Q^{(n+1)^2} h^{(n+1)^3} .\sum\limits_{n = - \infty }^ \infty (-1)^nZ^{n+1} Q^{(n+1)^2} h^{(n+1)^3}$$
$$\sum\limits_{n = - \infty }^ \infty (-1)^{n-1}Z^{2n} Q^{2n^2} h^{2n^3} A_{n-1}(Qh^3,h)=\sum\limits_{n = - \infty }^ \infty Z^{n} Q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^{n-1}Z^{n} Q^{n^2} h^{n^3}$$
$$\sum\limits_{n = - \infty }^ \infty (-1)^{n}z^{2n} q^{2n^2} h^{2n^3} A_{n-1}(qh^3,h)=\sum\limits_{n = - \infty }^ \infty z^{n} q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^{n}z^{n} q^{n^2} h^{n^3}$$
$$A_n(q,h)=A_{n-1}(qh^3,h)$$ $$A_{-1}(q,h)=A_{0}(qh^{-3},h)$$ $$A_1(q,h)=A_{0}(qh^3,h)$$
$$A_2(q,h)=A_{1}(qh^3,h)=A_{0}(qh^6,h)$$ $$A_{-2}(q,h)=A_{-1}(qh^{-3},h)=A_{0}(qh^{-6},h)$$ $$A_3(q,h)=A_{2}(qh^3,h)=A_{1}(qh^6,h)=A_{0}(qh^9,h)$$ $$A_n(q,h)=A_{0}(qh^{3n},h)$$
We need to find $A_{0}(q,h)$ to complete the proof:
We need to focus on $z^0$ terms
$$F(z)=\sum\limits_{n = - \infty }^ \infty z^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^nz^n q^{n^2} h^{n^3}$$
If we multiply terms by terms for $z^0$ terms
$$A_{0}(q,h)=\sum\limits_{m = - \infty }^ \infty (-1)^m q^{2m^2}$$ $$A_{n}(q,h)=A_{0}(qh^{3n},h)=\sum\limits_{m = - \infty }^ \infty (-1)^m q^{2m^2}h^{6nm^2}$$
Thus We can write that
$$\sum\limits_{n = - \infty }^ \infty (-1)^nz^{2n} q^{2n^2} h^{2n^3} \sum\limits_{m = - \infty }^ \infty (-1)^m q^{2m^2}h^{6nm^2}=\sum\limits_{n = - \infty }^ \infty z^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^nz^n q^{n^2} h^{n^3}$$
$$\sum\limits_{n = - \infty }^ \infty \sum\limits_{m = - \infty }^ \infty (-1)^{n+m}z^{2n} q^{2n^2+2m^2} h^{2n^3+6nm^2}=\sum\limits_{n = - \infty }^ \infty z^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^nz^n q^{n^2} h^{n^3}$$
$$\sum\limits_{n = - \infty }^ \infty (-1)^nz^{2n} q^{2n^2} h^{2n^3} \phi(q^2h^{6n}) =\sum\limits_{n = - \infty }^ \infty z^n q^{n^2} h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^nz^n q^{n^2} h^{n^3}$$
If we put that $z=1$ and $q=1$
$$\sum\limits_{n = - \infty }^ \infty (-1)^n h^{2n^3} \phi(h^{6n}) =\sum\limits_{n = - \infty }^ \infty h^{n^3} .\sum\limits_{n = - \infty }^ \infty (-1)^n h^{n^3} $$
Finally, the same method can be used for $\sum\limits_{n = - \infty }^ \infty z^{n^k}$ to get similiar identies where $k>3$