Coupled second-order differential equations I am trying to solve the following system of coupled ODEs:
\begin{align}
-x^2 f'' - 3xf' + (1-2a)f - (a+1)x^2g'' + (2-4a)xg' + (4a-2)g &= 0,\\
(a-1)x^2 f'' + (4a+2)xf' + (12-6a)f + 12xg' + (12a-24)g &= 0,
\end{align}
where $f$ and $g$ are function of $x$ and $a$ is a constant.
What method do you suggest for solving this system? Any suggestion will be appreciated!
Thanks!
 A: The idea:


*

*Following the @Michael Galuza suggestion you will get a system of two linear equations with constant coefficients wrt functions $F(t)=f(e^t)=f(x)$ and $G(t)=g(e^t)=g(x)$
$$
\left[\matrix{-1 & -(a+1)\\a-1 & 0}\right]\left[\matrix{F''\\G''}\right]+
\left[\matrix{-2 & 3(1-a)\\3(1+a) & 12}\right]\left[\matrix{F'\\G'}\right]+
\left[\matrix{1-2a & 4a-2\\12-6a & 12a-24}\right]\left[\matrix{F\\G}\right]=
\left[\matrix{0\\0}\right]
$$
that I write shortly as
$$
AH''+BH'+CH=\left[\matrix{0\\0}\right].
$$

*The matrix $A=\left[\matrix{-1 & -a-1\\a-1 & 0}\right]$ is non-singular for $a\ne\pm 1$. So you can left-multiply with $A^{-1}$ to get
$$
H''=DH'+EH,
$$
which with the notation $Z=\left[\matrix{H\\H'}\right]$ is equivalent to the linear system
$$
Z'=\left[\matrix{0 & I\\ E & D}\right]Z.
$$

*If $A$ is singular then the equations are dependent, and the system can be reduced to a first order DE after eliminating the coefficients for the second derivatives. For example, if $a=1$ we get the second equations as
$$
\left[\matrix{1 & 2}\right]H'+\left[\matrix{1 & -2}\right]H=0.\tag1
$$
Adding it to the first equation gives
$$
\left[\matrix{-1 & -2}\right]H''+\left[\matrix{-1 & 2}\right]H'=0\tag2
$$
which is, in fact, dependent (minus the derivative of $(1))$. Solving the equation $(1)$ gives you the dependence between $F$ and $G$ (the whole subspace of solutions)
$$
F'+2G'+F-2G=0\quad\Rightarrow\quad (F+2G)'+(F+2G)=4G\quad\Rightarrow\quad
$$
$$
\Rightarrow\quad F(t)+2G(t)=4\int e^{s-t}G(s)\,ds.
$$
Similar for $a=-1$.

A: For any equation of the form
$$
c_n x^n u^{(n)}(x) + \dots + c_1 x u'(x) + c_0 u(x) = 0,
$$
a natural trial solution will be $u = x^\alpha$ for some $\alpha$.  Let's see if this works by looking for a solution of the form
$$
\begin{bmatrix} f \\ g \end{bmatrix} = \begin{bmatrix} C  \\ D \end{bmatrix} x^\alpha. 
$$
Plugging these in to the above equations, we get
\begin{align}
(-\alpha(\alpha - 1) - 3 \alpha + (1 - 2a))C &&+ (-\alpha(\alpha - 1)(a + 1) + (2 - 4a) \alpha + (4a - 2)) D &= 0 \\
((a - 1)\alpha(\alpha - 1) + (4a + 2) \alpha + (12 - 6a))C &&+ (12 \alpha + (12a - 24)) D &= 0,
\end{align}
or
$$
\underbrace{\begin{bmatrix}-\alpha(\alpha - 1) - 3 \alpha + (1 - 2a) & -\alpha(\alpha - 1)(a + 1) + (2 - 4a) \alpha + (4a - 2) \\ (a - 1)\alpha(\alpha - 1) + (4a + 2) \alpha + (12 - 6a) & 12 \alpha + (12a - 24) \end{bmatrix}}_{\displaystyle{ \equiv M}} \begin{bmatrix} C \\ D \end{bmatrix} = 0
$$
There will be a non-trivial solution to this equation if and only if the determinant of this equation is zero;  working through this, the determinant is
$$
\det M = (a^2 - 1) \alpha (\alpha - 2) (\alpha + 3) (\alpha + 5).
$$
Thus, assuming $a \neq \pm 1$, our solutions will have either $\alpha = \{2,0,-3,-5\}.$  The relative values of $C$ and $D$ will then be found by plugging in each of the above values into the above matrix and finding its resulting null space.  For example, if $\alpha = 2$, we get
$$
\begin{bmatrix}-7 -2a & -6a \\ 14 + 4a & 12 a \end{bmatrix} \begin{bmatrix} C \\ D \end{bmatrix} = 0,
$$
which implies that $(7 + 2a) C + 6aD = 0$.  So one of the linearly independent solutions will be
$$
\begin{bmatrix} f \\ g \end{bmatrix} = \begin{bmatrix} 6a x^2 \\ -(7 + 2a)x^2 \end{bmatrix}.
$$
The other three linearly independent solutions will be found by applying the other three values of $\alpha$ to the matrix and finding the corresponding null spaces. 
Unfortunately, it's not immediately clear to me how to extend this method to the cases $a = \pm 1$.  I'll have to ruminate on this further, and I'll edit this post as necessary if I see a way to do it.
