My question comes from section 8.4 of Analysis of Numerical Methods by Isaacson and Keller. I am trying to get a better understanding of the formula stated in Theorem 3, p.409. First, some definitions and basic facts:
Definition 1. An $n$th order inhomogeneous linear difference equation is a system of equations of the form
$$ \hspace{2cm} \sum_{s=0}^{n} a_s u_{j+s} = c_{j+n}, \quad j = 0,1,2,\ldots \hspace{2cm} (1) $$
Here, the coefficients $a_0,a_1,\ldots,a_n$ and $c_n,c_{n+1},c_{n+2},\ldots$ are complex numbers with $a_n \neq 0$; the sequence $\{u_j\}_{j=0}^{\infty}$ is complex-valued, and is to be determined based on the given coefficients. In the special case where $c_{j+n} = 0 $ for all $j \geq 0$, i.e. when
$$ \hspace{2.3cm} \sum_{s=0}^{n} a_s u_{j+s} = 0, \quad j = 0,1,2,\ldots \hspace{2.3cm} (2) $$
the system is called homogeneous.
Proposition 1. The set of all complex-valued sequences $\{u_j \}_{j=0}^{\infty}$ satisfying the homogeneous equation (2) forms an $n$-dimensional vector space over $\mathbb{C}$.
Definition 2. A basis for the solution space of (2) is called a fundamental set of solutions.
Proposition 2. For each $(v_0,v_1,\ldots,v_{n-1}) \in \mathbb{C}^n$, the homogeneous initial value problem
$$
\hspace{1.8cm}
\begin{cases}
\displaystyle \sum_{s=0}^{n} a_s u_{j+s} = 0, \quad j = 0,1,2,\ldots \\[5pt]
u_0 = v_0, \; u_1 = v_1, \ldots, u_{n-1} = v_{n-1}
\end{cases}
\hspace{1.8cm} (3)
$$
has a unique solution.
Now here is the theorem which I am trying to better understand. I am also including the author's preceding paragraph because I think it is useful:
Let us now consider the inhomogeneous $n$th order linear difference equations (1). If $\{v_j\}$ is any particular solution of equation (1) and $\{u_j\}$ is a solution of the homogeneous system (2), then $\{u_j + v_j\}$ is a solution of (1). This solution of (1) can be made to satisfy any particular initial conditions by adjusting the $\{u_j\}$. We now develop a discrete analog of what is known as Duhamel's principle in the theory of differential equations (where integral representations of solutions are obtained).
Theorem 3. Let $\{u_j^{(0)}\}_{j=0}^{\infty}, \ldots, \{u_j^{(n-1)}\}_{j=0}^{\infty}$ be the fundamental set of solutions of the $n$th order homogeneous difference equation (2) that satisfy the initial conditions
$$ \hspace{2.2cm} u_i^{(k)} = \delta_{ik}, \quad i = 0,1,\ldots,n-1; \quad k = 0,1,\ldots,n-1. \hspace{2.2cm} (7) $$Then the solution of (1) subject to the initial conditions
$$u_{0} = v_0, \; u_1 = v_1, \; \ldots, \; u_{n-1} = v_{n-1}$$
is given by
$$
\hspace{2cm} u_j = \sum_{k=0}^{n-1} v_k u_j^{(k)} + \frac{1}{a_n} \sum_{k=0}^{j-n} c_{k+n} u_{j-k-1}^{(n-1)}, \qquad j = 0, 1,2, \ldots \hspace{2cm} (8)
$$
where we define $u_i^{(n-1)} := 0$ for all $i < 0$ and $c_j := 0$ for all $j < n$.
My question. How would we obtain the second summation in (8)?
I have gone through the proof and am convinced the formula is correct, but the proof isn't particularly enlightening. I am looking for the intuition behind the second summation. This is my current understanding: By the discussion preceding the theorem, we can express the solution of the IVP as the sum of a particular solution of (1) and a solution of the homogeneous system. The sequence defined by the first summation clearly satisfies the initial conditions $u_i = v_i, \, i =0,1,\ldots,n-1$, and it is clearly a solution of the homogeneous problem (because such solutions are closed under linear combinations). The sequence defined by the second summation, which I'll denote by $\{w_j \}_{j=0}^{\infty}$, satisfies $w_j = 0$ for $j = 0,1,\ldots,n-1$, and (with some algebra) can be shown to satisfy the inhomogeneous equations (1). So the formula works….but how would we intuitively come up with such a sequence $\{w_j\}$?
Best Answer
The idea behind Duhamel's principle is to convert a nonhomogeneous equation with homogeneous initial conditions into a homogeneous equation with non homogeneous initial conditions. I explained the physical interpretation in a previous answer of mine How to understand the Duhamel's principle. In short, the inhomogeneous part of the equation is like an external force, that you can decompose in localised impulses by superposition. Each impulse can be seen as causing a discontinuity in the response. If the response has homogeneous initial conditions, it starts from zero. Thus, viewing the discontinuity in the right limit, this amounts to specifying an initial condition at the impulse time.
To solve (1), by superposition with a solution of (3), you reduce to the case of homogeneous initial conditions, this corresponds to the first sum in (8). Next, by superposition, you can assume that $c$ is supported exactly at one index, i.e. you write: $$ c = \sum_{i=n}^\infty c_i\delta_i $$ with: $$ \delta_i=(\delta_{ij})_{j\in\mathbb N} $$ the sequence that is 1 at $i$ and zero elsewhere. This is the analogue of decomposing the external force into elementary impulses. You then just need to solve (with $D$ the finite difference operator): $$ Dg^{(i)} = \delta_i $$ and you can reconstruct the solution as: $$ u = \sum_{i=n}^\infty c_ig^{(i)} $$ the sum is well defined because the $g^{(i)}$ are supported on $[i,\infty)$ (by "causality"), so each term of $u$ only involves the summation of a finite number of terms.
Now, to solve for the elementary solutions $g$ (Green's functions), the core of Duhamel's principle is to convert it into the homogeneous equation with suitable initial conditions. Applying (1) to $j=i-n$, the only index where the RHS is non zero, gives: $$ \sum_{s=0}^na_sg_{i-n+s}^{(i)} = 1 $$ only, for $k<i$ $g^{(i)}_k = 0$ since it is solution to the homogeneous equation with zero initial conditions. The equation thus becomes: $$ a_n g_i^{(i)} = 1 $$ i.e. $g_j^{(i)}=0$ for $j=i-n+1, ... i-1$ and $g_i^{(i)} = 1/a_n$. After the impulse, $g$ satisfies the homogeneous equation again. To summarize: $$ \begin{align} g_{k+i-n+1}^{(i)} &= \frac{\delta_{k,n-1}}{a_n} && k=0 … n-1\\ \sum_{s=0}^na_s g_{j+s+i-n+1}^{(i)}&=0 && j\in\mathbb N \end{align} $$ so you recognize the definition of $u^{(n-1)}$ for the subsequence $(g_{k+i-n+1}^{(i)})_{k\in\mathbb N}$. The key property is the shift invariance of the difference equation. The rest if the values are just zero so: $$ g_j^{(i)} = \begin{cases} \frac{1}{a_n}u^{(n-1)}_{j-i+n-1} & j > i-n \\ 0 & j < i \end{cases} $$ Note the overlap region for $i-n<j<i$ where both give zero. You now have all you need to reconstruct the solution.
Explicitly: $$ \begin{align} u_j &= \sum_{k=0}^{n-1}v_ku^{(k)}_j+\sum_{i=n}^\infty c_ig_j^{(i)} \\ &= \sum_{k=0}^{n-1}v_ku^{(k)}_j+\sum_{i=n}^jc_ig_j^{(i)} && \text{(causality)}\\ &= \sum_{k=0}^{n-1}v_ku^{(k)}_j+\sum_{i=n}^jc_i\frac{1}{a_n}u^{(n-1)}_{j-i+n-1}\\ &= \sum_{k=0}^{n-1}v_ku^{(k)}_j+\frac{1}{a_n}\sum_{k=0}^{j-n}c_{k+n}u^{(n-1)}_{j-k-1} && k = i+n \end{align} $$
Hope this helps.