Autocorrelation function $R_{yy}(t_1,t_2)$?
Let $\sigma^2$ denote the common variance of the random variables comprising the zero-mean stationary Gaussian process $\{X(t)\}$. We want to find
$$R_{Y}(t,s) = E[Y(t)Y(s)] = E\big[(X(t))^2(X(s))^2\big].$$
Now,
given that $X(t)=x$, 
\begin{align}
E\big[(X(t))^2(X(s))^2\mid X(t)=x\big] &amp;=
E\big[x^2(X(s))^2\mid X(t)=x\big]\\
&amp;= x^2E\big[(X(s))^2\mid X(t)=x\big]
\end{align}
where we know (since $X(t)$ and $X(s)$ are jointly Gaussian random variables with covariance $R_X(t-s)$ and correlation coefficient $\rho = \sigma^{-2}R_X(t-s)$) that the conditional
distribution of $X(s)$ given that $X(t) = x$ is a Gaussian distribution
with mean $\rho x$ and variance $\sigma^2(1-\rho^2)$. It follows that
\begin{align}
E\big[(X(t))^2(X(s))^2\mid X(t)=x\big] &amp;= x^2E\big[(X(s))^2\mid X(t)=x\big]\\
&amp;= x^2 \left(\sigma^2(1-\rho^2) + \rho^2x^2\right).\end{align}
Hence, the random variable $E\big[(X(t))^2(X(s))^2\mid X(t)]$
equals $(X(t))^2 \left( \sigma^2(1-\rho^2) + \rho^2(X(t))^2\right)$ and the law of iterated expectation gives

\begin{align}
E\big[(X(t))^2(X(s))^2\big]&amp;=E\bigr[E\big[(X(t))^2(X(s))^2\mid X(t)\big]\bigr]\\
&amp;= E\big[(X(t))^2\left(\sigma^2(1-\rho^2) + \rho^2(X(t))^2\right)\big]\\
&amp;= \sigma^2(1-\rho^2)E[(X(t))^2] + \rho^2 E[(X(t))^4]\\
&amp;= \sigma^4\left((1 - \rho^2) + 3\rho^2\right)\\
&amp;= \sigma^4\left(1 + 2\rho^2\right)\\
&amp;= (\sigma^2)^2 + 2 (\rho\sigma^2)^2\\
&amp;= \left(R_X(0)\right)^2 + 2 \left(R_X(t-s)\right)^2\\
&amp;= E[(X(t))^2]E[(X(s))^2] + 2\left(E[X(t)X(s)]\right)^2.
\end{align}
This differs in subtle ways (not just the use of $t$ and $s$ in place of $t_1$ and $t_2$) from what the OP insists is the right answer as stated in several books that he claims to have read, which answer he has faithfully transcribed into his question above.