Peeter Joot's (OLD) Blog.

Math, physics, perl, and programming obscurity.

Posts Tagged ‘Ill conditioned LDEs’

PHY454H1S Continuum Mechanics. Lecture 20: Asymptotic solutions of ill conditioned equations. Taught by Prof. K. Das.

Posted by peeterjoot on March 30, 2012

[Click here for a PDF of this post with nicer formatting and figures if the post had any (especially if my latex to wordpress script has left FORMULA DOES NOT PARSE errors.)]


Peeter’s lecture notes from class. May not be entirely coherent.

Two ill conditioned LDEs.

We’ll consider two cases, both ones that we can solve exactly

  • With u(0) = 1, and letting \epsilon \rightarrow 0, we’ll look at solutions of the ill conditioned LDE

    \begin{aligned}\epsilon \frac{du}{dy} + u = y\end{aligned} \hspace{\stretch{1}}(2.1)

  • With u(0) = 0, u(1) = 2, and 0 < \epsilon \ll 1 we’ll look at the second order ill conditioned LDE

    \begin{aligned}\epsilon \frac{d^2u}{dy^2} + \frac{du}{dy} = 1\end{aligned} \hspace{\stretch{1}}(2.2)

The first order LDE.

Exact solution.

Our homogeneous equation is

\begin{aligned}\epsilon \frac{du}{dy} + u = 0,\end{aligned} \hspace{\stretch{1}}(2.3)

with solution

\begin{aligned}u \propto e^{-y/\epsilon}.\end{aligned} \hspace{\stretch{1}}(2.4)

Looking for a solution of the form

\begin{aligned}u = A(y) e^{-y/\epsilon},\end{aligned} \hspace{\stretch{1}}(2.5)

we find

\begin{aligned}\epsilon A' e^{-y/\epsilon} = y,\end{aligned} \hspace{\stretch{1}}(2.6)

and integrate to find

\begin{aligned}A(y) = (y - \epsilon) e^{x/\epsilon} + C.\end{aligned} \hspace{\stretch{1}}(2.7)

Application of the boundary value constraints give us

\begin{aligned}u = ( 1 + \epsilon ) e^{-y/\epsilon} + y - \epsilon.\end{aligned} \hspace{\stretch{1}}(2.8)

This is plotted in figure (1).

Figure 1: Plot of exact solution to simple first order ill conditioned LDE.

Limiting cases.

We want to consider the limiting case where

\begin{aligned}0 < \epsilon \ll 1,\end{aligned} \hspace{\stretch{1}}(2.9)

and we let \epsilon \rightarrow 0. If y = O(1), then we have

\begin{aligned}u \approx 1 \times 0 + y,\end{aligned} \hspace{\stretch{1}}(2.10)

or just

\begin{aligned}u = y.\end{aligned} \hspace{\stretch{1}}(2.11)

However, if y = O(\epsilon) then we have to be more careful constructing an approximation. When y is very small, but \epsilon is also of the same order of smallness we have

\begin{aligned}e^{-y/\epsilon} \ne 0.\end{aligned} \hspace{\stretch{1}}(2.12)

If \epsilon \rightarrow 0 and y \rightarrow O(\epsilon)

\begin{aligned}e^{-y/\epsilon} \rightarrow e^{-\epsilon O(1) /\epsilon} \rightarrow e^{-O(1)}\end{aligned} \hspace{\stretch{1}}(2.13)


\begin{aligned}u \approx e^{-y/\epsilon}\end{aligned} \hspace{\stretch{1}}(2.14)

Approximate solution in the inner region.

When y = O(\epsilon) define a new scale

\begin{aligned}Y = \frac{y}{\epsilon}\end{aligned} \hspace{\stretch{1}}(2.15)

\begin{aligned}y = O(1)\end{aligned} \hspace{\stretch{1}}(2.16)

so that our LDE takes the form

\begin{aligned}\frac{du}{dY} + u = \epsilon Y\end{aligned} \hspace{\stretch{1}}(2.17)

When \epsilon \rightarrow 0 we have

\begin{aligned}\frac{du}{dY} + u \approx 0.\end{aligned} \hspace{\stretch{1}}(2.18)

We have solution

\begin{aligned}\ln u = -Y + \ln C,\end{aligned} \hspace{\stretch{1}}(2.19)


\begin{aligned}u \propto e^{-Y} = e^{-y/\epsilon}.\end{aligned} \hspace{\stretch{1}}(2.20)

Question: Couldn’t we just Laplace transform.
Answer given: We’d still get into trouble when we take \epsilon \rightarrow 0. My comment: I don’t think that’s strictly true. In an example like this where we have an exact solution, a Laplace transform technique should also yield that solution. I think the real trouble will come when we attempt to incorporate the non-linear inertial terms of the Navier-Stokes equation.

Second order example.

Exact solution.

We saw above in the first order system that our specific solution was polynomial. While that was found by the method of variation of parameters, it seems obvious in retrospect. Let’s start by looking for such a solution, starting with a first order polynomial

\begin{aligned}u = A y + B.\end{aligned} \hspace{\stretch{1}}(2.21)

Application of our LDE operator on this produces

\begin{aligned}A &= 1 \\ B &= 0.\end{aligned}

Now let’s move on to find a solution to the homogeneous equation

\begin{aligned}\epsilon \frac{d^2u}{dy^2} + \frac{du}{dy} = 0\end{aligned} \hspace{\stretch{1}}(2.22)

As usual, we look for the characteristic equation by assuming a solution of the form u = e^{m y}. This gives us

\begin{aligned}\epsilon m^2 + m = (\epsilon m + 1) m = 0,\end{aligned} \hspace{\stretch{1}}(2.23)

with roots

\begin{aligned}m = 0, -1/\epsilon.\end{aligned} \hspace{\stretch{1}}(2.24)

So our homogeneous equation has the form

\begin{aligned}u(y) = A e^{-y/\epsilon} + B,\end{aligned} \hspace{\stretch{1}}(2.25)

and our full solution is

\begin{aligned}u(y) = A e^{-y/\epsilon} + B + y\end{aligned} \hspace{\stretch{1}}(2.26)

with the constants A and B to be determined from our boundary value conditions. We find

\begin{aligned}0 &= u(0) = A + B + 0 \\ 2 &= u(1) = A e^{-1/\epsilon} + B + 1.\end{aligned}

We’ve got B = -A and by subtracting

\begin{aligned}A ( e^{-1/\epsilon} -1 ) = 1.\end{aligned} \hspace{\stretch{1}}(2.27)

So the exact solution is

\begin{aligned}u = y + \frac{e^{-y/\epsilon} - 1}{e^{1/\epsilon} - 1}.\end{aligned} \hspace{\stretch{1}}(2.28)

This is plotted in figure (2).

Figure 2: Plot of our ill conditioned second order LDE

Solution in the regular region.

For small \epsilon relative to y our LDE is approximately

\begin{aligned}\frac{du}{dy} = 1,\end{aligned} \hspace{\stretch{1}}(2.29)

which has solution

\begin{aligned}u = y + C\end{aligned} \hspace{\stretch{1}}(2.30)

Our u(1) = 2 boundary value constraint gives us

\begin{aligned}C = 1.\end{aligned} \hspace{\stretch{1}}(2.31)

Our solution in the regular region where \epsilon \rightarrow 0 and y = O(1) is therefore just

\begin{aligned}u = y + 1.\end{aligned} \hspace{\stretch{1}}(2.32)

Solution in the ill conditioned region.

Now let’s consider the inner region. We’ll see below that when y = O(\epsilon), and we allow both \epsilon and y tend to zero independently, we have approximately

\begin{aligned}u \sim 1 - e^{-y/\epsilon}.\end{aligned} \hspace{\stretch{1}}(2.33)

We’ll now show this. We start with a helpful change of variables as we did in the first order case

\begin{aligned}Y = \frac{y}{\epsilon}.\end{aligned} \hspace{\stretch{1}}(2.34)

When y = O(\epsilon) and Y = O(1) we have

\begin{aligned}\not{{\epsilon}} \frac{d^2 u}{\not{{\epsilon}} \epsilon dY^2} + \frac{1}{{\epsilon}} \frac{du}{dY} = 1\end{aligned} \hspace{\stretch{1}}(2.35)


\begin{aligned}\frac{d^2 u}{dY^2} + \frac{du}{dY} = \epsilon.\end{aligned} \hspace{\stretch{1}}(2.36)

This puts the LDE into a non ill conditioned form, and allows us to let \epsilon \rightarrow 0. We have approximately

\begin{aligned}\frac{d^2 u}{dY^2} + \frac{du}{dY} = 0.\end{aligned} \hspace{\stretch{1}}(2.37)

We’ve solved this in our exact solution work above (in a slightly more general form), and thus in this case we have just

\begin{aligned}u = A + B e^{-Y}\end{aligned} \hspace{\stretch{1}}(2.38)

at Y = 0 we have

\begin{aligned}u = A + B = 0.\end{aligned} \hspace{\stretch{1}}(2.39)

so that

\begin{aligned}B = -A.\end{aligned} \hspace{\stretch{1}}(2.40)

and we find for the inner region

\begin{aligned}u = A (1 - e^{-Y}) = A( 1 - e^{-y/\epsilon} ) \sim 1 - e^{-y/\epsilon}.\end{aligned} \hspace{\stretch{1}}(2.41)

Taking these independent solutions for the inner and outer regions and putting them together into a coherent form (called matched asymptotic expansion) is a rich and tricky field. For info on that we’ve been referred to [1].


[1] EJ Hinch. Perturbation methods, volume 6. Cambridge Univ Pr, 1991.


Posted in Math and Physics Learning. | Tagged: , , , , | Leave a Comment »