solving $X^2 - 3X - A = 0$ where $A,X \in \mathbb{M_2(\mathbb{R})}$ Given $A = \begin{pmatrix} 7 & 3 \\ 3 & 7 \end{pmatrix}$ find a $2\times 2$ matrix $X$ s.t. $X^2 - 3X - A = 0$, in the previous parts I have diagonalised $A$ and got $P^{-1}AP = \begin
{pmatrix} 4 & 0 \\ 0 & 10 \end{pmatrix}$, but I'm not sure how this relates to finding $X$.
An idea was that because $A$ has $7,3$ repeating then $X = \begin{pmatrix} a & b \\ b & a \end{pmatrix}$ but this was only an idea and I'm not sure why $X$ would be in that form.
 A: $$X^2-3X-A=0$$
$$P(X^2-3X-A)P^{-1}=0$$
Which using the distruibitive property and the fact that $PX^2P^{-1}=(PXP^{-1})^2$, yields:
$$(PXP^{-1})^2-3(PXP^{-1})-PAP^{-1}=0$$
Now complete the square:
$$(PXP^{-1})^2-3(PXP^{-1})+\frac{9}{4}I_2=\frac{9}{4}I_2+PAP^{-1}$$
$$(PXP^{-1}-\frac{3}{2}I_2)^2=\frac{9}{4}I_2+PAP^{-1}= \begin{pmatrix} \frac{9}{4} & 0 \\ 0 & \frac{9}{4} \end{pmatrix}+\begin{pmatrix} 4 & 0 \\ 0 & 10 \end{pmatrix}$$
$$(PXP^{-1}-\frac{3}{2}I_2)^2= \begin{pmatrix} 6.25 & 0 \\ 0 & 12.25 \end{pmatrix}=\begin{pmatrix} 2.5 & 0 \\ 0 & 3.5 \end{pmatrix}^2$$
Thus, one of the solutions is obtained by solving 
$$PXP^{-1}-\frac{3}{2}I_2= \begin{pmatrix} 2.5 & 0 \\ 0 & 3.5 \end{pmatrix}$$
Which is easy to solve.
To determine all solutions you just need to find all square roots of the diagonal matrix:$$\begin{pmatrix} 6.25 & 0 \\ 0 & 12.25 \end{pmatrix}$$
I put one of four square roots I know, I am not sure if there are others. I think you can easily determine this
A: Let us deal with the general case first. Suppose $A=PDP^{-1}$ be a diagonalisation such that all eigenvalues of $A$ (diagonal entries of $D$) are distinct. If $p(X)=A$ for some matrix polynomial $p(X)$, then $p(P^{-1}XP)=D$. So, if we define $Y=P^{-1}XP$, we are essentially solving $p(Y)=D$.
$D$ is a polynomial in $Y$, so the two matrices commute, i.e. $DY=YD$. However, $D$ has distinct eigenvalues. Therefore $Y$ must be diagonal. It follows that if $D=\operatorname{diag}(d_1,\ldots,d_n)$ and $Y=\operatorname{diag}(y_1,\ldots,y_n)$, then all solutions of $p(Y)=D$ are obtained by solving $p(y_i)=d_i$ individually for each $i$.
As $p(y_i)=d_i$ and the $d_i$s are distinct, the $y_i$s in each $Y$ are distinct too. Therefore, if $q$ is the Lagrange interpolating polynomial for the data points $(y_1,d_1),\ldots,(y_n,d_n)$ ($y$ is abscissa and $d$ is ordinate), we have $Y=q(D)$ and $X=PYP^{-1}=q(PDP^{-1})=q(A)$. In other words, every solution $X$ is a polynomial in $A$.
Now, let's get to your specific case. The $A$ in your question has the form $\pmatrix{a&b\\ b&a}$. By symmetry, all powers of $A$ are also of the same form. Since $X$ is a polynomial in $A$, it is a linear combination of the powers in $A$. Hence $X$ must bear the form of $\pmatrix{a&b\\ b&a}$ too.
