$\begin{cases}x[n + 1] = -x[n] + 3 \\ y[n + 1] = -y[n] – e^{-n}\end{cases}$
Although I realize this system happens to be decoupled, I want to solve it using a general technique which doesn't depend on this contingency. I first solve the homogeneous system, which has coefficient matrix $\begin{bmatrix}-1 & 0 \\ 0 & -1\end{bmatrix}$, which in turn has eigenvalue $-1$ with all of $\mathbb{R}^2$ as its eigenspace, so the homogeneous solution should be $\vec{x_h} = A(-1)^n\begin{bmatrix}1 \\ 0\end{bmatrix} + B(-1)^n\begin{bmatrix}0 \\ 1\end{bmatrix}$.
Based on how undetermined coefficients works for systems of ODEs, I guess that $\vec{x_p} = \vec a + \vec b e^{-n}$ will be the form of the particular solution. Then the system becomes $\vec a + \vec b e^{-(n + 1)} = \begin{bmatrix}-1 & 0 \\ 0 & -1\end{bmatrix}(\vec a + \vec b e^{-n}) + \begin{bmatrix}3 \\ -e^{-n}\end{bmatrix}$. Expanding, $\begin{bmatrix}a_1 + \frac{b_1}{e} e^{-n} \\ a_2 + \frac{b_2}{e} e^{-n}\end{bmatrix}$ = $\begin{bmatrix}3 – a_1 – b_1 e^{-n} \\ -a_2 – b_2 e^{-n} – e^{-n}\end{bmatrix}$, leading to $\begin{cases} a_1 = 3 – a_1 \\ \frac{b_1}{e} = -b_1 \\ a_2 = -a_2 \\\frac{b_2}{e} = -b_2 – 1\end{cases}$ with solution $a_1 = \frac{3}{2}, b_1 = 0, a_2 = 0, b_2 = \frac{e}{e + 1}$, and finally the general solution $\vec x = A(-1)^n\begin{bmatrix}1 \\ 0\end{bmatrix} + B(-1)^n\begin{bmatrix}0 \\ 1\end{bmatrix} + \begin{bmatrix}\frac{3}{2} \\ \frac{1}{e^n + e^{n – 1}}\end{bmatrix}$.
I couldn't make this match Wolfram's answer exactly, and the graphs aren't identical. There isn't much information on this topic on Youtube, etc., so I just did what felt natural based on solving matrix ODEs.
Does this undetermined coefficients technique from matrix ODEs transfer over to matrix difference equations in the above manner, as one would expect? If so, did I make the correct particular solution guess (a sum of the forms one would guess for each individual equation in the system)? I wasn't sure if the $3$ forcing term would "resonate" with the $(-1)^n$ characteristic mode, requiring an extra term in the guess, but $(-1)^n$ technically isn't constant, and the coefficient system turned out consistent, so I suppose this should be unnecessary. Where did I go wrong?
Best Answer
$\def\A{{\bf A}} \def\B{{\bf B}} \def\C{{\bf C}} \def\M{{\bf M}} \def\I{{\bf I}} \def\O{{\bf 0}} \def\x{{\bf x}} \def\y{{\bf y}}$There is a minor sign error in your work. Note that $b_2 = -e/(1+e)$. Taking this into account, you will arrive at the correct solution. Below we consider the technique in a fairly general way.
Consider a vector of discrete functions $\y_n$ such that \begin{align*} \y_{n+1}-\y_n = \M\y_n,\tag{1} \end{align*} where $\M$ does not depend on $n$. For example, if $$\y_n = \left[\begin{array}{cccccc} 1&n&n^2&n^3&n^4&\cdots \end{array}\right]^T$$ then such an $\M$ exists, $$\M = \left[\begin{array}{cccccc} 0&0&0&0&0&\cdots \\ 1&0&0&0&0&\cdots \\ 1&2&0&0&0&\cdots \\ 1&3&3&0&0&\cdots \\ 1&4&6&4&0&\cdots \\ \vdots& & & & &\ddots \end{array}\right]$$ If $$\y'_n = \left[\begin{array}{cccccc} r^n&nr^n&n^2r^n&n^3r^n&n^4r^n&\cdots \end{array}\right]^T$$ then $$\M' = \left[\begin{array}{cccccc} r-1&0&0&0&0&\cdots \\ r&r-1&0&0&0&\cdots \\ r&2r&r-1&0&0&\cdots \\ r&3r&3r&r-1&0&\cdots \\ r&4r&6r&4r&r-1&\cdots \\ \vdots& & & & &\ddots \end{array}\right]$$ We can consider combinations of such discrete functions such as $$\y_n = \left[\begin{array}{cccc} 1&n&r^n&nr^n \\ \end{array}\right]^T$$ which results in a block diagonal form for $\M$, $$\M = \left[\begin{array}{cccccc} 0&0&0&0 \\ 1&0&0&0 \\ 0&0&r-1&0 \\ 0&0&r&r-1 \end{array}\right]$$
Suppose that we wish to solve \begin{align*} \x_{n+1} &= \A\x_{n} + \B\y_n\tag{2} \end{align*} where $\A,\B$ are constant matrices and where $\y_n$ is a discrete function satisfying (1). We search for homogeneous and particular solutions, $\x_n^h,\x_n^p$, such that \begin{align*} \x_n &= \x_n^h + \x_n^p. \end{align*} Clearly $\x_n^h = \A^n \x_0^h$. Assume that $\x_n^p = \C\y_n$ where $\C$ is independent of $n$. This implies that $$(\C\M+(\I-\A)\C-\B)\y_n = \O.$$ We assume that the elements of $\y_n$ are independent discrete functions of $n$. Thus, we must solve \begin{align*} \C\M+(\I-\A)\C&=\B.\tag{3} \end{align*} This can be solved using the vec trick, for example. A unique solution will exist if $\M$ and $\A-\I$ share no eigenvalues.
For the original question we have \begin{align*} \y_n &= \left[\begin{array}{c} 1 \\ e^{-n} \end{array}\right] \\ \M &= \left[\begin{array}{cc} 0&0\\ 0&\frac 1 e -1 \end{array}\right] \\ \A &= -\I \\ \B &= \left[\begin{array}{cc} 3&0\\ 0&-1 \end{array}\right] \end{align*} Equation (3) becomes \begin{align*} \left[\begin{array}{cc} 2c_{11} & \frac{1+e}{e}c_{12} \\ 2c_{21} & \frac{1+e}{e}c_{22} \end{array}\right] &= \left[\begin{array}{cc} 3&0\\ 0&-1 \end{array}\right] \end{align*} which is easily solved, with the result \begin{align*} \C &= \left[\begin{array}{cc} \frac32 & 0 \\ 0 & -\frac{e}{1+e} \end{array}\right] \end{align*} Thus, \begin{align*} \x_n &= \left[\begin{array}{cc} -1 & 0 \\ 0 & -1 \end{array}\right]^n \left[\begin{array}{c} x_0^h \\ y_0^h \end{array}\right] + \left[\begin{array}{cc} \frac32 & 0 \\ 0 & -\frac{e}{1+e} \end{array}\right] \left[\begin{array}{c} 1 \\ e^{-n} \end{array}\right] \\ &= \left[\begin{array}{c} (-1)^n x_0^h + \frac32 \\ (-1)^n y_0^h - \frac{e}{1+e}e^{-n} \end{array}\right] \end{align*}