ODE boundary value problem relaxation method reference request

Question

This is a somewhat basic question, I guess. Take the ODE boundary value problem $$ \frac1\lambda y''-y'=0, \qquad y(0) = 0, \quad y(1) = 1, $$ with the solution $$ y(x) = \frac{e^{x\lambda}-1}{e^\lambda -1}. $$ Discretize it on a uniform grid $x_k = kh$, $h=1/n$, $k=0,\ldots,n$, as follows: $$ \frac1\lambda \frac{y_{k-1}-2y_k+y_{k+1}}{h^2}-\frac{y_{k+1}-y_{k-1}}{2h} = 0, $$ and solve the resulting tridiagonal system, for, e.g., $n=100, \lambda=10^4$:

With increasing $n$, the scheme nevertheless converges.

a. Can you suggest a good reference that describes this kind of thing? I think that stiffness of initial-value problems is not quite the same thing, which is why I'm asking this question, but maybe I'm wrong.

b. Short of choosing a better method, is there a "filter" I can apply to the output $\mathbf{y}_{0:n}$ to recover a better approximation of the true solution (in particular, make the transformed output positive, as it should be).

Answer 1

This is a convection dominated problem. Since $1/\lambda$ is small, your equation is approximately $y'=0$. But the solution of this cannot satisfy your boundary conditions. The solution will rapidly change from 0 to 1 as you approach $x=1$, and a boundary layer is created at $x=1$. The term $y''/\lambda$ cannot be neglected in the boundary layer.

If you write

$$ y_k = \frac{1}{2}(1 + \lambda h/2) y_{k-1} + \frac{1}{2}(1 - \lambda h/2) y_{k+1} $$

the solution will be monotone if

$$ \lambda h < 2 $$

which is called the cell Peclet condition. So your mesh must be finer than

$$ h < \frac{2}{\lambda} $$

You can see that increasing $\lambda$ requires very fine mesh. If your mesh satisfies above condition you should get a good solution.

The other alternative is to use an upwind scheme. In your case it is (please check stability as above)

$$ \frac{1}{\lambda} \frac{y_{k-1} - 2 y_k + y_{k+1}}{h^2} - \frac{y_k - y_{k-1}}{h} = 0 $$

but this is only first order accurate.

Another way to analyze this is in terms of an M-matrix.

For more, you can refer to some texts, e.g.,

1. Strikwerda: Finite Difference Schemes and Partial Differential Equations
2. Wesseling: Principles of Computational fluid dynamics