Solving a system of differential equations I would like to get some help by solving the following problem:
$$ 
p'_1= \frac 1 x p_1 - p_2 + x$$
$$ p'_2=\frac 1{x^2}p_1+\frac 2 x p_2 - x^2 $$
with initial conditions $$p_1(1)=p_2(1)=0, x \gt0 $$
EDIT:
If I use Wolframalpha, I get


Where $u$ and $v$ are obviously $p_1$ and $p_2$.
Can anybody explain whats going on?

One approach is to express $p_2=-p_1'+\frac1xp_1+x$ from the first equation and substitute into the second to get:
$$-p_1''+\frac1xp_1'-\frac1{x^2}p_1+1=p_2'=\frac1{x^2}p_1+\frac2x\left(-p_1'+\frac1xp_1+x\right)-x^2$$
$$p_1''-\frac3xp_1'+\frac4{x^2}p_1=x^2-1$$
Multiplying by $x^2$ we get $x^2p_1''-3xp_1'+4p_1=x^4-x^2$, which is a Cauchy–Euler equation.
Solve the homogeneous equation $x^2p_1''-3xp_1'+4p_1=0$ first: $p_1=x^r$, hence 
$$x^r(r(r-1)-3r+4)=0\hspace{5pt}\Rightarrow\hspace{5pt} r^2-4r+4=(r-2)^2=0 \hspace{5pt}\Rightarrow\hspace{5pt} r=2$$
So the solution to the homogeneous equation is $C_1x^2\ln x+C_2x^2$.
Now we can use Green's function of the equation to find $p_1$: ($y_1(x)=x^2\ln x,\hspace{3pt} y_2(x)=x^2$)
$$\begin{align*}k(x,t)&=\frac{y_1(t)y_2(x)-y_1(x)y_2(t)}{y_1(t)y_2'(t)-y_2(t)y_1'(t)}=
\frac{t^2\ln t\cdot x^2-x^2\ln x\cdot t^2}{t^2\ln t\cdot 2t-t^2\cdot(2t\ln t+t)}=\frac{t^2\ln t\cdot x^2-x^2\ln x\cdot t^2}{-t^3}\\
&=\frac{x^2\ln x-x^2\ln t}{t}\end{align*}$$
Then ($b(x)$ is the in-homogeneous part, i.e. $b(x)=x^4-x^2$)
$$\begin{align*}p_1(x)&=\int k(x,t)b(t)dt=\int \frac{x^2\ln x-x^2\ln t}{t}t^2(t^2-1)dt\\
&=x^2\ln x\int t(t^2-1)dt-x^2\int t(t^2-1)\ln tdt\end{align*}$$
Compute the integral, find $p_1$ using you initial values and the substitute back to find $p_2$.
