%% 22. Time-dependent PDEs %% % % \setcounter{page}{278} % %% % % Many ODE BVPs that arise in applications % are related to time-dependent PDE\kern .5pt s.\ \ In the simplest case % the solution $y(x)$ of the BVP is a {\bf time-independent} % or {\bf steady solution} % of the PDE, that is, a solution $y(x,t)$ of the PDE that is % independent of $t$. Often a solution of % a PDE may converge to time-independent form as % $t\to\infty$. % %% % % \smallskip % {\em Laplace and heat equations.} % The archetypical pair of equations related in this way % are the Laplace equation, an ODE BVP when restricted to one % space dimension, and the heat equation, a PDE IBVP % ({\em initial boundary-value problem}). Specifically, % the {\bf 1D Laplace equation} is the equation % $y'' = 0$, that is, $d^{\kern .7pt 2}y/dx^2 = 0$. In this chapter, % since time derivatives will also come into play, % we will write this in the notation $y_{xx}^{} = 0$, where $y_{xx}^{}$ % is an abbreviation for $\partial^2y/\partial x^2.$\ \ For % example, here is the 1D Laplace equation with % homogeneous Dirichlet boundary conditions on $[-\pi,\pi]$: % $$y_{xx}^{} = 0, \quad x\in [-\pi,\pi], ~ y(\pm\pi) = 0. % \eqno (22.1)$$ % We hardly need a computer to see that the solution is % $y(x) = 0.$ % %% % % Now let us consider how (22.1) may arise in % describing steady solutions of the {\bf 1D heat equation,} % $u_t^{} = u_{xx}^{}$. Here % is an {\bf initial boundary-value problem (IBVP)} % for this equation with a particular choice of initial data, % a pulse centered at $x=1$: % $$u_t^{} = u_{xx}^{}, \quad x\in [-\pi,\pi], ~ u(\pm\pi,t) = 0, ~ % u(x,0) = \exp(-50(x-1)^4). \eqno (22.2)$$ % Physically, the zero solution just mentioned represents the % effect that eventually, all the heat % described by (22.2) will flow out the ends of the % interval, where the temperature is held at zero. % A plot of the initial condition together with % the solution at times $t = 0.01$, $0.1$, and $1$ shows the % beginning of this process.\footnote{As the experiments of this % chapter show, Chebfun can be very effective in solving % certain PDE\kern .5pt s. PDE\kern .5pt s are not Chebfun's main focus, however, % and its syntax is not as simple nor its capabilities as % comprehensive for PDE\kern .5pt s as for ODE\kern .5pt s.} % In this and other figures of this chapter displaying % PDE solutions at various times~$t$, we plot the curves in % orange, with the initial curve in brown. % u0 = chebfun('exp(-50*(x-1)^4)',[-pi,pi]); pdefun = @(t,x,u) diff(u,2); bc.left = @(t,u) u; bc.right = @(t,u) u; opts = pdeset('plot','off'); t = [0 .01 .1 1]'; [t,u] = pde15s(pdefun,t,u0,bc,opts); plot(u,CO,ibvp), hold on, plot(u0,CO,ibvp0), hold off title(['Fig.~22.1.~~Solution to (22.2) ' ... 'at $t=0, 0.01, 0.1, 1$'],FS,11), ylim([-.2 1.2]) %% % \vskip 1.02em %% % It is apparent in this figure that until the % last curve, the boundaries have had little effect. We can % confirm this by integrating the curves to determine the % total heat at each time, which doesn't diminish % much until $t = 1$. disp(' t heat'), disp([t sum(u)']) %% % % \noindent % As $t$ increases further, the left boundary begins to be % important as well as % the right one, as is clear from these images for % $t = 1$ (again) and $2,4,8,16.$ % u = u(:,end); t = [1 2 4 8 16]'; [t,u] = pde15s(pdefun,t,u,bc,opts); plot(u,CO,ibvp), hold on, plot(u(:,1),CO,ibvp0), hold off title('Fig.~22.2.~~Solution continued to $t=1,2,4,8,16$',FS,11) %% % \vskip 1.02em %% % % \noindent % Note how the pulse is now moving to the center, approaching its % asymptotic form $C \exp(-t/4) \cos(x/2)$, as determined % originally by Joseph Fourier around 1807. % The integrals confirm that the total heat content is % now diminishing speedily, with $98\%$ of the heat % gone by $t=16$. % disp(' t heat 0.76exp(-t/4)') disp([t sum(u)' 0.76*exp(-t/4)]) %% % % The Laplace and heat equations may be the archetypes, but % these simple linear equations % give only a hint of the rich relationship between % time-dependent PDE\kern .5pt s and associated ODE BVPs. % We will now look at some further % examples illustrating the phenomena of stability, % instability, and symmetry-breaking bifurcation. % %% % % \smallskip % {\em Stability and instability.}\ \ Let us add % an exponential term to (22.1) and (22.2) to make them % into ODE and PDE versions of the nonlinear % {\em Bratu equation}, considered % in Chapters~16 and~18.\ \ Here is (16.4) again, with the constant % on the nonlinear term % adjusted to compensate for the change in domain from % $[\kern .3pt 0,1]$ to $[-\pi,\pi]$. % $$y_{xx}^{} + {3e^y\over 4\pi^2} = 0 , % \quad x\in [-\pi,\pi], ~ y(\pm\pi) = 0 . \eqno (22.3)$$ % As shown in Chapter 16, there are two solutions to % this ODE BVP.\ \ We plot one of them as a dashed line, the % one we shall find is unstable. % N = chebop(-pi,pi); N.op = @(x,y) diff(y,2) + (.75/pi^2)*exp(y); N.lbc = 0; N.rbc = 0; y1 = N\0; plot(y1,CO,bvpnl) hold on, x = chebfun('x',[-pi pi]); N.init = .2*(pi^2-x^2); y2 = N\0; plot(y2,'--',CO,bvpnl), ylim([-.5 2.5]), hold off title(['Fig.~22.3.~~Two solutions of the ' ... 'Bratu equation (22.3)'],FS,11) %% % \vskip 1.02em %% % % \noindent % Now let us embed the problem in a time-dependent PDE with % initial condition $u(x,0) = u_0^{}(x)$, % $$u_t^{} = u_{xx}^{} + {3e^u\over 4\pi^2}, % \quad x\in [-\pi,\pi], ~ u(\pm\pi,t) = 0,~u(x,0) = u_0^{}(x) . \eqno (22.4)$$ % First of all, suppose we take $u_0^{}(x) = 0$. Then as $t\to\infty$, % the solution increases and % approaches a steady state, the lower solution of (22.3) just plotted. % u0 = 0*x; pdefun = @(t,x,u) diff(u,2)+(.75/pi^2)*exp(u); t = (0:2:12)'; [t,u] = pde15s(pdefun,t,u0,bc,opts); plot(u,CO,ibvp), hold on, plot(u0,CO,ibvp0) plot(y1,CO,bvpnl), plot(y2,'--',CO,bvpnl), hold off, ylim([-.5 2.5]) title('Fig.~22.4.~~Solution to (22.4) at $t=0, 2,\dots,12$',FS,11) %% % \vskip 1.02em %% % % \noindent % Many other initial conditions also lead to the same steady % state, which is a stable solution of this PDE.\ \ For example, % in this plot we choose $u_0^{}$ to have % a few wiggles and include extra values $t = 0.05, 0.2$ % to show how transient their influence is. % u0 = sin(x)+.3*sin(3*x)+.1*sin(13*x); t = [0 0.05 0.2 2:12]'; [t,u] = pde15s(pdefun,t,u0,bc,opts); plot(u,CO,ibvp), hold on, plot(u0,CO,ibvp0), ylim([-1.5 2.5]) plot(y1,CO,bvpnl), plot(y2,'--',CO,bvpnl), hold off title(['Fig.22.5.~~Solution to (22.4) at ' ... '$t=0,.05,.2$ and $2,4,\dots,20$'],FS,11) %% % \vskip 1.02em %% % % As another example, next % we take $u_0^{}$ to be $0.9$ times the {\em upper\/} % solution of the BVP (22.3). Despite that starting point, % the PDE converges to the other, lower solution. The % upper solution of the BVP is an unstable solution of the PDE. % u0 = 0.9*y2; t = (0:4:32)'; [t,u] = pde15s(pdefun,t,u0,bc,opts); plot(u,CO,ibvp), hold on, plot(u0,CO,ibvp0), ylim([-.5 2.5]) plot(y1,CO,bvpnl), plot(y2,'--',CO,bvpnl), hold off title(['Fig.~22.6.~~Solution to (22.4) at ' ... '$t=0, 4,\dots,32$'],FS,11) %% % \vskip 1.02em %% % % \noindent % If we take $u_0^{}$ to be a little larger than the upper % solution of the BVP, the same instability is manifested as % divergence to $\infty$ in finite time, at $t\approx 6.914$ % for this particular initial condition. % u0 = 1.1*y2; t = [0 4 6 6.5 6.8 6.9 6.913]'; [t,u] = pde15s(pdefun,t,u0,bc,opts); plot(u,CO,ibvp), hold on, plot(u0,CO,ibvp0) plot(y1,CO,bvpnl), plot(y2,'--',CO,bvpnl), hold off title(['Fig.~22.7.~~Solution to (22.4) at ' ... '$t=0,4,6,6.5,6.8,6.9,6.913$'],FS,11) %% % \vskip 1.02em %% % % As we have discussed at many points in this book, phenomena % of stability and instability are of pervasive importance in % the mathematical sciences. Generally speaking, a stable % solution is likely to be observable in an experiment and % an unstable one is not. Also, generally speaking, a stable % solution is likely to emerge as a steady state of a % time-dependent process and an unstable one is not. % Still, a study of unstable solutions may be % important to a full understanding of a problem, and they % can often be determined by solving time-independent BVPs. % %% % % \smallskip % {\em Symmetry-breaking bifurcation.} % We have seen in Chapters 16--18 that the number of solutions % of an ODE BVP may change as a parameter varies. % Often when new solutions appear, this corresponds % to a loss of stability of an old solution, and % the new solutions are less symmetric % than the original one.\footnote{Symmetry-breaking % bifurcations are at the root of {\em phase transitions} in % physics. Physicists believe that some particularly important % phase transitions occurred soon after the Big Bang, when a % single unified force separated into the different fundamental % forces recognized today such as the electromagnetic and weak forces.} % %% % % To illustrate this effect, let us look at the % time-dependent {\bf Allen--Cahn equation}, % $$u_t^{} = \varepsilon\kern .3pt u_{xx}^{} + u - u^3, \quad % x\in [-1,1], ~ u(\pm 1,t) = 0, ~u(x,0) = u_0^{}(x). \eqno (22.5)$$ % For $\varepsilon > \varepsilon_1^{} \approx 0.406$, % the unique steady solution % of this problem is $u\equiv 0$, as we illustrate here with % $\varepsilon = 0.6$. % The associated ODE BVP, which appeared in Exercise~18.1, is % $$\varepsilon\kern .3pt y_{xx}^{} + y - y^3 = 0, \quad % x\in [-1,1], ~ y(\pm 1) = 0, \eqno (22.6)$$ % and we show its (zero) solution superimposed on the plot. % x = chebfun('x'); ep = 0.5; pdefun = @(t,x,u) ep*diff(u,2)+u-u.^3; u0 = 0.5*cos(pi*x/2)+.7*sin(pi*x); t = (.5:.5:4)'; [t,u] = pde15s(pdefun,t,u0,bc,opts); plot(u,CO,ibvp), hold on, plot(u0,CO,ibvp0), ylim([-1.3 1.3]) title(['Fig.~22.8.~~Convergence to 0 for Allen--Cahn ' ... 'eq.~(22.5) with $\varepsilon=0.6$'],FS,11) y1 = 0*x; plot(y1,CO,bvpnl), hold off %% % \vskip 1.02em %% % % Reducing $\varepsilon$ to $0.15$, on the other hand, illustrates % a different situation. Here, the zero solution is still % valid, but it is unstable. With the same initial condition % as before, the curve approaches a positive form as $t\to\infty$. % By symmetry, if the initial condition were negated, the % steady solution as $t\to\infty$ would be negated too. % ep = 0.15; pdefun = @(t,x,u) ep*diff(u,2)+u-u.^3; [t,u] = pde15s(pdefun,t,u0,bc,opts); plot(u,CO,ibvp), hold on, plot(u0,CO,ibvp0), ylim([-1.3 1.3]) title(['Fig.~22.9.~~New solutions and broken symmetry ' ... 'with $\varepsilon=0.15$'],FS,11) N = chebop(-1,1); N.lbc = 0; N.rbc = 0; N.op = @(x,y) ep*diff(y,2) + y - y^3; N.init = u0; y2 = N\0; y3 = -y2; plot(y1,CO,bvpnl,'--'), plot([y2 y3],CO,bvpnl), hold off %% % \vskip 1.02em %% % What if one starts from an initial condition that is an odd function % of $x$, % hence equally close to the upper and lower solutions? For this % value of $\varepsilon$, the solution quickly converges to zero. u0 = .7*sin(pi*x); [t,u] = pde15s(pdefun,t,u0,bc,opts); plot(u,CO,ibvp), hold on, plot(u0,CO,ibvp0), ylim([-1.3 1.3]) title(['Fig.~22.10.~~Odd initial function: ' ... 'convergence to a saddle point'],FS,11) plot(y1,CO,bvpnl,'--'), plot([y2 y3],CO,bvpnl), hold off %% % \vskip 1.02em %% % % \noindent % The zero solution is not truly stable, however. It is a saddle % point, and if it is perturbed by a function that is not odd, the % perturbations will generally grow. % %% % If $\varepsilon$ is reduced further for this problem, % below another bifurcation point $\varepsilon_2^{}\approx 0.1$, a further % set of solutions appears. ep = 0.02; pdefun = @(t,x,u) ep*diff(u,2)+u-u.^3; u0 = .3*sin(pi*x); [t,u] = pde15s(pdefun,t,u0,bc,opts); plot(u,CO,ibvp), hold on, plot(u0,CO,ibvp0), ylim([-1.3 1.3]) title(['Fig.~22.11.~~Two more unstable steady states ' ... 'with $\varepsilon=0.02$'],FS,11) N.op = @(x,y) ep*diff(y,2) + y - y^3; N.init = cos(pi*x/2); y2 = N\0; y3 = -y2; N.init = sin(pi*x); y4 = N\0; y5 = -y4; plot(y1,CO,bvpnl,'--'), plot([y2 y3],CO,bvpnl) plot([y4 y5],'--',CO,bvpnl), hold off %% % \vskip 1.02em %% % % \noindent % We have marked the new solutions with dashed curves, because % in fact these again are saddle points, not stable points. % If the initial condition were not an odd function, the solution % would eventually converge to one of the stable steady states % of a single sign. % %% % % \smallskip % {\em Time-harmonic solutions and resonance.} % Up to now we have looked at PDE\kern .5pt s of first order in $t$.\ \ An % ODE BVP may also be relevant to time-oscillatory % solutions of a PDE that is of second order in $t$. % (We considered such problems in Exercise~6.3 and % in the Application of Chapter~20, % Why is New York hotter than San Francisco?'') % Suppose we have a linear PDE % $$u_{tt}^{} = Lu + e^{i\kern .3pt\omega\kern .3pt t} f , \eqno (22.7)$$ % where $L$ is a linear operator acting on functions of $x$, % $\omega$ is frequency, % and $f$ is a function of $x$. % A {\bf time-harmonic} solution of (22.7) is a solution % of the form % $$u(x,t) = e^{i\kern .3pt\omega\kern .3pt t} y(x) \eqno (22.8)$$ % for some function $y$. Inserting % (22.8) in (22.7) gives % $$-\omega^2 u = Lu + e^{i\kern .3pt\omega\kern .3pt t} f ,$$ % or after dividing by $e^{i\kern .3pt\omega\kern .3pt t}$, % $$-\omega^2 y = Ly + f , \eqno (22.9)$$ % known as the {\bf reduced equation} associated with (22.7). % One of the fundamental techniques of mathematical physics % is to solve reduced equations to find time-harmonic solutions % of time-dependent equations.\footnote{In quantum mechanics, % the technique is applied to the Schr\"odinger equation, % which differs from (22.7) in being of first order in $t$. We saw % this in the Application of Chapter~6, Eigenstates of % the Schr\"odinger equation.''} % %% % % \smallskip % {\em Further links between ODE\kern .5pt s and PDE\kern .5pt s.} % The difference between ODE\kern .5pt s and PDE\kern .5pt s is nothing more than % the number of independent variables, so it is hardly surprising % that connections between the two arise in many ways. In this % chapter we have emphasized steady-state forms that a solution % of a PDE may settle down to as $t\to\infty$. % A related situation is that % sometimes as $t\to\infty$ the solution of a PDE approaches not % a fixed function but a {\em traveling wave,} a function % $u(x,t) = y(x-c\kern .7pt t)$ for some constant wave velocity $c$ % that in most cases is not known a priori. Here again ODE\kern .5pt s play % a natural role. Another way ODE\kern .5pt s arise from PDE\kern .5pt s is in % {\em separation of variables,} and the time-harmonic solutions % just looked at fall in this category. % Separation of variables applies also in purely spatial % problems, for example, if a radially-symmetric elliptic PDE is % factored into an ODE in $r$ and an ODE in $\theta$\kern .5pt: this is the % most familiar way in which one encounters the Bessel equation (Exercise~22.1). % Another PDE-ODE link arises for problems that have {\em similarity % solutions.} And if we move up one or more dimensions, it is worth % noting that just as a PDE in one space dimension may reduce to % an ODE in % some limit, so a PDE in two or more space dimensions may reduce % to a PDE in one dimension less. For example, steady solutions of % the heat equation in a 2D spatial domain will be described by solutions % of the Laplace equation, a PDE, in the same domain. % %% % % \begin{center} % \hrulefill\\[1pt] % {\sc Application: solitons and the KdV equation} \\[-3pt] % \hrulefill % \end{center} % %% % % The simplest of all time-dependent PDE\kern .5pt s is the {\em first-order % 1D linear wave equation} (or {\em advection equation}), % $$u_t^{} = -u_x^{}, \eqno (22.10)$$ % whose solutions consist of an arbitrary wave form $y(x)$ traveling % rightward at speed~$1$ without changing shape, % $$u(x,t) = y(x-t). \eqno (22.11)$$ % For example, here we start with a Gaussian of amplitude % $50$ centered at $x=-1$ and run to $t=2$ using % the Chebfun {\tt spin} command (stiff PDE integrator''). % The wave slides along to $x=1$ without changing shape. % dom = [-5,5]; x = chebfun('x',dom); u0 = 50*exp(-30*(x+1)^2); S = spinop(dom,[0 2]); S.lin = @(u) -diff(u); S.init = u0; u1 = spin(S,160,.01,'plot','off'); plot(u0,'r',u1,'k',LW,1), ylim([-20 130]) text(-1,57,'$t=0$',FS,9,HA,CT), text(1,57,'$t=2$',FS,9,HA,CT) title(['Fig.~22.12.~~Advection with the linear ' ... 'wave equation (22.10)'],FS,11) %% % \vskip 1.02em %% % % For a nonlinear equation, on the other % hand, change of shape will be the % rule. In this application we look at % the most famous of all nonlinear wave % equations, the {\em KdV} ({\em Korteweg--de Vries\/}) equation, % which we write in the form % $$u_t^{} = -0.03(u^2)_x^{} - 0.01 u_{xxx}^{}. \eqno (22.12)$$ % The equation blends first-order nonlinear {\em advection}, involving % $(u^2)_x^{}$ or equivalently $2\kern .3pt uu_x^{}$, with third-order % linear {\em dispersion}, involving $u_{xxx}^{}$. % The constants $0.03$ and $0.01$ have no great significance % and have been chosen to make the behavior approximately match % that of $u_t^{} = -u_x^{}$ for our initial condition. % In the next figure we solve (22.12) to $t=2$ as before with {\tt spin}, % which is made for nonlinear problems like this. % S = spinop(dom,[0 2]); S.lin = @(u) -0.01*diff(u,3); S.nonlin = @(u) -0.03*diff(u.^2); S.init = u0; u1 = spin(S,160,.01,'plot','off'); plot(u0,'r',u1,'k',LW,1), ylim([-20 130]) text(-1,57,'$t=0$',FS,9,HA,CT), text(.78,50,'$t=2$',FS,9,HA,CT) title(['Fig.~22.13.~~Nonlinear propagation with ' ... 'the KdV equation (22.12)'],FS,11) %% % \vskip 1.02em %% % % To first approximation, Figure~22.13 is like Figure~22.12. % However, the peak is moving a little more slowly, and it is losing % amplitude. Meanwhile a tail of dispersive % oscillations is appearing behind. (The small oscillations % ahead of the main pulse are due to the use of periodic % boundary conditions in all {\tt spin} simulations.) A complicated mix of effects % like this is what one must expect to see, in general, with a nonlinear PDE. % %% % % But the KdV equation has special properties that % this experiment does not reveal. It has special solutions, % called {\em solitary waves\/} or {\em solitons}, % that travel at fixed speed with exactly % uniform shape. For any constant $c>0$, the associated soliton % profile is given by the formula % $$y(x) = 50\kern .7pt c\kern 1.3pt \hbox{sech}(\sqrt{100\kern .3pt c} % \kern .8pt (x+1)/2)^2 , \eqno (22.13)$$ % where $\hbox{sech}(x) = 1/\cosh(x)$, the reciprocal of the hyperbolic % cosine. % Moreover, the speed at which the soliton moves % is equal to $c$, which means that tall solitons move % faster than short ones. Here we see this behavior with % a simulation for $t\in [\kern .3pt 0,2\kern .2pt ]$ and $c=1$. Note that the % $\hbox{sech}^2$ function is much like a Gaussian, the difference % in shape being too little to notice by eye. % c = 1; soliton1 = 50*c*sech(sqrt(100*c)*(x+1)/2)^2; S.init = soliton1; u1 = spin(S,160,.01,'plot','off'); plot(S.init,'r',u1,'k',LW,1), ylim([-20 130]) text(-1,57,'$t=0$',FS,9,HA,CT), text(1,57,'$t=2$',FS,9,HA,CT) title(['Fig.~22.14.~~If the initial condition is a soliton, ' ... 'it translates unchanged'],FS,11) %% % \vskip 1.02em %% % % \noindent % Since we took $c=1$, this soliton has traveled at speed $1$. % If we double the height, we double the speed. % c = 2; soliton2 = 50*c*sech(sqrt(100*c)*(x+2)/2)^2; S.init = soliton2; u1 = spin(S,160,.005,'plot','off'); plot(S.init,'r',u1,'k',LW,1), ylim([-20 130]) text(-2,107,'$t=0$',FS,9,HA,CT), text(2,107,'$t=2$',FS,9,HA,CT) title('Fig.~22.15.~~A taller soliton moves faster',FS,11) %% % \vskip 1.02em %% % % Another property of the KdV equation is truly extraordinary. % {\em Solitons can pass through one another.} For example, if % we start with a solution consisting of two well separated solitons, % the rear one taller than the front one, then the rear % one will catch up, pull % ahead, and eventually return to the same shape and % speed it had before. Here is an illustration. % S.init = soliton1 + soliton2; u1 = spin(S,160,.005,'plot','off'); plot(S.init,'r',u1,'k',LW,1), ylim([-20 130]) text(-1.3,64,'$t=0$',FS,9,HA,CT), text(1.4,64,'$t=2$',FS,9,HA,CT) title('Fig.~22.16.~~A fast soliton passes through a slower one',FS,11) %% % \vskip 1.02em %% % % \noindent % Note the locations of the peaks in this figure % for the curve at $t=2$. The speeds have returned to % essentially their initial values 2 and 1, but % the tall soliton is a bit further % along than $x=2$, and the short soliton is not as far along % as $x=1$. This is because a shift of each soliton is % introduced during the time of their interaction around $t=1$. % The details of these nonlinear interactions are subtle and fascinating, % as the reader can easily explore with further computations. % A waterfall plot gives the idea. % spin(S,160,.005,'plot','waterfall'); set(gca,FS,9) title('Fig.~22.17.~~Two solitons',FS,11), colormap([0 0 0]) xlabel('$x$',FS,10), ylabel('$t$',FS,10), zlabel(' ') %% % \vskip 1.02em %% % % A soliton is an example of a traveling wave, which we % defined earlier as a solution of % a time-dependent PDE of the form % $$u(x,t) = y(x-c\kern .7pt t) \eqno (22.14)$$ % for some function $y(x)$. % Of course we know that in the case of the KdV % equation, the formula for $y$ must be (22.13), but % let us suppose we do not know this. % Inserting (22.14) in (22.12) gives % $-c\kern .7pt y_x^{} = -0.03 \kern .5pt (y^2)_x^{} - 0.01 y_{xxx}^{},$ % that is, % $$y_{xxx}^{} = 100\kern .7pt c % \kern .7pt y_x^{} -3\kern .5pt (y^2)_x^{}.$$ % This looks like a nonlinear 3rd-order ODE.\ \ However, % all three terms are differentiated, so we % can integrate it once to get the second-order equation % $$y_{xx}^{} = 100\kern .7pt c\kern .7pt y-3\kern .5pt y^2 + C$$ % for some constant $C$. Now let us use this ODE to find the % shape of a traveling wave solution. If we assume that $y$ and % its derivatives approach 0 as $|x|\to \infty$, then we have $C = 0$, % and the equation becomes % $$y_{xx}^{} = 100\kern .7pt c\kern .7pt y-3\kern .5pt y^2. % \eqno (22.15)$$ % This could be dealt with analytically, or one % could try various things on the computer. For example, suppose % we regard (22.15) as an IVP on the interval $x\in [-2,2\kern .2pt ]$, with % $c=1$. With initial conditions $y(-2)=0$ and % $y'(-2)= 10^{-4}$ we get a soliton of just the expected amplitude % $50$.\ \ The horizontal location has no particular significance; it would be % different with a different value of $y'(-2)$. % N = chebop(-2,2); N.op = @(y) diff(y,2) - 100*y + 3*y^2; N.lbc = [0; 0.0001]; y = N\0; plot(y,CO,ivpnl), axis([-5 5 -20 130]) title(['Fig.~22.18.~~Soliton travelling wave profile ' ... 'determined from an ODE'],FS,11) %% % \vskip 1.02em %% % % We have not mentioned the scientific uses of the KdV equation. % These began with water waves, and % the story is often told of John Scott Russell's chase on horseback of a % solitary wave in % the Union Canal in Scotland in 1834. Two centuries later, the % important applications are in optics. Light pulses travelling through % optical fibers experience nonlinear effects combined with dispersion, and % the KdV equation and its relatives are fundamental nowadays in % optical technology. % %% % % \smallskip % {\sc History.} There are few equations with more impact in % the history of mathematics than the KdV equation. % The understanding of solitons and their interactions came in the 1960s % and led to sensational new discoveries related to so-called % completely integrable dynamical systems; some of the key names % were Zabusky, Kruskal, and Lax. % The technical details are exhilarating, and so are the wider % implications, for it was work around the KdV equation that led to % some of our deepest understanding of the significance of nonlinearity. % All this was sparked by computers. Though the KdV equation goes % back to the 19th century, the great advances began % with a 1955 paper by Fermi, Pasta, and Ulam, % who worked on the Maniac I computer at Los Alamos. % \smallskip % %% % % {\sc Our favorite reference.} % {\em The PDE Coffee Table Book} was a group % project at Oxford in the late 1990s coordinated by % Trefethen together with Kristine Embree. % The idea was to present 100 PDE\kern .5pt s in the simplest % possible terms, each getting exactly a two-page spread with % just the right information about mathematics, applications, and history, % and always with a colorful computed illustration. Alas, the % project halted when it was only $34\%$ complete, but that % $34\%$ is worth checking out: see % {\em The (Unfinished) PDE Coffee Table Book,} % \verb|https://people.maths.ox.ac.uk/trefethen/pdectb.html|. % \smallskip % %% % % \begin{displaymath} % \framebox[4.7in][c]{\parbox{4.5in}{\vspace{2pt}\sl % {\sc Summary of Chapter 22.} Many ODE BVPs are % related to associated time-dependent PDE\kern .5pt s in one space dimension. % The ODE\kern .5pt s may arise as limits\/ $t\to\infty$ in which solutions % approach a steady state, a steady oscillation, or a traveling wave. % ODE\kern .5pt s may also arise from purely spatial PDE\kern .5pt s by the process % of separation of variables. % \vspace{2pt}}} % \end{displaymath} % %% % \smallskip\small\parskip=1pt\parindent=0pt % {\em Exercise $22.1$. Bessel equation for a circular drum.} % A 2D membrane with a fixed boundary $\Gamma$ oscillates according % to the wave equation $u_{tt}^{} = \Delta u$, where $\Delta$ is the % Laplacian operator, $\Delta u = u_{xx}^{} + u_{yy}^{}$. This is a PDE % in two space variables and time. % {\em (a)}~Assume that $u$ takes the form $u(x,y,t) = % \exp(i\kern .5pt \omega\kern .3pt t) v(x,y)$ for a fixed frequency~$\omega$ and function % $v$. Show that this assumption % reduces the equation to a PDE eigenvalue equation in two space variables. % {\em (\kern .7pt b)} Now assume that $\Gamma$ is the unit circle and % that~$v$ takes the polar coordinates % form $v(r,\theta) = \exp(i k \theta) y(r)$ for % a fixed integer wave number~$k$. Using the polar coordinates % representation of the Laplacian, % $\Delta v = v_{rr}^{} + r^{-1} v_r^{} + r^{-2} v_{\theta\theta}^{}$, % show that this reduces the problem % to an ODE eigenvalue problem involving Bessel's equation. % \par % {\em \underline{Exercise $22.2$}. % Asymptotic behavior of advection-diffusion problem.} % {\em (a)} Execute {\tt chebgui} and run the demo labeled % Advection-diffusion equation 2 under the tab Demos/PDE-scalar. % The IBVP here is $u_t^{} = 0.3\kern .2pt u_{xx}^{} + 10\kern .3pt u_x^{}$ for % $x\in[-1,1]$ with initial condition $u(x,0) = \exp(-10x^4/(1-x^2))$ % and boundary conditions $u(\pm 1) = 0$. % Describe in words how the solution behaves as $t\to\infty$. % {\em (\kern .7pt b)} This behavior can be analyzed by looking % for solutions of this equation of the form $u(x,t) = e^{\lambda t} % y(x)$, where~$y$ is a fixed function of $x$. Write down the % formulas showing that $y$ and $\lambda$ are solutions of an % eigenvalue problem. Determine $y$ and $\lambda$, either % analytically or numerically, for the value of $\lambda$ corresponding to % the slowest decay of $u(x,t)$ as $t\to\infty$. Show % that~$y$ has a shape corresponding to the observations of part {\em (a)}. % {\em (c)} By pressing the Export solution'' button in Chebgui and % then typing \verb|max(u(:,end))|, find the value of % $\max_{x\in[-1,1]} u(x,0.25)$ (since the final time of this simulation % is $t=0.25$). Combine this figure with the result of part {\em (\kern .7pt b)} % to estimate $\max_{x\in[-1,1]} u(x,5)$ to an accuracy of within % $1\%$. {\em (d)} The asymptotic rate of decay is different if the advection % term $10\kern .3pt u_x^{}$ is removed from the equation. (Physically, % heat transfer is faster with convection than just conduction --- the % oceanic effect that we found makes San Francisco cooler than New York in % the summer.) % Estimate $\max_{x\in[-1,1]} u(x,5)$ in this case. % \par % {\em \underline{Exercise $22.3$}. Metastability of Allen--Cahn equation.} % {\em (a)} Execute {\tt chebgui} or {\tt pde15s}, as you % prefer, for the Allen--Cahn equation % $u_t^{} = 0.05\kern .3pt u_{xx}^{} + u = u^3$ for $x\in [-5,5]$ and % $t\in [\kern .3pt 0,50\kern .3pt ]$ with boundary % conditions $u(-5) = u(5) = -1$ and initial condition % $u(x,0) = -1+2\exp(-x^2)$. Show how the initial function quickly changes % to a stalagmite'' form that looks like a steady state. % {\em (\kern .7pt b)} Now run the calculation for % $t\in [\kern .3pt 0,1000\kern .3pt ]$ and show that the apparent steady % state was not truly steady, only metastable. Estimate the time at % which the stalagmite is extinguished --- say, the time $t$ at which % $\max_{x\in[-5,5]} u(x,t)$ becomes negative. How does this critical time change % if the diffusion coefficient $0.5$ is increased to $0.6$ or $0.7\kern .3pt$? % (With values of the diffusion coefficient much less than $0.5$ the % time soon becomes nearly infinite in practice --- too long to be % measured by a Chebfun computation --- though it is always finite in % principle.) % {\em (c)}~Write down the ODE in the variable $y$ for a time-independent % solution of the Allen--Cahn equation. One solution of the ODE is % $y(x) = -1$. The waveform observed for $t=50$ in part {\em (a)} % is a pseudo-solution'': it nearly satisfies the ODE, though not % exactly. Insert this % solution in the ODE to quantify nearly.'' %