Is MLE of $\theta$ asymptotically normal when $(X,Y)\sim e^{-(x/\theta+\theta y)}\mathbf1_{x,y>0}$? Suppose $(X,Y)$ has the pdf
$$f_{\theta}(x,y)=e^{-(x/\theta+\theta y)}\mathbf1_{x>0,y>0}\quad,\,\theta>0$$
Density of the sample $(\mathbf X,\mathbf Y)=(X_i,Y_i)_{1\le i\le n}$ drawn from this population is therefore
\begin{align}
g_{\theta}(\mathbf x,\mathbf y)&=\prod_{i=1}^n f_{\theta}(x_i,y_i)
\\&=\exp\left[{-\sum_{i=1}^n\left(\frac{x_i}{\theta}+\theta y_i\right)}\right]\mathbf1_{x_1,\ldots,x_n,y_1,\ldots,y_n>0}
\\&=\exp\left[-\frac{n\bar x}{\theta}-\theta n\bar y\right]\mathbf1_{x_{(1)},y_{(1)}>0}\quad,\,\theta>0
\end{align}
The maximum likelihood estimator of $\theta$ can be derived as 
$$\hat\theta(\mathbf X,\mathbf Y)=\sqrt\frac{\overline X}{\overline Y}$$

I wish to know whether the limiting distribution of this MLE is normal or not.

It is clear that a sufficient statistic for $\theta$ based on the sample is $(\overline X,\overline Y)$.
Now I would have said that the MLE is asymptotically normal without a doubt if it was a member of the regular one-parameter exponential family. I don't think that is the case, partly because we have a two-dimensional sufficient statistic for a one-dimensional parameter (as in $N(\theta,\theta^2)$ distribution, for example).
Using the fact that $X$ and $Y$ are in fact independent Exponential variables, I can show that the exact distribution of $\hat\theta$ is such that
$$\frac{\hat\theta}{\theta}\stackrel{d}{=} \sqrt F\quad,\text{ where }F\sim F_{2n,2n}$$
I cannot possibly proceed to find the limiting distribution from here.
Instead I can argue by WLLN that $\overline X\stackrel{P}\longrightarrow\theta$ and $\overline Y\stackrel{P}\longrightarrow 1/\theta$, so that $\hat\theta\stackrel{P}\longrightarrow\theta$.
This tells me that $\hat\theta$ converges in distribution to $\theta$. But this does not come as a surprise, since $\hat\theta$ is a 'good' estimator of $\theta$. And this result is not strong enough to conclude whether something like $\sqrt n(\hat\theta-\theta)$ is asymptotically normal or not. I could not come up with a reasonable argument using CLT either.

So a question remains whether the parent distribution here satisfies the regularity conditions for the limiting distribution of MLE to be normal.

 A: A direct proof for asymptotic normality:
The log-likelihood here is
$$L = -\frac {n \bar x}{\theta} - \theta n \bar y$$
The first and second derivatives are
$$\frac {\partial L}{\partial \theta } = \frac {n \bar x}{\theta^2} - n\bar y,\;\;\;\frac {\partial^2 L}{\partial \theta^2 } = -\frac {2n \bar x}{\theta^3} $$
The MLE $\hat \theta_n$ satisfies
$$\frac {\partial L(\hat \theta_n)}{\partial \theta }=0$$
Applying a mean value expansion around the true value $\theta_0 $ we have
$$\frac {\partial L(\hat \theta_n)}{\partial \theta } = \frac {\partial L(\theta_0)}{\partial \theta } + \frac {\partial^2 L(\tilde \theta_n)}{\partial \theta^2 }(\hat \theta_n - \theta_0) =0$$
for some $\tilde \theta_n$ in between $\hat \theta_n$ and $\theta_0$. Re-arranging we have,
$$(\hat \theta_n - \theta_0) = -\left(\frac {\partial^2 L(\tilde \theta_n)}{\partial \theta^2 }\right)^{-1}\frac {\partial L(\theta_0)}{\partial \theta }$$
But in our single-parameter case, the inverse is just the reciprocal, so, inserting also the specific expressions of the derivatives,
$$(\hat \theta_n - \theta_0) = \frac {\tilde \theta^3_n}{2n\bar x}\left(\frac {n \bar x}{\theta^2_0} - n\bar y\right)$$
$$\implies \sqrt{n}(\hat \theta_n - \theta_0) = \frac {\tilde \theta^3_n}{2\bar x \theta_0^2}\sqrt{n}\cdot\left(\bar x - \theta_0^2\bar y  \right)$$
$$\implies \sqrt{n}(\hat \theta_n - \theta_0) = \frac {\tilde \theta^3_n}{2\bar x \theta_0^2}\cdot\left (n^{-1/2}\sum_{i=1}^n(x_i-\theta_0^2 y_i)\right)$$
The variance of the sum is
$$\text{Var}\left(\sum_{i=1}^n(x_i-\theta_0^2 y_i)\right) = 2n\theta_0^2 $$
Manipulating the expression we can write, using $S_n$ for the sum of the i.i.d. elements,
$$\sqrt{n}(\hat \theta_n - \theta_0) = \left(\frac {\tilde \theta^3_n}{\sqrt{2}\bar x \theta_0}\right)\cdot\frac {\sum_{i=1}^n(x_i-\theta_0^2 y_i)}{\sqrt{n}\sqrt{2}\theta_0} $$
$$\sqrt{n}(\hat \theta_n - \theta_0) = \left(\frac {\tilde \theta^3_n}{\sqrt{2}\bar x \theta_0}\right)\cdot\frac {S_n}{\sqrt{\text{Var}(S_n)}}$$
More over, we have that  $E(x_i-\theta_0^2 y_i) = 0$, so $E(S_n)=0$. So we have the subject matter of a classical CLT, and one can verify that the Lindeberg condition is satisfied. It follows that
$$\frac {S_n}{\sqrt{\text{Var}(S_n)}} \to_d N(0,1)$$
Due to the consistency of the estimator, we also have 
$$\left(\frac {\tilde \theta^3_n}{\sqrt{2}\bar x \theta_0}\right) \to_p \frac{\theta_0}{\sqrt{2}}$$
and by Slutsky's Theorem we arrive at
$$\sqrt{n}(\hat \theta_n - \theta_0) \to_d N (0, \theta_0^2/2)$$
Nice. Double the information, half the variance (compared to the case where we would estimate $\theta_0$ based on a sample from a single random variable).
PS: The fact that in the above expressions $\theta_0$ appears in the denominator, points towards @whuber's comment that MLE's asymptotic normality needs the unknown parameter to be away from the boundary of the parameter space (in our case, away from zero). 
