A friend came up with this problem, and we and a few others tried to solve it. It turned out to be really hard, so one of us asked his professor. I came with him, and it took me, him and the professor about an hour to say anything interesting about it.
We figured out that for positive $x$, assuming $f$ exists and is differentiable, $f$ is monotonically increasing. (Differentiating both sides gives $f'(x)*[\text{positive stuff}]=2x$). So $f$ is invertible there. We also figured out that f becomes arbitrarily large, and we guessed that it grows faster than any linear function. Plugging in $f{-1}(x)$ for $x$ gives $x+f(x)=[f^{-1}(x)]^2$. Since $f(x)$ grows faster than $x$, $f^{-1}$ grows slower and therefore $f(x)=[f^{-1}(x)]^2-x\le x^2$.
Unfortunately, that's about all we know... No one knew how to deal with the $f(f(x))$ term. We don't even know if the equation has a solution. How can you solve this problem, and how do you deal with repeated function applications in general?

Here is something: Iterating
$$f_0(x):=0, \quad f_{n+1}(x):= x^2-f_n\bigl(f_n(x)\bigr)$$
one obtains a sequence of polynomials whose lowest terms stabilize at
$$x^2-x^4+2x^6-4x^8+8x^{10}-16 x^{12}+32 x^{14}-65 x^{16}+\ldots-17316 x^{28} +\ldots\ .$$
The sequence of coefficients so obtained (apart from the sign) is listed here: OEIS, but I can't make anything out of the explanation given there.