Find all polynomials $P$ such that $P(x^2+1)=P(x)^2+1$ Find all polynomials $P$ such that 
$$P(x^2+1)=P(x)^2+1$$ 
 A: Let $P(y)=\sum_{0\le r\le n}a_ry^r$
So, $$P(1+x^2)=\sum_{0\le r\le n}a_r(1+x^2)^r=a_0+a_1(1+x^2)+a_2(1+\binom 21x^2+x^4)+\cdots
+a_{n-1}(1+\binom {n-1}1x^2+\binom {n-1}2x^4+\cdots+\binom {n-1}{n-2}x^{2(n-2)}+\binom {n-1}{n-1}x^{2(n-1)})
+a_n(1+\binom n1x^2+\binom n2x^4+\cdots+\binom n{n-1}x^{2(n-1)}+\binom  n nx^{2n})$$
$$=x^{2n}a_n+x^{2n-2}(a_n\binom n{n-1}+a_{n-1})+x^{2n-4}(a_n\binom n{n-2}+a_{n-1}\binom {n-1}{n-2}+a_{n-2})+x^2(a_n\binom n1+a_{n-1}\binom{n-1}1+\cdots+a_2\binom21+a_1)+\sum_{0\le r\le n}a_r$$
and $$\{P(x)\}^2+1=\{\sum_{0\le r\le n}a_rx^r\}^2+1$$
$$=a_n^2x^{2n}+x^{2n-1}2a_na_{n-1}
+x^{2n-2}(a_{n-1}^2+2a_na_{n-2})+x^{2n-3}2(a_na_{n-3}+a_{n-1}a_{n-2})
+x^{2n-4}(a_{n-2}^2+2a_na_{n-4}+2a_{n-1}a_{n-3})+\cdots+x^2(a_1^2+2a_0a_2)+\sum_{0\le r\le n}a_r^2+1$$
Comparing the coefficients of the different powers of $x$
$r=n\implies a_n=a_n^2\implies a_n=1$ as $a_n\ne0$
$r=n-1\implies 2a_na_{n-1}=0\implies a_{n-1}=0$
$r=n-2\implies a_n\binom n{n-1}+a_{n-1}=a_{n-1}^2+2a_na_{n-2}\implies a_{n-2}=\frac n2$
$r=n-3\implies 2(a_na_{n-3}+a_{n-1}a_{n-2})=0\implies a_{n-3}=0$
$r=n-4\implies a_n\binom n{n-2}+a_{n-1}\binom {n-1}{n-2}+a_{n-2}=a_{n-2}^2+2a_na_{n-4}+2a_{n-1}a_{n-3}\implies a_{n-4}=\frac {n^2}8=\frac1{2!}\left(\frac n2\right)^2$
$r=n-5\implies 2(a_na_{n-5}+a_{n-1}a_{n-4}+a_{n-2}a_{n-3})=0\implies a_{n-5}=0$
$r=n-6\implies 2(a_na_{n-6}+a_{n-1}a_{n-5}+a_{n-2}a_{n-4})+a_{n-3}^2=a_n\binom n{n-3}+a_{n-1}\binom{n-1}{n-3}+a_{n-2}\binom{n-2}{n-3}+a_{n-3}\implies a_{n-6}=\frac1{3!}\left(\frac n2\right)^3-\frac n3$
A: From @lab bhattacharjee post we see that the formulae for $P(1+x^2)$ and $P(x)^2+1$ are
$$
\begin{eqnarray}
P(1+x^2)&=&\sum_{k=0}^{n}{\left(\sum_{i=0}^{n-k}{a_{n-i}\,{{n-i}\choose{k}}}
 \right)\,x^{2\,k}}  \tag{1} \\
P(x)^2+1&=&\sum_{k=0}^{n}{\left(\sum_{i=0}^{k}{a_{n-i}\,a_{n-k+i}}\right)\,x^{
 2\,n-k}} 
+\sum_{k=0}^{n-1}{\left(\sum_{i=0}^{k}{a_{i}\,a_{k-i}}
 \right)\,x^{k}}+1  \tag{2} \\
\end{eqnarray}
$$
The coefficients of the monomials with odd exponents in $(1)$ are $0$. For the coefficients of the exponent $2n$ we get the equation 
$$ a_{n}^2=a_{n}$$
for $n>0$ and therefore 
$$a_{n}=1$$
If $a_{n}=0$ then  it would be $degree(g)<n$. For $n=0$ one gets the equation $a_{0}^2+1=a_{0}$ instead.
The remaining $a_i$ can calculated successively. After calculation $a_n \ldots a_{n-(k-1)}$ for even $k=2m$ the value $a_{n-k}$  can be calculated by the coefficient 
of $x^{2n-2m}$
$$
\sum_{i=0}^{2m}{a_{n-i}\,a_{n-2m+i}}=
\sum_{i=0}^{m}{a_{n-i}\,{{n-i}\choose{n-m}}}
$$
which gives 
$$
a_{n-2m}=\frac{\sum_{i=0}^{m}{a_{n-i}\,{{n-i}\choose{n-m}}}-\sum_{i=1}^{2m-1}{a_{n-i}\,a_{n-2m+i}}}{2a_n} \tag{3}
$$
indexes of the $a$ on the right hand side of this equation are all larger than $2n-2m$ and therefore already known. For odd $k=2m+1$ we got the equation
$$
\sum_{i=0}^{2m+1}{a_{n-i}\,a_{n-2m-1+i}}=0
$$
and therefore
$$
a_{n-2m-1}=-\frac{\sum_{i=1}^{2m}{a_{n-i}\,a_{n-2m-1+i}}}{2a_n} \tag{4}
$$
This shows that $a_{n-2m-1}$ must be $0$. This can be proven by inductions starting with $m=0$ whch gives $a_{n-1}=0$ from $(4)$. In equation $(4)$ one factor $a_j$ of 
each summand has an odd index because 2m-1 is odd  and so $a_i=0$ and also the whole summand by induction. So the whole sum is $0$.
So $0=a_{n-1}=a_{n-3}=...$. For $n$ odd this means especially $a_0=0$.
With the formulae $(3)$ and $(4)$ we can calculate the sequence $a_n, a_{n-1}, \ldots, a_0$ from the coefficients of $x^{2n}, \ldots, x^{n}$. But it is possible that 
$a_i$ calculated do not fullfil the equations for the coefficients of $x^{0}, \ldots, x^{n-1}$. These coefficients give raise to the following equations 
$$
\sum_{i=0}^{2m}{a_{i}\,a_{2m-i}}=\sum_{i=0}^{n-m}{a_{n-i}\,{{n-i}\choose{m}}} \tag{5}
$$
or 
$$ a_0^2+1=\sum_{i=0}^{n} a_i $$
if $k=m=0$
$$
\sum_{i=0}^{2m+1}{a_{i}\,a_{2m+1-i}}=0 \tag{6}
$$
So we can only condlude that there is at most one polynomial for each degree.
The polynomials from the the post of @Hurkyl shows that there exists polynomial with $degree(P)=2^n$ for each $n$.
From @Hurkyl We take the definition of the polynomial 
$$Q(x)=x^2+1$$
and we know that the problem can be restated as
Find all polynomial $P$ that 
$$P(Q(x))=Q(P(x))$$ 
or 
$$ P \circ Q = Q \circ P $$
therefore
$$ P \circ Q^n = Q^n \circ P$$
and @Hurkyl showed that all
$$P=Q^n$$ 
are solutions
I checked for $degree(P)$ from $0$ to $16$ that there are only the following polynomials:
$$
-\frac{\sqrt{3}i-1}{2} \\
\frac{\sqrt{3}i+1}{2} \\
x \\
x^2+1 \\
x^4+2x^2+2 \\
x^8+4x^6+8x^4+8x^2+5 \\
x^{16}+8x^{14}+32x^{12}+80x^{10}+138x^8+168x^6+144x^4+80x^2+26 \\
$$
Besides the constatn solutions the only solutions that exists where the solutions @Hurkyl found.
Let $P$ be a polynomial that fullfills the equation $P * Q = Q *P$. We showed that $a_0$=0 if $degree(P)$ is odd. Therefore $P(0) =0$ and  also
$$P(Q^n(0))=Q^n(P(0))=Q^n(0)$$
So $Q^n(0), n=0,1,2,3,\ldots$ is a strictly increasing and therefore infinite sequence with $P(x)=x$ and therefore $P(x)-x=0$. But if a polynomial $P(x)-x$ is $0$ for 
infinite many values $x$ the the polynomial es equal to $0$.
So $P=id$. THis means that $x$ is the only polynomial with odd degree that satisfies our functional equation.
The remaining problem: Show that $degree(P)=2^n$ if $degree(P)$ is even.
Edit:
The solution fo the remaining problem is  now already included in my other post
A: One solution is $P(x) = x^2 + 1$:
$$ P(x^2 + 1) = (x^2 + 1)^2 + 1 = P(x)^2 + 1$$
Another solution is $P(x) = x^4 + 2x^2 + 2$:
$$ \begin{align} P(x^2 + 1) &= (x^2 + 1)^4 + 2(x^2 + 1)^2 + 2
\\&= x^8 + 4 x^6 + 8 x^4 + 8 x^2 + 5
\\&= (x^4 + 2x^2 + 2)^2 + 1
\\&= P(x)^2 + 1 \end{align}$$
How did I find these? Turns out it's obvious after rewriting the equation! Let $Q(x) = x^2+1$. Then the equation is
$$ P(Q(x)) = Q(P(x)) $$
or more succinctly, $P \circ Q = Q \circ P$. It's now clear that $P=Q$ is one solution, $P = Q \circ Q$ is another solution, $P = Q \circ Q \circ Q$ is yet another, and so forth.
Note that $P(x) = x$ is in this family as well, as the identity function is the empty repeated composition (just like the empty sum is $0$ and the empty product is $1$).
This argument is easily adapted to show the set of solutions is a monoid under composition.
I haven't worked out the complete solution though. Alas I don't know much about the monoid of polynomials under composition.

This answer to a similar question cites a 1922 theorem of Ritt that contains a characterization of polynomials that commute under composition; in particular, we can conclude the repeated compositions of $Q$ are indeed the entire solution space, excluding nonconstant polynomials. The remaining solutions are thus the two functions $P(x) = \beta$ where $\beta$ is a solution to $\beta = Q(\beta)$.
A: It depends on whether the coefficients of the polynomial come from a field (or ring) of characteristic $0$. [ http://en.wikipedia.org/wiki/Characteristic_%28algebra%29 ]
There are extra solutions in characteristic $p$, such as $P(x) = x^p$. It could be a hard problem to determine whether solutions exist that are not generated by compositions of $x^p$ and $x^2+1$.  The general problem of determining all commuting pairs of polynomials over a finite field is well-known and unsolved.
In characteristic $0$, solutions of $P(Q(x))=Q(P(x))$ are classified. The case $Q(x)=x^2+1$ does not fall into one of the families of nontrivial solutions (it is not linearly conjugate to $x^2$ or to the degree 2 Chebyshev polynomial), so that the only possibililty is $P(x) = Q^{\circ \hskip0.7pt n}$ for some $n$.
A: A full characterization of solutions which is quite short and sweet.
Plugging in $x = \pm a$, we see that $P(x) = \pm P(-x)$ for each $x$. Clearly one sign holds infinitely often. By looking at the roots of the polynomial $P(x) - P(-x)$ or $P(x) + P(-x)$, we easily see that $P(x)$ is either odd or even.
If $P(x)$ is odd, then note that $P(0) = 0$. Then $P(1) = 1$, $P(2) = 2$, $P(5) = 5$, etc. and by induction we can get infinitely many values $a$ such that $P(a) = a$. By looking at the roots of the polynomial $P(x) - x$ we see that $P(x) = x$ for all $x$ then.
Now suppose $P$ is even. This means $P(x) = Q(x^2)$ for some polynomial $Q$. Let $R(x) = x^2 + 1$. Remark there clearly exists a polynomial $S(x)$ such that $P(x) = S(x^2+1)$, just shift $Q$. Now remark that $P \circ R = R \circ P$. Hence $S \circ R \circ R = R \circ S \circ R \implies S \circ R = R \circ S$, so $S$ is a solution to the functional equation. By an easy induction on degree, it follows the only even solutions to the equations are iterations of $R$, so we are done. (EDIT : Also there is the constant solution $P(x) = c$ where $c = c^2 + 1$, because the induction starts on degree 1 which I forgot)
A: There is another method (besides that posted here) to prove  that solutions have the structure $\sum_{k=0}^n\,a_{2k}x^{2k}$ without constructing formulas to 
calculate the  coefficients as in my other answer.
We have 
$$P(x)^2-P(-x)^2=(P(x^2+1)-1)-(P((-x)^2+1)-1)=0$$
At least one of $P(x)-P(-x)$ and $P(x)+P(-x)$ must have infinitely many zeros and therefore must be identical to $0$.
Asume that $P(-x)=-P(x)$. This means that $P$ is an odd function and $P(0)=0$. Let us define
$$Q(x)=x^2+1$$
then $P(Q^n(0)=Q^n(0)$ and so 
$$ P(x)-x=0, \quad x=0,\,Q(0),\,Q^2(0),\ldots$$
for infintely many $x$, so $P(x)=x$.
$P$ must be an even polynomial if not $P(x)=x$. But the even polynomials are exactly the polynomial that contain only even powers of x.
The polynomials $Q^n$ are even polynomials that satisfy the functional equation. They have a degree of $2^n$
In contrast to my other answer I was not able to prove with this method that there is at most one polynomial of a certain degree that satisfy the functional equation.

Proof
a polynomial is even $\Longleftrightarrow$ the polynomial contains only even powers of $x$
The even polynomial $P(x)$ can be expressed as sum 
$f(x)+g(x)$ were $f(x)$ is a polynomial that only contains even powers of $x$  and $g(x)$ is a polynomial that contains only odd powers of $x$.
We have 
$$P(x)-f(x)=g(x)$$
The left side is an even polynomial, the right side an odd one, so $P(x)=f(x)$. So polynomials that contain only even powers of $x$ are exactly the even polynomials. 
EDIT:
The following completes the proofs:
Lemma $Q(x)=x^2+1$,$P(Q(x)=Q(P(x))$ and $P(x)$ is even. Then there is a polynomial $T(x)$ with $P=T \circ Q = Q \circ T$
So a polynomial $P$ of degree $n$ can be reduced to a polynomial with degree $\frac{n}{2}$ that also satisfies the functional equation. This process can repeated until we 
arrive at a polynomial with odd degree. But the only polynomial with odd degree that satisfies the functional equation is $T(x)=x$.
Therefore the functional equation 
$$Q \circ P = P \circ Q $$
has exactly the following polynomial solutions:
$$
-\frac{\sqrt{3}i-1}{2} \\
\frac{\sqrt{3}i+1}{2} \\
x \\
Q(x)\\
Q^2(x) \\
Q^3(x) \\
\ldots
$$
Proof of the Lemma
We define $R(x)=\sqrt[+]{x-1}$, then
$$R(Q(x))=Q(R(x))=x,\quad \forall x>1$$
and
$$R(x)^2=x-1, \quad \forall x>1$$
We have 
$$ Q \circ P = P \circ Q$$
and therefore 
$$ Q \circ P \circ R = P \circ Q \circ R = P = P \circ R \circ Q , \quad \forall  x>1$$ 
For $x> 1$ the polynomial 
$$T(x)=P(R(x))=\sum_{k=0}^{n}a_{2n}(x-1)^{n}$$
satisfies 
$$P=T \circ Q = Q \circ T$$
 But if a polynomial equation is satisfied for infinite many $x$ it is satisfied for all $x$.
The coefficients of the polynomials can be efficiently calculated by the formulae $(3)$ and $(4)$ of my other post
