A question about a formal power series manipulation  I want to find a function $f(x,y)$ which can satisfy the following equation,
$\prod _{n=1} ^{\infty} \frac{1+x^n}{(1-x^{n/2}y^{n/2})(1-x^{n/2}y^{-n/2})} = exp [ \sum _{n=1} ^\infty \frac{f(x^n,y^n)}{n(1-x^{2n})}]$


*

*I would like to know how this is solved.


(..though I landed into this through a different route (calculating Witten Index!), such expressions also occur in finite dimensional representation theory where a generating function for the character of the (anti)symmetric powers of a representation (the LHS) is written as a Plethystic exponential (the RHS) of the original (generally fundamental) representation...)
One can perturbatively check that the following function satisfies the above equation, 


*

*$f(x,y) = \sqrt{x} (\sqrt{y} + 1/\sqrt{y}) + x (1 + y + 1/y) + x^{3/2} (y^{3/2} + 1/y^{3/2}) + x^2 (y^2 + 1/y^2)$
$\quad\quad\quad\quad\quad + \frac{(x y)^{5/2} }{(1 - \sqrt{x y})}(1 - 1/y^2) + \frac{(x/y)^{5/2}}{(1 - \sqrt{x/y})} (1 - y^2) 
$
The paper doesn't state any proof or explanation for how this was obtained but the above is order-by-order in $x$ checkable to be right after truncating the original equation at any finite value of $n$. (...I don't know how to check this keeping the full sum/product over $n$..)

Now I tried to do something obvious but it didn't work!
$\prod_{n=1}^{\infty} \frac{ (1+x^n) }{1+x^n -x^{\frac{n}{2}} \left(y^{\frac{n}{2}}  + y^{-\frac{n}{2}}\right) } = \exp \left[ \sum_{n=1}^{\infty}   \frac{ I_{ST}(x^n,y^n) } {n (1-x^{2n}) } \right] $
$\Rightarrow \sum_{n=1}^{\infty} \left[ \ln (1+x^n) - \ln(1-(\sqrt{xy})^n) 
- \ln\left(1- \left(\sqrt{\frac{x}{y}}\right)^n\right) \right]  = \sum_{n=1}^\infty  \frac{I_{ST}(x^n,y^n)} {n(1-x^{2n})} $
Now we expand the logarithms and we have, 
$\sum_{n=1}^{\infty} \left[  \sum_{a=1}^{\infty} (-1)^{a+1} \frac{x^{na}}{a} + 
\sum_{b=1}^{\infty} \frac{ (\sqrt{xy})^{nb} } {b} + \sum_{c=1}^{\infty} 
\frac{ (\sqrt{x/y})^{nc} }{c}   \right] =  \sum _{n=1}^\infty  
\frac{f(x^n,y^n)} {n(1-x^{2n})}$
$\Rightarrow \sum _{a=1} ^{\infty} \frac{1}{a} \left[  \sum _{n=1} ^{\infty} \left( (-1)^{a+1}x^{na} + (xy)^{\frac{na}{2}} + \left(\frac{x}{y}\right)^{\frac{na}{2}} \right)  \right] =  \sum _{n=1} ^\infty  \frac{f(x^n,y^n)} {n(1-x^{2n})}$ 
By exchanging $a$ and $n$ (relabeling on the LHS) and matching the patterns on both sides and picking out the $n=1$ term one sees that one way this equality can hold is if,
2.
$f(x,y) = (1-x^2) \sum _{a=1} ^{\infty} \left[ x^a + (xy)^{\frac{a}{2}} + (\frac{x}{y})^{\frac{a}{2}}  \right]$
$\Rightarrow f (x,y) = (1-x^2) \left(-1 + \frac{1}{1-x} -1 + \frac{1}{1-\sqrt{xy}} - 1 + \frac{1}{1-\sqrt{\frac{x}{y}} } \right)$
But this solution is neither the one above which could be perturbatively checked to be true nor does it satisfy the original equation! Why?  
After doing a series expansion of the above (using Series on Mathematica) one sees that this above derived equation (2) differs from (1) in having just one extra term of $x^2$. (...I would like to know what is wrong in the derivation that gives (2) this one extra term compared to the non-derivable but perturbatively checked correct answer (1)...) 
 A: 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})$.
