%% 20. Boundary and interior layers % % \setcounter{page}{253} % %% % % We have seen examples of boundary and interior layers % already in Chapters~4, 7, and~18. This is a topic where % the use of a tool like Chebfun is particularly compelling. % With a computation or two, one quickly sees where % layers lie and how their thickness depends on parameters. % %% % % To start the discussion, here is % a reprise of the example of equation (5.3), % $$y'' = y , \quad x \in [\kern .3pt 0,40\kern .3pt], % ~y(0) = y(40) = 1. \eqno (20.1)$$ % ODEformats L = chebop(0,40); L.op = @(x,y) diff(y,2) - y; L.lbc = 1; L.rbc = 1; y = L\0; plot(y,CO,bvp), axis([0 40 -1.5 1.5]) title('Fig.~20.1.~~Simplest example of boundary layers (20.1)',FS,11) %% % \vskip 1.02em %% % % \noindent % Just by looking at the image, one sees immediately at some % level what is going on: the solution wants to be zero'' and % more or less achieves this aim, apart from regions % near the boundaries, % where it goes through rapid transitions to satisfy the % boundary conditions. What does it mean that it wants to be % zero''? This ODE is easy to analyze: the general solution % is $A\exp(-x) + B \exp(x)$. This implies % that all components grow exponentially % either as $x$ increases or as $x$ decreases. Therefore, the % only way $y$ can deviate much from zero without deviating % {\em hugely\/} from zero is for this to happen near a boundary. % %% % % Problems with boundary or interior layers usually contain % a large or a small parameter. In (20.1), the parameter is the length % of the interval, 40. To move to a more standardized % formulation, let us % transplant the equation to the interval $[-1,1]$ % with a small constant $\varepsilon$ multiplying % the highest-order derivative: % $$\varepsilon \kern .3pt y'' = y , \quad x \in [-1,1], % ~y(-1) = y(1) = 1. \eqno (20.2)$$ % The rescaling multiplies % $y''$ by $(20)^2$, so the right choice if we wish % to match (20.1) is $\kern .5pt \varepsilon = (20)^{-2}$. % As mentioned in Chapter~7, a problem like this with a small % parameter multiplying the highest derivative is called a % {\bf singular perturbation problem.} % Here is the solution to (20.2), looking the % same as for (20.1). % L = chebop(-1,1); ep = (1/20)^2; L.op = @(x,y) ep*diff(y,2) - y; L.lbc = 1; L.rbc = 1; y = L\0; plot(y,CO,bvp), axis([-1 1 -1.5 1.5]) title('Fig.~20.2.~~Same problem rescaled to $[-1,1]$, (20.2)',FS,11) %% % \vskip 1.02em %% %% % This example is so simple that we can quickly work out % exact formulas. % The general solution to the rescaled problem (20.2) is % $A\exp(-\varepsilon^{-1/2}x) + % B \exp(\varepsilon^{-1/2}x)$, or % as we may rewrite it with redefinitions of $A$ and $B$, % $$y(x) = A\exp(-\varepsilon^{-1/2}(1+x)) % + B \exp(-\varepsilon^{-1/2}(1-x)). \eqno (20.3)$$ % The boundary conditions correspond to the equations % $$A + B\exp(-2\kern .5pt\varepsilon^{-1/2}) = % B + A\exp(-2\kern .5pt\varepsilon^{-1/2}) = 1,$$ % with solution % $$A = B = {1-\exp(-2\kern .5pt\varepsilon^{-1/2})\over % 1-\exp(-4\kern .5pt\varepsilon^{-1/2})} = % 1 + O(\exp(-2\kern .5pt\varepsilon^{-1/2})).$$ % Therefore (20.3) can be written % $$y(x) = \exp(-\varepsilon^{-1/2}(1+x)) + % \exp(-\varepsilon^{-1/2}(1-x)) + O(\exp(-3\kern.3pt\varepsilon^{-1/2})). % \eqno (20.4)$$ % At $x=0$, for example, we have % $$y(0) = 2 \exp(-\varepsilon^{-1/2}) % + O(\exp(-3\kern .5pt \varepsilon^{-1/2})),$$ % which for $\varepsilon = (20)^{-2}$ becomes % $$y(0) = 2\exp(-20) + O(\exp(-60)).$$ % We confirm this result numerically: y(0), 2*exp(-20) %% % % Although this example is very simple, it illustrates a general idea. % Consider the approximation (20.4), a sum of two terms. % The first term consists of an approximation % at the left boundary, and the second % at the right. Together, they provide an % accurate approximation to the solution throughout $[-1,1]$. % Such combinations are the business % of {\bf boundary layer analysis,} which constructs % approximate solutions to all kinds of problems % by combining different approximations in different % regions. This is an advanced art, % applicable for nonlinear as well as linear problems. % In the remainder of this chapter, we will take some % introductory steps into this subject for linear equations. % %% % % The starting point of boundary layer analysis is this principle: % \begin{enumerate} % \item {\em Outside a boundary or interior % layer, terms involving\/ $\varepsilon$ are negligible.} % \end{enumerate} % For example, for the equation (20.2), outside any boundary or interior layers, % we can expect % $$y = 0$$ % to hold approximately. This equation applicable outside % layers is called the {\bf outer equation}. % %% % % If (20.2) is modified to % include a further term, the outer equation will adjust accordingly. For % example, suppose we consider % $$\varepsilon \kern .3pt y'' = y + 2x , \quad x \in [-1,1], % ~y(-1) = y(1) = 1. \eqno (20.5)$$ % Now the outer equation is % $$y = -2x,$$ % and the solution changes shape to look like $y = -2x$ away % from $x = \pm 1$. % L.op = @(x,y) ep*diff(y,2) - y - 2*x; y = L\0; plot(y,CO,bvp), axis([-1 1 -2 2]) title('Fig.~20.3.~~Modified problem (20.5)',FS,11) %% % \vskip 1.02em %% %% % % To proceed further with boundary layer analysis we add % a pair of principles applicable inside a transition layer % (i.e., a boundary or interior layer): % \begin{enumerate} % \setcounter{enumi}{1} % \item % {\em Inside a layer, a low order derivative % is negligible compared with a higher order one.} % \item % {\em Inside a layer, a variable coefficient can be % approximated by a constant, at least if it is locally nonzero.} % \end{enumerate} % Let us see how these ideas work out in a % slightly more complicated example, equation (7.8), % $$\varepsilon \kern .3pt y'' + x\kern .1pt y' + x\kern .1pt y = 0, % \quad x\in [-2,2\kern .2pt ], ~ y(-2) = -4, ~ % y(2) = 2. \eqno (20.6)$$ % With $\varepsilon = 0.001$, the % solution shows an interior layer at $x=0$, as we saw already % in Fig.~7.6. % L = chebop(-2,2); L.op = @(x,y) .001*diff(y,2) + x*diff(y) + x*y; L.lbc = -4; L.rbc = 2; y = L\0; plot(y,CO,bvp), ylim([-6 18]) title('Fig.~20.4.~~Interior layer (20.6)',FS,11), hold on %% % \vskip 1.02em %% % % \noindent % For a boundary layer analysis of this % result, applying principle (1), we replace % $\varepsilon$ by zero in (20.6) and get the outer equation % $$x\kern .1pt y' + x\kern .1pt y = 0 . \eqno (20.7)$$ % Thus we immediately have a prediction as to the behavior % of a solution to (20.6). In every region, % (i) $y''$ will % be very large (of order $\varepsilon^{-1}$), % or (ii) $x$ will be close to % zero, or (iii) $y(x)$ will approximate $Ce^{-x}$ for some $C$. % %% % From the picture it is clear that case (iii) applies for % $|x|\gg 0$, with outer solutions $-4\exp(-(x+2))$ on the left % and $2\exp(-(x-2))$ on the right. Here are thick dashed lines % marking these approximations: fleft = chebfun('-4*exp(-(x+2))',[-2,.2]); h1 = plot(fleft,'--r',LW,2); fright = chebfun('2*exp(-(x-2))',[-.1,2]); h2 = plot(fright,'--r',LW,2); ylim([-6 18]) title('Fig.~20.5.~~Outer solutions on left and right',FS,11) %% % \vskip 1.02em %% % % \noindent % To explain the rest of the solution, we % apply principle (2). In a transition layer, % we can expect $y''$ to be much larger than $y'$ and $y'$ to be % much larger than $y$. An approximation of (20.6) in a layer % is accordingly % $$\varepsilon \kern .3pt y'' + x\kern .1pt y' = 0. \eqno (20.8)$$ % First, let us use this equation to see why there is no % boundary layer at $x=2$. Applying % principle (3), we approximate (20.8) further by % $$\varepsilon \kern .3pt y'' + 2\kern .3pt y' = 0. \eqno (20.9)$$ % The general solution to this equation is % $A + B\exp(-2\kern .3pt \varepsilon^{-1}x)$, and the % minus sign in the second of these terms implies that % it decreases rather than increases with increasing $x$. % Similarly, at $x=-2,$ % $$\varepsilon \kern .3pt y'' - 2\kern .3pt y' = 0, \eqno (20.10)$$ % with general solution $A + B\exp(2\kern .3pt \varepsilon^{-1}x)$, % and again the sign is such that no boundary layer is possible. % %% % % Only $x=0$ remains a candidate for a transition layer, matching % the shape of the solution in the plot. % Let us accordingly consider (20.8) for $x\approx 0$. Here we % rescale variables to see what is happening. We define % $$s = x\kern .5pt \varepsilon^{-1/2}, \quad u(s) = y(x),$$ % which leads to % $$y''(x) = \varepsilon^{-1} u''(s), \quad % y'(x) = \varepsilon^{-1/2} u'(s),$$ % whereupon (20.8) becomes % $$u'' + s u' = 0, \quad u(-\infty) = % -4\kern .3pt e^{-2}, ~ u(\infty) = 2\kern .3pt e^2 . \eqno (20.11)$$ % One solution of the ODE part of this equation % is $u(s) = A$ for any constant~$A$, and the other can % be obtained by setting $w=u'$, yielding the equation % $w' + sw = 0$. Separation of variables gives % $w(s) = B\exp(-s^2/2)$ for a constant $B$, or after % integration and redefinition of $B$, % $u(s) = B \hbox{\kern 1pt erf} (s/\sqrt 2)$, where % $\kern .5pt\hbox{erf}\kern .5pt$ is the error function. Thus we have % $$u(s) = A + B\hbox{\kern 1pt erf}\left({s\over \sqrt 2 \kern 1pt}\right ) , \quad % A = e^2 - 2e^{-2}, ~ B = e^2+2e^{-2}.$$ % Here we superimpose on the plot a new dashed line corresponding to % this approximation. % set(h1,CO,[1 .7 .7]), set(h2,CO,[1 .7 .7]) A = exp(2)-2*exp(-2); B = exp(2)+2*exp(-2); ep = .001; s = @(x) x/sqrt(ep); fmiddle = chebfun(@(x) A + B*erf(s(x)/sqrt(2)),[-.8 .8]); plot(fmiddle,'--r',LW,2), hold off title('Fig.~20.6.~~Approximation to the transition layer',FS,11) %% % \vskip .2em %% % % This completes our study of (20.6). Now let us look % at the same equation, except with a sign change on the % second derivative term: % $$-\varepsilon \kern .3pt y'' + x\kern .1pt y' + x\kern .1pt y = 0, % \quad x\in [-2,2\kern .2pt ], ~ y(-2) = -4, ~ % y(2) = 2. \eqno (20.12)$$ % The solution looks completely different. % Instead of an interior layer, we % have boundary layers at each end, but they are so thin % as to be almost invisible in the plot. % L.op = @(x,y) -.001*diff(y,2) + x*diff(y) + x*y; y = L\0; plot(y,CO,bvp), ylim([-5,5]) title('Fig.~20.7.~~Sign change (20.12)',FS,11) %% % \vskip 1.02em %% % % \noindent % An experienced eye immediately % conjectures that the boundary layers are probably of width % $O(\varepsilon)$ rather than $O(\varepsilon^{1/2})$, as a % closeup at the left boundary confirms. % plot(y{-2,-1.98},CO,bvp), ylim([-5,5]) title('Fig.~20.8.~~Closeup near left boundary',FS,11) %% % \vskip 1.02em %% % The scaling of the widths of boundary and interior layers % is of fundamental importance in applications. % Let us illustrate this matter further with a constant-coefficient % equation, an advection-diffusion problem, % $$\varepsilon \kern .3pt y'' + y' + y = 0 , \quad % x \in [\kern .3pt 0,1],~ y(0)=y(1) = -1. \eqno (20.13)$$ % Here are the solutions for $\varepsilon = 0.1,$ $0.02,$ % and $0.004$. L = @(ep) chebop(@(x,y) ep*diff(y,2)+diff(y)+y,[0,1],'dirichlet'); y = @(ep) L(ep)\(-1); for ep = [.1 .02 .004] plot(y(ep),CO,bvp), ylim([0 2]), hold on end title(['Fig.~20.9.~~Solutions to (20.13) with ' ... '$\varepsilon = 0.1, 0.02, 0.004$'],FS,11), hold off %% % \vskip 1.02em %% % % \noindent % Let us informally define the width of a boundary layer to be % the distance from~$0$ to the value of~$x$ at which % the solution first reaches $y(x)=0.5$ (compare Exercise~18.2). % The widths for four values of $\varepsilon$ % show convincing $O(\varepsilon)$ dependence. % width = @(ep) min(roots(y(ep)-0.5)); disp(' epsilon boundary layer width') for ep = 10.^(-1:-1:-4) fprintf('%10.5f %10.7f\n',ep,width(ep)) end %% % % Here is another variation on the theme of (20.6) and (20.12), % with the left-hand boundary condition made positive to enhance the % effect: % $$0.001y'' + x\kern .1pt y' - y = 0, % \quad x\in [-2,2\kern .2pt], ~ y(-2) = 4, ~ % y(2) = 2. \eqno (20.14)$$ % The transition here is called a {\em corner layer}. % L = chebop(-2,2); L.lbc = 4; L.rbc = 2; L.op = @(x,y) .001*diff(y,2) + x*diff(y) - y; y = L\0; plot(y,CO,bvp), ylim([-1 5]) title('Fig.~20.10.~~Corner layer (20.14)',FS,11) %% % \vskip 1.02em %% % % \noindent % Almost anything is possible with interior and boundary % layers, and we give this example as just another illustration. % %% % % \begin{center} % \hrulefill\\[1pt] % {\sc Application: why is New York hotter than San Francisco?} \\[-3pt] % \hrulefill % \end{center} % %% % % New York is hot in July, around 9 degrees hotter than % San Francisco ($25^\circ$C vs.\ $16^\circ$ are typical % 24-hour averages). % Yet in January, New York is 9 degrees colder than % San Francisco ($1^\circ$ vs.\ $10^\circ$). % What's going on? Why is the weather so much more moderate % on the west coast of the USA than the east coast? % We shall explain the effect in words first, % then explore the mathematics and the connection with boundary layers. % %% % The first key fact is that in the middle latitudes of Earth, prevailing % winds blow from west to east. % This means that San Francisco's weather comes from % over the ocean, whereas New York's comes from over the land. %% % % The second key fact is that ocean is liquid, and land is solid. % When the summer sun beats down on the land, it heats up, but only % down to a depth of a few meters, because the heat transfer is limited % by conduction, which is slow. This lets the surface % get very hot. Water, on the other hand, transfers heat % by the much faster mechanism of convection --- fluid motion. % The stirring isn't enough to bring % the summer heat down to the great depths, but it certainly brings it deeper % than a few meters. The heat is spread over a much greater % volume, so that the surface of the ocean never gets % anywhere near as hot as the surface of the land, and San Franciscans % can get by without air conditioning. % %% % % To model this effect mathematically, we shall be very rough, % since the point is just to understand the main mechanism. % Let us imagine that the earth is 100 meters deep, with a temperature % $T$ fixed at $13^\circ$C at the bottom (the average of the figures % reported above for both New York and San Francisco). % Our depth variable will be $x$, in units of meters, % running from $-100$ to $0$. % At the earth's surface, let us imagine that heat flows in at % the rate of 100 watts per square meter in July and flows out % at the same rate in January, with a sinusoidal dependence on $t$ in-between, % $$\hbox{incoming heat in } W/m^2 = 100 \exp(2\pi it) ,$$ % where $t$ is measured in years. % (We use a complex exponential for the usual reason of convenience % explained in Chapter~4; the physical temperature is the real part.) % To get to an ODE, % we start with a PDE, the heat equation, % $$\rho \kern .5pt c \kern .5pt {\partial T\over \partial t} = % \left({\hbox{year}\over\hbox{sec}}\right) k % \kern .7pt {\partial^2 T\over \partial x^2}, \eqno (20.15)$$ % where $\rho$ denotes density, $c$ is heat capacity, and $k$ is % conductivity. % (PDE\kern .5pt s are considered systematically in Chapter 22.) % We measure $\rho$ in kilograms per cubic meter, % $c$ in Joules per kilogram per degree, and $k$ in watts per meter per degree. % The quotient $(\hbox{year}/\hbox{sec})$, which is equal to % about $3.1\times 10^7$, is included to reconcile our % choice of units of years for $t$ and watts (Joules per second) % in the definition of $k$. % With the not too unreasonable choices % $$\rho = 3000, \quad c = 2000 ,$$ % (20.15) reduces to % $${\partial T\over \partial t} = % 5.3 \kern .5pt k \kern 1pt {\partial^2 T\over \partial x^2}. % \eqno (20.16)$$ % Meanwhile the boundary condition at $x = -100$ is % $$T(-100,t) = 13 ,$$ % and our assumption concerning % the heat influx gives the boundary condition at $x=0$, % $${\partial T\over \partial x}(0,t) = 100 \kern .5pt % k^{-1}\exp(2\pi it).$$ % We now reduce the problem to an ODE by separating variables, assuming % the solution takes the form % $$T(x,t) = 13 + e^{2\pi it} (y(x)-13)$$ % for some function $y$. The PDE (20.16) becomes % $$2\pi i (y-13) = 5.3 \kern .3pt k \kern .5pt y'', \eqno (20.17)$$ % and the boundary conditions become % $$y(-100) = 13, \quad y'(0) = 100\kern .5pt k ^{-1}. \eqno (20.18)$$ % %% % % We are ready to compute a solution. Taking $k=3$ as a % reasonable value of the conductivity, here is the % temperature distribution we get from (20.17)--(20.18). % L = chebop(-100,0); k = 3; L.op = @(y) 5.3*k*diff(y,2) -2i*pi*(y-13); L.lbc = 13; L.rbc = @(y) diff(y)-100/k; y = L\0; plot(real(y),CO,bvp), grid on, ylim([0 60]) title('Fig.~20.11.~~July temperature in the earth',FS,11) xlabel('depth $x$',FS,10), ylabel('temperature (${}^\circ$C)',FS,9) %% % \vskip 1.02em %% % % \noindent % Note how hot it is in the boundary layer at the surface! New York's summers % begin to make sense. Note also the dip below $13^\circ$C in the % temperature at a depth of about 3 meters. This effect is genuine, % reflecting a bit of wavelike behavior in the conduction of % heat.\footnote{The thick stone walls of Italian % villas exploit this effect to keep cool by day and warm by night. % A thickness of $3\hbox{m}/20 = 15$cm ought to be of the right % order of magnitude, since % days are about $20^2$ times shorter than years and the % effect in question involves a square root.} % %% % As for mild San % Francisco, we are not about to go into details of the physics of % convection in the upper ocean. To get the general idea, however, let's % imagine that ocean is like land, but with 100 times greater % conductivity. The heat now reaches much deeper, % giving a much more moderate surface temperature. k = 300; L.op = @(y) 5.3*k*diff(y,2) -2i*pi*(y-13); L.rbc = @(y) diff(y)-100/k; y = L\0; plot(real(y),CO,bvp), grid on, ylim([0 60]) title('Fig.~20.12.~~July temperature in the ocean',FS,11) xlabel('depth $x$',FS,10), ylabel('temperature (${}^\circ$C)',FS,9) %% % \vskip 1.02em %% % % The details of our models have hardly been precise, % but the main message of % Figs.~20.11 and~20.12 is genuine. New York is hotter than % San Francisco in the summer because its weather is controlled by % a thinner boundary layer. % %% % % \smallskip % {\sc History.} % Boundary-layer analysis began with an epochal % contribution by the German fluid mechanician % Ludwig Prandtl. An important feature of fluid mechanics is viscosity, % that is, friction, but viscous effects are often insignificant % away from solid boundaries. % Prandtl realized that one could exploit this effect mathematically % by dropping certain terms from the Navier--Stokes partial differential % equations of fluid mechanics near a wall and dropping % other terms away from the wall. These simplifications led to % major advances in engineering, making % it possible to solve problems that had previously been % inaccessible. Prandtl presented these methods % in 1904, and ever since then, % boundary-layer analysis has been recognized as one of the % fundamental tools of applied mathematics. % %% % % \smallskip % {\sc Our favorite reference.} Generations of % applied mathematics graduate % students at MIT have learned asymptotics from % the exceptionally rich book by Bender and % Orszag, {\em Advanced Mathematical Methods for Scientists % and Engineers,} McGraw-Hill, 1978. % \smallskip % %% % % \begin{displaymath} % \framebox[4.7in][c]{\parbox{4.5in}{\vspace{2pt}\sl % {\sc Summary of Chapter 20.} Singular perturbation problems, where % a small parameter\/ $\varepsilon$ % multiplies the highest derivative, typically have % solutions featuring rapid transitions known as boundary % or interior layers. The development of approximate solutions % for such problems starts from three principles. % (1) Outside a layer, terms involving\/ $\varepsilon$ are negligible. % (2) Inside a layer, a low order derivative % is negligible compared with a higher order one. % (3) Inside a layer, a variable coefficient can be % approximated by a constant, at least if it is locally nonzero. % \vspace{2pt}}} % \end{displaymath} % %% % \smallskip\small\parskip=1pt\parindent=0pt % \newpage % {\em \underline{Exercise $20.1$}. Corner layer.} % {\em (a)} Plot the solution of $\varepsilon \kern .1pt y'' % + x\kern .1pt y' - y = 0$, $y(-1) = 1$, $y(1) = 1$, $\varepsilon = 0.001$ % (compare Fig.~20.10). % Which two terms are balanced to determine the % outer solution, outside the transition region $x\approx 0\kern .5pt$? % Write down the outer equation and solve it to explain % the behavior of $y$ in this region. % {\em (\kern .7pt b)} Which two terms are balanced for $x\approx 0\kern .5pt$? % Write down the inner equation and show how it can capture the shape % of $y(x)$ in this region. % {\em (c)} % Plot $y$ again, but now on a semilogy scale and zooming % in with {\tt axis([-.1 .1 1e-4 1])}. % On the same plot, superimpose curves corresponding to the other % two terms of the ODE, $\varepsilon\kern.1pt y''$ and $x\kern .1pt y'$. % For what range of values of $x$ are the two terms that are balanced for % the outer solution at least 10 times bigger than the omitted term? % What is the answer to the same question for the inner solution? % \par % {\em \underline{Exercise $20.2$}. Cusp.} % {\em (a)} If we change the equation of the previous exercise by just % one coefficient, $\varepsilon \kern .1pt y'' % + x\kern .1pt y' - y/2 = 0$, the solution changes in an interesting way. % Plot the solution for $\varepsilon = 10^{-3}$. % Find the outer equation and explain the shape of the solution. % {\em (\kern .7pt b)} Use {\tt tic} and {\tt toc} to measure the % computing time for this problem with $\varepsilon = 10^{-3}$, % $10^{-4}$, and $10^{-5}$. The reason for the slow-down is that % Chebfun is utilizing a Chebyshev grid of thousands of points, requiring % the solution of a matrix problem of dimension in the thousands. % Change the chebop domain to {\tt L.domain = [-1 0 1]} % and repeat the same timings; also plot the solution in the case % $\varepsilon = 10^{-6}$. Now Chebfun has introduced a breakpoint % at $x=0$, so that separate Chebyshev grids are utilized on either side, % and this makes the matrices smaller. (Although the user can % introduce breakpoints like this at fixed locations, Chebfun % does not offer adaptive gridding capabilities for ODE BVPs.) % \par % {\em \underline{Exercise $20.3$}. Variable oscillations.} % Solve $\varepsilon \kern .3pt y'' + (x^2-1) y = 0$, $x\in[-2,2\kern .2pt]$ % with $y(-2) = 1$ and $y(2) = 2$ for % $\varepsilon =10^{-5}$ and plot $y$. % Explain the nature of this solution and derive an approximate % formula for the wavelength of the oscillations in the oscillatory % regions. (Compare Fig.~7.5.) % Note that for this problem, principle (1) of our general method of % boundary layer analysis does not apply because if the % term involving $\varepsilon$ is dropped, there is just one % other term remaining, with nothing to balance it. % \par % {\em \underline{Exercise $20.4$}. A double boundary layer.} % {\em (a)} Solve $\varepsilon^3 \kern .1pt y'' + x^3y' + % (x^3-\varepsilon) y = 0$, $x\in[\kern .3pt 0,1]$ % with $y(0) = 2$ and $y(1) = 1$ for % $\varepsilon =0.005$ and plot the resulting solution. % This problem has two boundary layers at $x=0$, one of width % $O(\varepsilon)$ and the other of width $O(\varepsilon^{1/2})$. % Find a way to show this graphically based on solutions computed % for various values of~$\varepsilon$. % {\em (\kern .7pt b)} A good global approximation for $y$ is % $y(x) \approx % 2\kern .2pt \exp(-x/\varepsilon) + e\kern .3pt % [\kern .3pt \exp(-\varepsilon/2\kern .2pt x^2) % + \exp(-x) - 1\kern .3pt]$. What is the maximum error % of this approximation for $\varepsilon = 10^{-1},$ % $10^{-2}$, and $10^{-3}\kern.5pt$? % \par % {\em \underline{Exercise $20.5$}. Exponential % ill-conditioning and pseudo-nonuniqueness.} % {\em (a)} % Solve $\varepsilon\kern .1pt y'' - x\kern .1pt y' +y= 1$, $x\in[-1,1]$ % with $y(-1) = 0$ and $y(1) = 0$ for % $\varepsilon = 1/8, 1/16, 1/32, \dots$ until you get a solution % that looks completely unlike the others. Mathematically, % this solution is not correct. However, it is a % {\em pseudo-solution} in the sense that the residual % $\varepsilon\kern .1pt y'' - x\kern .1pt y' + y-1$ % is very small. Confirm this by computing % $\max_{x\in[-1,1]}^{} |\kern .3pt \varepsilon \kern .1pt % y''-x\kern .1pt y' + y - 1|$ for each of % your values of\/ $\varepsilon$. % {\em (\kern .7pt b)} % Explain this effect by boundary layer analysis % as follows. Find a formula for outer solutions to this problem and % note that it matches the observed incorrect solution. % Now find formulas for inner % solutions corresponding to boundary layers at % $x\approx \pm 1$ and note that % these provide enough boundary conditions to match any choice of % the outer solution. % Thus although the linear operator $\cal L$ that % maps right-hand side functions to solutions of this BVP is % mathematically nonsingular, hence has no null function, it % has {\em pseudo null functions} with exponentially small residual. % This ODE BVP is exponentially ill-conditioned, and we may think of it as % exponentially close to being underdetermined. % (The reason this has happened is because of the sign change in the % coefficient of $y'$ in this equation at $x=0$. This is related % to the theories of {\em exponential dichotomy} and {\em pseudospectra}. % See Chapters 10 and 11 of Trefethen and Embree, {\em Spectra and % Pseudospectra,} Princeton, 2005.) % \par % {\em \underline{Exercise $20.6$}. Exponential ill-conditioning % of the adjoint problem.} % Show that the adjoint of the operator $\cal L$ % of Exercise 20.5 is ${\cal L}^* \kern -3pt: v\mapsto % \varepsilon\kern .1pt v'' + x\kern .1pt v' + 2v$ with boundary % conditions $v(\pm 1) = 0$. % Solve the problem ${\cal L}^*\kern -1pt v = 1$ numerically % for $\varepsilon = 1/8, 1/16, 1/32, \dots.$ What goes wrong now when % when $\varepsilon$ gets sufficiently small? %