Examples of systems of linear differential equations:
$$ \begin{align} \left. \begin{aligned} 2x'+3x+5y'-y &= e^t \nonumber \\ x'-x+3y'+y &= \sin t \nonumber \end{aligned} \right \} \hskip1em & \text{(a)} \\ \left. \begin{aligned} (D^2 + 3D -1)x + y &= 6 + t^2 \nonumber \\ (D+2)x -(D^2 + D) y &= t \nonumber \end{aligned} \right \} \hskip1em & \text{(b)} \end{align} $$Systems of linear DE with constant coefficients leads readily to solutions by means of operator $\boldsymbol D$.
The above system (a) can be rewritten into
$$ \begin{align} (2D+3)x + (5D-1)y &= e^t, \nonumber \\ (D-1)x + (3D+1)y &= \sin t. \nonumber \end{align} $$The pair of equations
$$ \begin{align} f_1(D)x + g_1(D)x &= h_1(t), \nonumber \\ f_2(D)x + g_2(D)x &= h_2(t)\nonumber, \end{align} $$where $D$ is operator $d/dt$, and coefficients in front of $x$ and $y$ are linear differential operators, is called a system of two linear DE.
We have learnt already that linear differential operator obey all the rules of algebra. So we will use the rules and methods which are used in solving an algebraic system of equations, e. g.:
The above system can be solved by elimination and we will use such approach also to systems of DE, which are written using the operator $D$.
Solve the system of differential equations
$$ \begin{align} \frac{dx}{dt} &= y \nonumber \\ \frac{dy}{dt} &= -x + 2y \nonumber \end{align} $$First, we rewrite the above equations by means of operator $D$:
$$ \begin{align} Dx - y &= 0 \label{ref:pair1a} \\ Dy + x - 2y &= 0 \label{ref:pair1b} \\[2mm] Dx - y &= 0 \hskip2em \nonumber / \cdot(D-2) \\ y(D-2) + x &= 0 \nonumber \nonumber \\[2mm] D(D-2)x - (D-2)y &= 0 \nonumber \\ x + (D-2)y &= 0 \nonumber \nonumber \\[2mm] \end{align} $$Let us add both equations to eliminate $(D-2)y$ and we get
$$ \begin{align} (D^2-2D)x +x &= 0 \nonumber \\ D^2x -2Dx +x &= 0 \nonumber \\ x'' -2x' +x &= 0 \nonumber \\ m^2 -2m +1 &= 0 \nonumber \implies m = \{1, 1\}\\ \end{align} $$So far, we are able to express $x$:
$$ \underline{\underline{x = c_1 e^t + c_2 t e^t}} $$Now we may substitute $x$ into $(\ref{ref:pair1a})$ to get $y$:
$$ \begin{align} & D(c_1 e^t + c_2 t e^t) - y = 0 \nonumber \\ & \underline{\underline{y = c_1 e^t + c_2 t e^t + c_2 e^t}} \nonumber \end{align} $$If we substitute $x,\ y,\ y'$ into $(\ref{ref:pair1b})$ we will find that the equality holds, so above determined $x$ and $y$ are complete solution.
Solve the system of differential equations
$$ \begin{align} 3\frac{dx}{dt} +3x +2y&= e^t \label{ref:pair2a} \\ 4x - 3\frac{dy}{dt} +3y&= 3t \label{ref:pair2b} \end{align} $$First, we rewrite the above equations by means of operator $D$:
$$ \begin{align} 3Dx +3x +2y&= e^t \nonumber\\ 4x - 3Dy +3y&= 3t \nonumber\\[2mm] x(3D+3) + 2y & = e^t \nonumber \hskip2em / \cdot(3D-3)\\ 4x - y(3D-3) &= 3t \nonumber \hskip2em / \cdot 2\\[2mm] x(3D+3)(3D-3) + 2y(3D-3) & = e^t(3D-3) \nonumber \\ 8x - 2y(3D-3) &= 6t \nonumber \\[2mm] \end{align} $$Let us add both last equations to eliminate $2y(3D-3)$:
$$ \begin{align} x(3D+3)(3D-3) + 8x &= 6t \nonumber \\ x(9D^2 -9) + 8x &= 6t \nonumber \\ 9x'' -9x + 8x &= 6t \nonumber \\ 9x'' -x &= 6t \nonumber \\ 9m^2 -1 &= 0 \implies m = \{\frac 1 3, -\frac 1 3\}\nonumber \\ \end{align} $$So far, we are able to collect $x$:
$$ \underline{\underline{x = c_1 e^{1/3t} + c_2 e^{-1/3t} - 6t} } $$Substituting $x'$ into $(\ref{ref:pair2a})$ brings $y$:
$$ \underline{\underline{y = -2c_1 e^{1/3t} -c_2 e^{-1/3t} + 9 + \frac 1 2 e^t + 9t}} $$If we substitute $x,\ y,\ y'$ into $(\ref{ref:pair2b})$ we will observe that the equality holds, so above determined $x$ and $y$ are complete solution.
Solve the system of differential equations
$$ \begin{align} 2\frac{dx}{dt} +\frac{dy}{dt} -x &= e^t \label{ref:pair3a} \\ 3\frac{dx}{dt} + 2\frac{dy}{dt} +y&= t \label{ref:pair3b} \end{align} $$First, we rewrite the above equations by means of operator $D$:
$$ \begin{align} 2Dx + Dy -x &= e^t \nonumber \\ 3Dx + 2Dy + y &= t \nonumber \\[2mm] 2Dx -x + Dy &= e^t \nonumber \hskip2em /\cdot(2D+1) \\ 3Dx + (2D+1)y &= t \nonumber \hskip2em / \cdot(-D)\\[2mm] 2D(2D+1)x - (2D+1)x + D(2D+1)y &= e^t(2D+1) \nonumber \\ -3D^2x - D(2D+1)y &= t(-D)\nonumber \\[2mm] \end{align} $$Let us add both equations to eliminate $D(2D+1)y$.
$$ \begin{align} (4D^2 + 2D)x -(2D+1)x-3D^2x &= e^t(2D+1) -1 \nonumber \\ 4D^2x +2Dx -2Dx -x-3D^2x &= 2De^t +e^t-1 \nonumber \\ D^2x - x &= 2e^t +e^t - 1 \nonumber \\ x''-x &= 3e^t - 1 \nonumber \\ \end{align} $$The solution for $x$ is
$$ \underline{\underline{x(t) = c_1e^t + c_2e^{-t} + \frac 3 2 e^t t + 1}}. $$Let us solve $x''-x = 3e^t - 1$ by annihilator approach:
$$ \begin{align} (D^2-1)x &= 3e^t - 1 \nonumber \\ (D^2-1)x(D-1)D &= (3e^t -1)(D-1)D \nonumber \\ (D^2-1)x(D-1)D &= 0 \nonumber \\ (m^2-1)(m-1)m &= 0 \implies m= \{ 1, -1, 1, 0\} \nonumber \end{align} $$The first two roots belong to complementary function the next two are to form particular solution. So we rewrite constants $c_3,\ c_4$ to more convenient $A,\ B$ and solve them.
$$ \begin{align} x(t) &= c_1e^t + c_2e^{-t} + c_3e^t t + c_4e^0 \nonumber \\ x_p(t) &= Ae^t t + B \nonumber \end{align} $$The constants $A,\ B$ of particular solution have to be found by means of method of undetermined coefficients. For that purpose we need $x_p'(t),\ x_p''(t)$
$$ \begin{align} x_p'(t) &= A(e^t \cdot 1 + e^t t) = Ae^t + Ae^t t \nonumber \\ x_p''(t) &= Ae^t +Ae^t t + Ae^t = 2Ae^t + Ae^t t \nonumber \\ \end{align} $$Now we can substitute $x_p(t),\ x_p''(t)$ into DE $x''-x = 3e^t - 1$ to determine $A,\ B$.
$$ 2Ae^t + Ae^t t - Ae^t t - B = 3e^t - 1 \nonumber \\ e^t(2A) - B = 3e^t - 1 \nonumber \implies 2A =3,\ B = 1 $$Finaly we may collect $x(t)$:
$$ {\underline{x(t) = c_1e^t + c_2 e^{-t} + \frac 3 2 e^t t + 1}}. $$To determine $y(t)$ we can use $(\ref{ref:pair3a})$. For that purpose we have to differentiate $x(t)$ to get $x'(t)$, then from $(\ref{ref:pair3a})$ express $y'(t)$ and integrate $y'(t)$.
$$ \begin{align} x'(t) &= c_1e^t - c_2e^{-t} + \frac 3 2 (e^t t + e^t) \nonumber \\ y'(t) &= e^t(-c_1-2) + e^{-t}(3c_2) - e^t t(\frac 3 2) + 1 \nonumber \\ y(t) &= e^t(-c_1-\frac 1 2) - e^{-t}(3c_2) - \frac 3 2 e^t t + t + \color{red}{d_1} \nonumber \end{align} $$The solutions $x(t),\ y(t)$ have to comply both equations $(\ref{ref:pair3a}),\ (\ref{ref:pair3b})$. The constant $d_1$ can be determined again by method of undetermined coefficients from $(\ref{ref:pair3b})$.
$$ \underline{\underline{y(t) = e^t(-c_1-\frac 1 2) - e^{-t}(3c_2) - \frac 3 2 e^t t + t -2}} $$If we substitute $x,\ x',\ y,\ y'$ into system of DE, it will confirm that the equalities hold.
Solve the system of differential equations
$$ \begin{align} \frac{dx}{dt} &= 3e^{-t} \nonumber \\ \frac{dy}{dt} &= x + y \nonumber \end{align} $$First, rewrite the above equations by means of operator $D$ and then we solve as an basic algebra task (to eliminate $x$).
$$ \begin{align} Dx &= 3e^{-t} \nonumber \\ Dy &= x + y \nonumber\\[2mm] Dx &= 3e^{-t} \nonumber \\ Dy - x - y &= 0 \hskip2em /\cdot D \label{ref:pair40} \\[2mm] Dx &= 3e^{-t} \nonumber \\ D(Dy - x- y) &= 0 \nonumber \\[2mm] Dx &= 3e^{-t} \label{ref:pair4a} \\ D^2y - Dx - Dy &= 0 \label{ref:pair4b} \\[2mm] D^2y - Dy &= 3e^{-t} \nonumber \hskip2em \text{Note: (\ref{ref:pair4a}) + (\ref{ref:pair4b})} \end{align} $$Finally we have to solve $y'' - y' = 3e^{-t}$ which is solved by means of auxiliary equation $(m^2 -m) = 0 \implies m = \{0, 1\}$, so we can express $y_c$
$$ y_c = c_1 + c_2e^t $$and $y_p$ has to be determined. That can be accomplished by the method of undetermined coefficients, because considering the members on the right side ($e^{-t}$) we can expect particular solution containing $e^{-t}$ and its derivatives:
$$ \left. \begin{aligned} y_p &= Ae^{-t} \nonumber \\ y_p' &= -Ae^{-t} \nonumber \\ y_p'' &= Ae^{-t} \nonumber \\ \end{aligned} \right\} \hskip2em Ae^{-t}+Ae^{-t} = 3e^{-t} \implies A = \frac 3 2 $$Then
$$ \underline{\underline{y = c_1 + c_2e^t + \frac 3 2 e^{-t}}}. $$Now there are more ways how to find $x$. We can use $(\ref{ref:pair40})$. For that purpose we have to differentiate $y$.
$$ \begin{align} &Dy = c_2 e^t - \frac 3 2 e^{-t} \nonumber \\ &c_2e^t - \frac 3 2 e^{-t} - x - (c_1 + c_2e^t + \frac 3 2 e^{-t}) = 0 \nonumber \\ &e^t(c_2 - c_2) + e^{-t}(-\frac 3 2 - \frac 3 2) - c_1 = x \nonumber \\ \end{align} $$When simplified we have $x$
$$ \underline{\underline{x = -c_1 - 3e^{-t}}} $$Since we are solving system of DE, it is good idea to test the solutions against both given DE to ensure that no part of solution has been lost.