What would be the distribution of the following equation:
$$y = a^2 + 2ad + d^2$$
where $a$ and $d$ are independent non-central chi-square random variables with $2 \textbf{M}$ degrees of freedom.
OBS.: The r.v.'s generating both $a$ and $d$ have $\mu = 0$ and $\sigma^2 \neq 1$, let's say $\sigma^2 = c$.
Best Answer
If $a, d\sim\chi^2_{2M}$ are independent, then $X=a+d$ will have $\chi^2_{4M}$ distribution. Since $X$ is non-negative, CDF of $Y=a^2+2ad+d^2=(a+d)^2=X^2$ can be found by noting $$F_Y(y)=P(Y\leq y)=P(X^2\leq y)=P(X\leq \sqrt{y})=F_X(\sqrt{y}).$$ Therefore, $$f_Y(y)=\frac{1}{2\sqrt{y}}f_X(\sqrt{y})=\frac{1}{2^{2M+1}\Gamma(2M)}y^{M-1}e^{-\sqrt{y}/2}.$$
If $a$ and $d$ are correlated then things are much more intricate. See for example N. H. Gordon & P. F. Ramig's Cumulative distribution function of the sum of correlated chi-squared random variables (1983) for a definition of multivariate chi-squared and distribution of its sum.
If $\mu\neq 2M$ then you are dealing with non-central chi-squared so the above will no longer be valid. This post may provide some insight.
EDIT: Based on the new information it seems $a$ and $d$ are formed by summing up normal r.v. with non-unit variance. Recall if $Z\sim N(0, 1)$ then $\sqrt{c}Z\sim N(0, c)$. Since now $$a=c\sum_{i=1}^{2M}Z_i^2=d,$$ both $a,d$ will have chi-squared distribution scaled by $c$, i.e. $\Gamma(M, 2c)$ distribution. In this case $X=a+d$ will be $\Gamma(2M, 2c)$ distributed. As a result, for $Y=X^2$ we have $$f_Y(y)=\frac{1}{2(2c)^{2M}\Gamma(2M)}y^{M-1}e^{-\sqrt{y}/2c}.$$