\(\displaystyle x' = 4x - 13y + 2e^{2t}\cos(3t) \\

y' = x\)

using undetermined coefficients. For the associated homogeneous solution, I get \(\displaystyle \vec{x_h} = Ae^{2t}\begin{bmatrix}13\cos(3t) \\ 2\cos(3t) + 3\sin(3t)\end{bmatrix} + Be^{2t}\begin{bmatrix}13\sin(3t) \\ 2\sin(3t) - 3\cos(3t)\end{bmatrix}\).

My first guess for the particular solution was \(\displaystyle \vec{x_p} = \begin{bmatrix}a_1 \\ a_2\end{bmatrix}e^{2t}\cos(3t) + \begin{bmatrix}b_1 \\ b_2\end{bmatrix}e^{2t}\sin(3t)\).

This was a risky guess, since it closely resembles parts of \(\displaystyle \vec{x_h}\), but I don't know the exact conditions under which the default guess fails when working with systems. Sure enough, I ended up with an inconsistent system for \(\displaystyle a_1\), \(\displaystyle a_2\), \(\displaystyle b_1\), and \(\displaystyle b_2\) as expected.

I looked at the math I used to solve the corresponding second-order equation, \(\displaystyle y'' - 4y' + 13y = 2e^{2t}\cos(3t)\), and saw that the guess \(\displaystyle y_p = ate^{2t}\cos(3t) + bte^{2t}\sin(3t)\) was successful, so my next guess was \(\displaystyle \vec{x_p} = \begin{bmatrix}a_1 \\ a_2\end{bmatrix}te^{2t}\cos(3t) + \begin{bmatrix}b_1 \\ b_2\end{bmatrix}te^{2t}\sin(3t)\).

I used the exponential shift rule (which I'm assuming/hoping is valid for systems?) to get \(\displaystyle e^{2t}(D + 2)\left[\begin{bmatrix}a_1 \\ a_2\end{bmatrix}t\cos(3t) + \begin{bmatrix}b_1 \\ b_2\end{bmatrix}t\sin(3t)\right] = \begin{bmatrix}4 & -13 \\ 1 & 0\end{bmatrix}\left(\begin{bmatrix}a_1 \\ a_2\end{bmatrix}te^{2t}\cos(3t) + \begin{bmatrix}b_1 \\ b_2\end{bmatrix}te^{2t}\sin(3t)\right) + \begin{bmatrix}2e^{2t}\cos(3t) \\ 0\end{bmatrix},\) but this once again led to an inconsistent algebraic system, namely, \(\displaystyle \begin{cases}3b_1 + 2a_1 = 4a_1 - 13a_2 \\ -3a_1 + 2b_1 = 4b_1 - 13b_2 \\ a_1 = 2 \\ b_1 = 0 \\ 3b_2 +2a_2 = a_1 \\ -3a_2 + 2b_2 = b_1 \\ a_2 = 0 \\ b_2 = 0\end{cases}\).

After computing the solution with Wolfram, it appears that the required guess was the sum of the two guesses above, \(\displaystyle \vec{x_p} = \begin{bmatrix}a_1 \\ a_2\end{bmatrix}te^{2t}\cos(3t) + \begin{bmatrix}b_1 \\ b_2\end{bmatrix}te^{2t}\sin(3t) + \begin{bmatrix}c_1 \\ c_2\end{bmatrix}e^{2t}\cos(3t) + \begin{bmatrix}d_1 \\ d_2\end{bmatrix}e^{2t}\sin(3t)\), but this form seems to break the general pattern. For example, if the forcing term in a second-order linear ODE is \(\displaystyle e^t\), and \(\displaystyle y_h = e^t\) is an associated homogeneous solution, then the guess \(\displaystyle y_p = ate^t\) is recommended; you don't need \(\displaystyle y_p = ate^t + be^t\). The corresponding second-order equation above serves as another example. According to this answer, we can conclude that the guess for systems where all but one component of the forcing term are \(\displaystyle 0\) should be the same as if that non-zero component were the forcing term of a single linear ODE.

Why does this system, by contrast, require the original terms to be kept when the new terms with an extra factor of $t$ are introduced? Additionally, why does the corresponding second-order equation, \(\displaystyle y'' - 4y' + 13y = 2e^{2t}\cos(3t)\), provide misleading information about what the guess should be?