Solution to non homogeneous system of ODE. This is the system:
$$ y_1^{'} = 2y_1+y_2+xe^x \\ y_2^{'}=3y_1 + 4y_2 + 2e^x$$
Solution to this system is sum of solution to the corresponding homogeneous system and one particular solution of this system. So when i find solution of homogeneous system i can use variation of parameters method (and i am supposed to, since text of the problem requires it even though i'm not very smooth with this method), but the problem is, when i try to solve homogeneous system:
$$ y_1^{'} = 2y_1+y_2\\ y_2^{'}=3y_1 + 4y_2 $$
i have the following:
firstly i find characteristic polynomial and it's $$\lambda^2-6\lambda+5$$ so i have characteristic numbers $$\lambda_1=1$$ and $$\lambda_2=5$$
Now, for $\lambda_1$ i end up with the following system 
$$ y_1^{'} = y_1+y_2\\ y_2^{'}=3y_1 + 3y_2 $$ and since solutions are following type of functions $$y=ae^x$$ where $a$ is constant i end up with this algebraic system of equations $$a_1=a_1+a_2 \\ a_2=3a_2+3a_1$$ and it has only trivial solution, and similar thing happens for the second characteristic number, which made me stuck here, even before what i thought that would be the hard part of the problem, and that is the method of variations of parameters. Any ideas on how to solve this?
 A: If I am understanding, you are trying to solve the homogeneous part first.
For $\lambda_1 = 1$, the RREF of $[A - \lambda_1 I] v_1 = [A - I] v_1 = 0$ gives
$$\begin{bmatrix}1 & 1 \\0 & 0\end{bmatrix} v_1 = 0 \implies v_1 = \begin{bmatrix} - 1 \\ 1 \end{bmatrix}$$
RREF details for $\lambda_1$
$$[A - I] v_1 = \begin{bmatrix}1 & 1 \\3 & 3\end{bmatrix} v_1 = 0$$


*

*Divide row $2$ by $3$:


$$\begin{bmatrix}1 & 1 \\1 & 1\end{bmatrix} v_1 = 0$$


*

*Set row $2 =$ row $2 - $ row $1$:


$$\begin{bmatrix}1 & 1 \\0 & 0\end{bmatrix} v_1 = 0$$


*

*We can now solve for $v_1$:


$$v_1 = \begin{bmatrix} - 1 \\ 1 \end{bmatrix}$$
Repeating this process for $\lambda_2 = 5$, the RREF of $[A - \lambda_2 I] v_2 = [A - 5I] v_2 = 0$ gives
$$\begin{bmatrix}1 & -\dfrac{1}{3} \\0 & 0\end{bmatrix} v_2 = 0 \implies v_2 = \begin{bmatrix} 1 \\ 3\end{bmatrix}$$
I will assume you can continue here to write the homogeneous solution.
A: Since exp(x) appears in both forcing terms, it makes sense to factor it out:
y = g exp(x)
Under which the equations become
g1' = g1 + g2 + x
g2'= 3( g1 + g2 ) +2
Add them up, with g= g1+g2
g' = 4 g +x +2
Solution is simple... 
g= C exp(4x) - x/4 - 9/16
Then get g2 by integration
g2' = 3 g + 2
Or g2 = 3C/4 + 5/16 x - 3/8 x^2 + D
And g1= g-g2
If the purpose of the problem was to get you to use the general formalism of characteristic equation etc, this doesn't help much. More often than not I find that formalism gets people into trouble on text book problems.
A: From the second equation, one has
$$  y_1=\frac13y_2^{'}-\frac43y_2-\frac23e^x. \tag{1}$$
Differentiating the second equation, one has
$$ y_2''=3y_1' + 4y_2' + 2e^x.\tag{2}$$
Using (1) and the first equation, one has
$$ y_2''=3(2y_1+y_2+xe^x)+4y_2' + 2e^x=6(\frac13y_2^{'}-\frac43y_2-\frac23e^x)+3y_2+3xe^x+4y_2' + 2e^x$$
or
$$ y_2''-6y_2'+5y_2=(3x-2)e^x.\tag{3}$$
Note that $y_2''-6y_2'+5y_2=0$ has two independent solutions $e^x,e^{5x}$ with the Wronskian
$$ W(e^x,e^{5x})=\left|\begin{matrix}e^x&e^{5x}\\e^x&5e^{5x}\end{matrix}\right|=4e^{6x}$$
and that (3) has a particular solution
$$ y_p=-e^x\int\frac{e^{5x}(3x-2)e^x}{W(e^x,e^{5x})}dx+e^{5x}\int\frac{e^{x}(3x-2)e^x}{W(e^x,e^{5x})}dx$$
which gives
$$ y_2=\frac{1}{64}(-24x^2+20x+5)e^x.$$
Thus the solution of (3) is
$$ y_2=\frac{1}{64}(-24x^2+20x+5)e^x+C_1e^x+C_2e^{5x} $$ 
and hence one has from (1)
$$ y_1=\frac{1}{192}(72x^2-108x-123)e^x-C_1e^x+\frac13e^{5x}.$$
