Choosing random variables $X \sim U(-1,1)$ and $Y = X^2$, it is possible to show that X and Y are uncorrelated, yet not independent.
I was wondering what the probability distribution of $Y$ now looks like. We should have $$P_Y(y) = \int_{-1}^{1}dx~~\frac{1}{2} \delta(x^2-y),$$ where $\delta$ denotes Dirac's delta distribution.
My first intuition would be to go for a substitution via $u = x^2$, but this is bad for the domain of integration. So to circumvent this inconvenience, I would take the integration domain into an indicator function $\chi_{[-1,1]}(x)$. Concatenating those ideas, I have $$ \frac{1}{2}\int \frac{du}{2 \sqrt{u}}~~\biggr(\chi_{[-1,0]}(- \sqrt{u}) + \chi_{[0,1]}( \sqrt{u})\biggr) \delta(u-y)~=~\frac{1}{4\sqrt{y}}\biggr(\chi_{[-1,0]}(- \sqrt{y}) + \chi_{[0,1]}( \sqrt{y})\biggr).$$
Are parts of this solution correct? How should one proceed to arrive at the correct solution?
You forgot to do the change of variables $du= 2x dx$.
In any case, working with the CDF of $Y$ avoids dealing with the Dirac delta and the degenerate joint distribution.
One can show $$P(Y \le y) = P(X^2 \le y) = \begin{cases} 0 & y \le 0 \\ \sqrt{y} & 0 < y < 1 \\ 1 & y \ge 1\end{cases}$$ since in the case $y \in (0,1)$ we have $P(X^2 \le y) = P(-\sqrt{y} \le X \le \sqrt{y})$.
Now that you have the CDF of $Y$, you can find the PDF by differentiating.