%% 19. Periodic ODEs %% % % \setcounter{page}{241} % %% % % Suppose we have an ODE defined for % $t\in (-\infty,\infty)$ whose coefficients are $2\pi$-periodic, % or in short, a {\bf periodic ODE.}\ \ Will the solutions % be $2\pi$-periodic too? % %% % % Almost any experiment will show you that % they don't {\em have} to be. For example, % consider the problem % $$y' + (1+\cos(t)) y = 1, \quad y(0) = 0. \eqno (19.1)$$ % By Theorem 11.1, a unique solution exists for all $t$, negative % and positive. This solution is % not periodic, however, as we see by plotting it % for $t\in [\kern .3pt 0,10\pi]$. % ODEformats L = chebop(0,10*pi); L.op = @(t,y) diff(y) + (1+cos(t))*y; L.lbc = 0; y = L\1; plot(y,CO,ivp), ylim([0 3]) title('Fig.~19.1.~~Nonperiodic solution to (19.1)',FS,11) pilabels = {'0','2\pi','4\pi','6\pi','8\pi','10\pi'}; set(gca,FS,8,XT,pi*(0:2:10),XTL,pilabels), set(gca,YT,0:3) %% % \vskip 1.02em %% % % \noindent % Although this solution isn't periodic, it is very close to % periodic after an initial phase. It would seem % that there exists a truly periodic solution to this ODE that % will oscillate between extremes of about $0.5$ and $2.4$. % (From now on in this chapter, for simplicity, {\bf periodic} usually means % $2\pi$-periodic.) % Chebfun will find this solution if you % specify the boundary condition \verb|'periodic'|. % L = chebop(0,2*pi); L.op = @(t,y) diff(y) + (1+cos(t))*y; L.bc = 'periodic'; y = L\1; plot(y,CO,bvp), ylim([0 3]) title('Fig.~19.2.~~Periodic solution to (19.1)',FS,11) pilabels2 = {'0','\pi/2','\pi','3\pi/2','2\pi'}; set(gca,FS,8,XT,pi*(0:.5:2),XTL,pilabels2), set(gca,YT,0:3) %% % \vskip 1.02em %% % % \noindent % The minimum and maximum are as expected, % ymin = min(y), ymax = max(y) %% % % \noindent % and the value at $t=0$ is a little bit above the minimum, % y0 = y(0) %% % % \noindent % If we solve the IVP on $[\kern .3pt 0,10\pi]$ with this value % specified as the initial condition, the periodic behavior is seen. % L = chebop(0,10*pi); L.op = @(t,y) diff(y) + (1+cos(t))*y; L.lbc = y0; y10pi = L\1; plot(y10pi,CO,ivp), ylim([0 3]) title('Fig.~19.3.~~More periods of the same solution',FS,11) set(gca,FS,8,XT,pi*(0:2:10),XTL,pilabels), set(gca,YT,0:3) %% % \vskip 1.02em %% % % Encouraged by this experiment, let us try a second-order periodic equation, % $$y'' + (1+\cos(t)) y = 1, \quad y(0) = y'(0) = 0. \eqno (19.2)$$ % Again the first solution we find is nonperiodic. % L = chebop(0,10*pi); L.op = @(t,y) diff(y,2) + (1+cos(t))*y; L.lbc = [4;0]; y = L\1; plot(y,CO,ivp), ylim([-50 50]) title('Fig.~19.4.~~Nonperiodic solution to (19.2)',FS,11) set(gca,FS,8,XT,pi*(0:2:10),XTL,pilabels) %% % \vskip 1.02em %% % % \noindent % This time, the curve is not settling down to a periodic form, but % oscillating and growing exponentially. % Nevertheless, again there exists a periodic solution. % L = chebop(0,2*pi); L.op = @(t,y) diff(y,2) + (1+cos(t))*y; L.bc = 'periodic'; y = L\1; plot(y,CO,bvp), ylim([-4 4]) title('Fig.~19.5.~~Periodic solution to (19.2)',FS,11) set(gca,FS,8,XT,pi*(0:.5:2),XTL,pilabels2) %% % \vskip 1.02em %% % % \noindent % As before, we can recover periodic % oscillations over $[\kern .3pt 0,10\pi]$ if we set the % initial conditions appropriately. % %% % % Nonlinear periodic ODE\kern .5pt s often have periodic solutions too. % For example, the equation % $$y'' + (1+\cos(t)) (y + (y/4)^3) = 1 \eqno (19.3)$$ % is a nonlinear variant of (19.2). Here we compute a periodic % solution and immediately plot it on $[\kern .3pt 0,10\pi]$. % L = chebop(0,2*pi); L.op = @(t,y) diff(y,2) + (1+cos(t))*(y+(y/4)^3); L.bc = 'periodic'; y = L\1; y10pi = chebfun(@(t) y(t),[0 10*pi],'trig'); plot(y10pi,CO,bvpnl) title('Fig.~19.6.~~Periodic solution to (19.3)',FS,11) set(gca,FS,8,XT,pi*(0:2:10),XTL,pilabels), ylim([-4 4]) %% % \vskip 1.02em %% % % Of course, as usual, the behavior of nonlinear problems can vary % greatly. If $(y/4)^3$ is changed to $(y/3)^3$ in (19.3), experiments % indicate that there is no solution. % For another example, in Chapter 1 we considered the van der Pol % equation (1.2), % $$0.3y''-(1-y^2)y'+y = 0.$$ % Like all autonomous equations, this ODE is not just $2\pi$-periodic, % but $T$-periodic for any period $T$. For any $T$, it has the $T$-periodic % solution $y = 0$. For $T=2\pi$, and for most other values of $T$, % this is the only $T$-periodic solution. % For the special value $T\approx 4.0725$, however, % there is another nontrivial $T$-periodic solution % corresponding to the limit cycle plotted in Chapter~1. % In fact, this solution is one of an infinite family of % distinct $T$-periodic solutions, since it could be shifted in~$t$ by % any constant $\Delta t$. % %% % % Let us take stock of the situation. We have seen that % a periodic ODE may or may not have periodic solutions, and % that it will almost certainly have nonperiodic ones. % The remainder of this chapter is organized around three questions for % periodic ODE\kern .5pt s: % \par\medskip % 1. Does there exist a periodic solution? % \par\smallskip % 2. If so, is it unique? % \par\smallskip % 3. What can we say about nonperiodic solutions? % \label{question3} % \par\medskip % \noindent % These questions have rich histories going back to % the nineteenth century, and % as usual with ODE\kern .5pt s, there are applications involving % both time and space. For example, % periodicity in time is associated % with the dynamics of the solar % system or of rotating or oscillating % machinery, and periodicity in space is associated % with propagation of sound and light waves in crystals. % We will focus mainly on linear problems, with brief comments % on nonlinear ones at the end. With regret, we shall mainly % look just at questions (1) and (2), though the Application is % of type (3). % %% % % Before turning to second-order periodic equations, the traditional % focus in this subject area, let us consider the first-order case. % Following the notation of Theorem 2.3, % we start with the scalar linear problem % $$y' - a(t) y = g(t) \eqno (19.4)$$ % where $a$ and $g$ are assumed to be % continuous and $2\pi$-periodic. We % can answer questions (1)--(2) above without much difficulty. % First we note that by Theorem~10.1, there is a unique solution % $y(t)$ of (19.4) associated with any initial % condition $y(0) = y_0$. If $y(2\pi)\ne y(0)$, then by % definition, $y(t)$ is not $2\pi$-periodic. If $y(2\pi) = y(0)$, on % the other hand, then since the coefficients are periodic, the % solution on $[2\pi,4\pi]$ will be the same as on $[\kern .3pt 0,2\pi]$, % the solution on $[4\pi,6\pi]$ will be the same as on $[2\pi,4\pi]$, and % so on, implying that the solution on all % of $(-\infty,\infty)$ is $2\pi$-periodic. Thus the question % of $2\pi$-periodicity reduces to the question of whether % $y(2\pi) = y(0)$. % %% % One special case of (19.4) is the situation where % $a(t) = 0$, % $$y' = g(t). \eqno (19.5)$$ % Here the ODE reduces to an integral, % $$y(t) = y(0) + \int_0^t g(s)\kern .5pt ds , \eqno (19.6)$$ % implying that the solution for any initial condition % $y(0) = y_0$ will be periodic if and only if % $$\int_0^{2\pi} g(s)\kern .5pt ds = 0 . \eqno (19.7)$$ % Thus the answer to question (1) is that there is a periodic % solution if and only if (19.7) holds, and the answer to (2) % is that it is never unique since % we could always add a constant. %% % % The more substantial special case of (19.4) is the homogeneous % case $g(t) = 0$, % $$y' -a(t) y = 0. \eqno (19.8)$$ % Here, following Theorem 2.2 (separation of % variables), the solution is given by % $$y(t) = e^{h(t)} y(0), \eqno (19.9)$$ % where $h(t)$ is defined by % $$h(t) = \int_0^t a(s)\kern .5pt ds . \eqno (19.10)$$ % In particular, we have % $$y(2\pi ) = e^{h(2\pi )} y(0), \eqno (19.11)$$ % which implies that $y(2\pi) = y(0)$ if and only if % $$\left(1-e^{h(2\pi)}\right) y(0) = 0.$$ % This equation holds if and only if % $y(0)$ is zero or $\exp(h(2\pi))=1$, that is, % $$h(2\pi) = 2\pi i \times \hbox{integer}. \eqno (19.12)$$ % If (19.12) holds, (19.8) is said to be {\bf critical}; % otherwise it is {\bf noncritical}. % Thus we see that if % (19.8) is noncritical, then (19.8) has the unique periodic solution % $y(t) = 0\kern 1pt$: the % answers to (1) and (2) are yes and yes. If % it is critical, then any value $y(0)$ in (19.9) % gives a $2\pi$-periodic solution, so the answers are yes and no. % %% % % We shall now see that the behavior of the inhomogeneous problem (19.4) % depends on the same condition (19.12) of criticality. % By Theorem 2.3, the solution can be written % $$y(t) = e^{h(t)} y(0) + \kern -2pt \int_0^t e^{h(t)-h(s)} g(s)\kern .5pt ds, % \eqno (19.13)$$ % where $h(t)$ is still defined by (19.10). % Thus the value at $t = 2\pi$ is % $$y(2\pi) = e^{h(2\pi)} y(0) + \kern -2pt \int_0^{2\pi} % e^{h(2\pi)-h(s)} g(s)\kern .5pt ds, \eqno (19.14)$$ % so the solution will be $2\pi$-periodic if and only if % $$\left(1 - e^{h(2\pi)}\right) y(0) = % \int_0^{2\pi} e^{h(2\pi)-h(s)} g(s)\kern .5pt ds. \eqno (19.15)$$ % This leads to the following theorem encompassing the % two special cases just considered. % %% % % \smallskip % {\em % {\bf Theorem 19.1. First-order scalar linear periodic ODE\kern .5pt s % (\textsf{\textbf{FLaShi}}).} % Consider the equation $$y' - a(t) y = g(t), \eqno (19.16)$$ % where $a(t)$ and $g(t)$ are continuous and $2\pi$-periodic. % If the equation is % noncritical in the sense that\/ $(19.12)$ does not hold, % then there exists a unique $2\pi$-periodic solution~$y$. % If it is critical, then there exist infinitely many % $2\pi$-periodic solutions if the right-hand % side of $(19.15)$ is zero, and no $2\pi$-periodic solutions if % it is nonzero. % } % \smallskip % %% % % \def\qed{\vrule height4pt width3pt depth4pt} % {\em Proof.} If the equation is noncritical, we may divide % (19.15) by $1-e^{h(2\pi)}$ to get % $$y(0) = \left.\int_0^{2\pi} e^{h(2\pi)-h(s)} g(s)\kern .5pt ds\right/ % \kern -2pt \left( 1 - e^{h(2\pi)}\right) ,$$ % producing a $2\pi$-periodic solution (19.13) as claimed. % Uniqueness follows since this condition on $y(0)$ determines % $y(t)$ uniquely. % If the equation is critical, then the left-hand side of (19.15) is % zero regardless of the value of $y(0)$. If the right-hand side % is zero, then each choice of $y(0)$ gives a distinct $2\pi$-periodic % solution, whereas if the right-hand side is nonzero, no choice of % $y(0)$ gives a periodic solution. \qed % \smallskip % %% % % The last few pages have considered scalar problems. For % a system of equations, the pattern is similar. Extending % Chapter~14, consider the equation % $${\bf y}' - {\bf A}(t) {\bf y} = {\bf g}(t), \eqno (19.17)$$ % where ${\bf y}(t)$ for each $t$ is an $n$-vector, % ${\bf g}(t)$ is a given $2\pi$-periodic continuous % $n$-vector function of $t$, % and ${\bf A}(t)$ is a given $2\pi$-periodic continuous $n\times n$ matrix % function of $t$. % The general solution of (19.17) is % $${\bf y}(t)={\bf Y}(t){\bf c}+{\bf Y}(t)\int_0^t{\bf Y}^{-1}(s)\kern .5pt {\bf g}(s)\kern .5pt ds % \eqno (19.18)$$ % for $n$-vectors ${\bf c}$, where % ${\bf Y}(t)$ is the {\bf fundamental matrix}, the $n\times n$ % matrix function whose columns are the linearly independent solutions % of the homogeneous problem % ${\bf y}' - {\bf A}(t) {\bf y} = {\bf 0}$, % ${\bf y}^{(1)}(t),\dots,{\bf y}^{(n)}(t)$ % with ${\bf Y}(0) = {\bf I}$, the $n\times n$ identity. % The values at $t=0$ and $t = 2\pi$ are % $${\bf y}(0)={\bf c}$$ % and % $${\bf y}(2\pi)={\bf Y}(2\pi){\bf c}+{\bf Y}(2\pi)\int_0^{2\pi} % {\bf Y}^{-1}(s)\kern .5pt {\bf g}(s)\kern .5pt ds,$$ % so the solution will be $2\pi$-periodic if and only if % $$\left({\bf I}-{\bf Y}(2\pi)\right) {\bf c} = % {\bf Y}(2\pi)\int_0^{2\pi} {\bf Y}^{-1}(s)\kern .5pt {\bf g}(s)\kern .5pt ds. % \eqno (19.19)$$ % The criticality condition will be the condition that % ${\bf I}-{\bf Y}(2\pi)$ is singular, % or equivalently, % $${\bf Y}(2\pi) \hbox{ has an eigenvalue equal to\/ } 1. \eqno (19.20)$$ % We say that (19.17) is {\bf critical} if this condition % holds and {\bf noncritical} otherwise, i.e., all eigenvalues % of ${\bf Y}(2\pi)$ are different from $1$. % Theorem 19.1 generalizes as follows. % %% % % \smallskip % {\em % {\bf Theorem 19.2. First-order linear periodic systems of ODE\kern .5pt s % (\textsf{\textbf{FLashi}}).} % Consider the $n$-dimensional system % $${\bf y}' - {\bf A}(t) {\bf y} = {\bf g}(t), \eqno (19.21)$$ % where ${\bf A}(t)$ and ${\bf g}(t)$ are continuous and $2\pi$-periodic. % If the equation is noncritical as defined by $(19.20)$, % then there exists a unique $2\pi$-periodic solution~$\bf y$. % If it is critical, then there exist infinitely many % $2\pi$-periodic solutions if % $\int_0^{2\pi} {\bf Y}^{-1}(s)\kern .5pt {\bf g}(s)\kern .5pt ds$ % is in the range of the matrix % ${\bf I}-({\bf Y}(2\pi))^{-1}$, and no $2\pi$-periodic solutions if % it is not in this range. % } % \smallskip % %% % % \def\qed{\vrule height4pt width3pt depth4pt} % {\em Proof.} If (19.21) is noncritical, we may % multiply (19.19) on the left by $({\bf I}-{\bf Y}(2\pi))^{-1}$ to get % $${\bf c} = ({\bf I} - {\bf Y}(2\pi))^{-1}\kern 1pt {\bf Y}(2\pi) % \int_0^{2\pi} {\bf Y}^{-1}(s)\kern .5pt {\bf g}(s)\kern .5pt ds,$$ % that is, % $${\bf c} = (({\bf Y}(2\pi))^{-1}-{\bf I})^{-1} % \int_0^{2\pi} {\bf Y}^{-1}(s)\kern .5pt {\bf g}(s)\kern .5pt ds,$$ % and this vector corresponds to a $2\pi$-periodic solution % (19.18) as claimed. Uniqueness follows since % a given $\bf c$ corresponds to one and only one solution % ${\bf y}(t) = {\bf Y}(t) {\bf c}$. % If the equation is critical, then (19.19) is a singular system % of equations, having infinitely many solutions if the % right-hand side is in the range of % ${\bf I}-{\bf Y}(2\pi)$ and no solutions if it is not in that range. For the % right-hand side of (19.19) to be in the range of ${\bf I}-{\bf Y}(2\pi)$ is the % same as for $\int_0^{2\pi} {\bf Y}^{-1}(s)\kern .5pt {\bf g}(s)\kern .5pt ds$ % to be in the range of $({\bf Y}(2\pi))^{-1} % - {\bf I}$. \qed % \smallskip % %% % % Having written down these results for a first-order system % of ODE\kern .5pt s, let us immediately specialize them to the % case that has had the most attention. % {\bf Hill's equation} is the scalar, homogeneous, second-order ODE % $$y'' + f(t)y = 0 , \eqno (19.22)$$ % where $f$ is real and $2\pi$-periodic.\footnote{In the literature % the period is often taken to be $\pi$ instead of $2\pi$.} % Although one could introduce a nonzero right-hand side, % the equation is conventionally considered in this homogeneous form. % Note that there is no first-order term, so there is no damping. % {\bf Mathieu's equation} is the special case % $$y'' + (a + b\cos(t))y = 0 , \eqno (19.23)$$ % which is the simplest instance of Hill's equation in which % the periodic coefficient is not simply a constant. It corresponds to % taking the first nonconstant term in a Fourier series of % $f$ and is thus a natural first step in the study of (19.22). % Again the equation is usually considered in this homogeneous form. % %% % % As always, (19.22) is equivalent to a first-order system of % two variables via the identification % ${\bf y} = (y_1^{},y_2^{})^T = (y,y')^T$, which transforms it to % the equations % $$y_1' = y_2^{}, \quad y_2' + f(t) y_1^{} = 0. \eqno (19.24)$$ % For this system of ODE\kern .5pt s, % ${\bf Y}(t)$ takes the form % {\bf Y}(t) = \pmatrix{u(t) & v(t) \cr\noalign{\vskip 3pt} % u'(t) & v'(t)} , % where $u$ and $v$ are the unique solutions of the IVPs % $$u'' + f(t)u = 0 , \quad u(0) = 1, ~ u'(0) = 0$$ % and % $$v'' + f(t)v = 0 , \quad v(0) = 0, ~ v'(0) = 1.$$ % If we write the vector % $\bf c$ as ${\bf c} = (a,b)^T$, the solution (19.18) (here % with ${\bf g}={\bf 0}$) becomes $y(t)=au(t) + bv(t)$. % The solution will be $2\pi$-periodic if and only if % $y(2\pi) = y(0)$ and $y'(2\pi) = y'(0)$, that is, % {\bf D} \pmatrix{a\cr\noalign{\vskip 3pt}b} = % \pmatrix{a\cr\noalign{\vskip 3pt}b} , % with ${\bf D}$ defined by % {\bf D} = % \pmatrix{u(2\pi) & v(2\pi) \cr\noalign{\vskip 3pt} u'(2\pi) & v'(2\pi)} . % \eqno (19.25) % The criticality condition will be % $${\bf D} \hbox{ has an eigenvalue equal to\/ } 1. \eqno (19.26)$$ % By the usual reasoning we see that if ${\bf D}$ has % no eigenvalue equal to $1$, then (19.22) has a unique periodic solution, % namely $y=0$, corresponding to $a=b=0$, % whereas if ${\bf D}$ has an eigenvalue equal to $1$, then % (19.22) has infinitely many periodic solutions, corresponding to a % nonzero eigenvector $(a,b)^T$. % %% % We summarize these results in a theorem. %% % % \smallskip % {\em % {\bf Theorem 19.3. Hill's equation % (\textsf{\textbf{fLaSHi}}).} % Consider the Hill equation % $$y'' + f(t)y = 0 , \eqno (19.27)$$ % where $f(t)$ is real, continuous, and $2\pi$-periodic. % If the equation is noncritical as defined by $(19.25)$--$\kern 1pt (19.26)$, % then there exists a unique $2\pi$-periodic solution $y$, namely % $y(t) = 0$. % If it is critical, then there exist infinitely many % $2\pi$-periodic solutions. % } % \smallskip % %% % % The equations (19.22)--(19.27) are homogeneous, meaning % that in the generic case, the only solution % is the zero solution. One reason such equations are % interesting nonetheless is that we are effectively speaking here % of eigenvalue problems. In the Mathieu equation % (19.23), for example, $-a$ plays the % role of an eigenvalue of the operator $y \mapsto y'' + b\cos(t) % \kern .5pt y$ with periodic boundary conditions, % For $b=0$, the nonzero eigenvalues are $-1, -4, -9, \dots,$ and each % is a double eigenvalue because of the translational symmetry. % For nonzero $b$ the symmetry is broken and the eigenvalues become % simple. This has dynamical consequences, a hint of which appears % in Exercise 19.10. % %% % % We close this chapter with a few more remarks concerning % periodic solutions of autonomous nonlinear ODE\kern .5pt s, which are % $T$-periodic for any $T$. % For a second-order scalar problem of this type, or a system of % two first-order equations, a phase plane analysis may % quickly give an understanding of what periodic solutions may exist. % For example, the phase plane portrait of the damped linear % oscillator (8.4) makes it clear that $y=0$ will be the % unique $T$-periodic solution of this equation for any $T$. For % another example, the phase plane of the nonlinear pendulum % equation (9.7), $y'' =-\sin(y),$ explains why this problem % has nonzero $T$-periodic solutions for all $T$ in the range % $(2\pi, \infty)$ (Exercise 19.5). % Similarly, the Lotka--Volterra system (10.3) % has nonzero $T$-periodic solutions for a range of % values of $T$ (Exercise 19.8). The van der Pol oscillator (1.2), % with its unique limit cycle, has nonzero $T$-periodic % solutions only when $T$ takes exactly the right value, as mentioned % below (19.3) (or an integer multiple of this value). % %% % % Suppose we have an autonomous problem like the van der Pol % oscillator with a limit cycle, corresponding to a nonzero % $T$-periodic solution for a particular value of $T$. % As mentioned earlier, this solution will not be unique, since % $t$ can be shifted by any amount $\Delta t$, but one can make % it unique, or make its multiplicity finite rather than infinite, % by imposing one or more {\bf anchor conditions.} % For example, here is the van der Pol equation with % the anchor condition $y(0) = 1$ and an unknown period $T$: % $$0.3y''-(1-y^2)y'+y = 0, ~~ y(0) = 1, ~~ % y \hbox{ T-periodic}. \eqno (19.28)$$ % From a phase diagram such as Fig.~9.11, it is % clear that there will be two solutions to this problem (for % minimal $T$), depending % on where in the oscillation the anchor condition is reached, and % the two solutions will have $y'(0)>0$ and $y'(0) < 0$. % %% % % Although (19.28) has the form of a $T$-periodic equation for % unknown $T$, it can be converted to a $2\pi$-periodic problem % if we transform the independent variable to % $s = 2\pi t/T\in [\kern .3pt 0,2\pi]$: % $$0.3{T^2\over 4\pi^2}{d^{\kern .7pt 2} y\over d s^2} - % {T\over 2\pi}(1-y^2){dy\over ds} +y = 0, ~~ % y ~ 2\pi\hbox{-periodic with } y(0) = 1. \eqno (19.29)$$ % This is a $2\pi$-periodic problem with an unknown parameter, but % a trick that is sometimes used is to think of $T$ not as a number % but as a solution to the trivial differential equation % $${dT\over ds} = 0, ~~ % T \hbox{ is } 2\pi\hbox{-periodic}. \eqno (19.30)$$ % The unknown parameter has become an unknown constant value of % one component of a solution of a system of ODE\kern .5pt s, namely % (19.29)--(19.30)! This formulation % is not really simpler in any genuine way, but it has the % advantage that one can now attempt to solve the problem simply % by applying a software tool % for solving nonlinear ODE\kern .5pt s. % %% % % On p.~\pageref{question3} the question (3) was raised of nonperiodic % solutions of periodic ODE\kern .5pt s. This is a very important subject which is % dealt with by {\em Floquet theory.} For a scalar problem, the main % result asserts that every solution can be written as a complex % exponential $\exp(iax)$ for some $a$, real or complex, times a % periodic function, a pattern readily seen in Fig.~19.4. % For systems of equations $a$ becomes a matrix $A$. % We shall not discuss this further except in the following application. % %% % % \begin{center} % \hrulefill\\[1pt] % {\sc Application: band gaps and forbidden frequencies}\\[-3pt] % \hrulefill % \end{center} % %% % % The material world is made of atomic and molecular structures % that are often periodic, at least on % a microscopic scale, and electromagnetic % waves propagate through these structures in distinctive ways. % The resulting field of X-ray crystallography has been the basic tool % by which structure of all kinds of molecules are determined, % as recognized by Nobel Prizes to von Laue % (1914), Bragg and Bragg (1915), Watson and Crick (1962), % Hodgkin (1964), and others.\footnote{William % and Lawrence Bragg were father and son, and when the % prize was announced, 25-year-old Lawrence was serving in the trenches % in World War I.\ \ Incredibly, he remained in the trenches for % a further year after winning the award. % Nobel Prize winners get better treatment nowadays.} With % the advent of quantum mechanics it was seen that electrons % too are associated with waves, whose distinctive modes of propagation % though periodic structures determine whether a material is % an insulator, a conductor, or a semiconductor. % %% % % Here we shall give an indication of a fundamental % phenomenon of wave propagation in crystals that was investigated % by Felix Bloch, another Nobel Prize winner, in the % beginning of the quantum era. % Let us return to the Schr\"odinger equation as given in % (6.11), now (just for convenience) % with periodic boundary conditions, % $$-y''/2 + V(x) y = \lambda y, \quad x\in [\kern .3pt 0,d\kern .5pt ], ~ % y(0) = y(d\kern .5pt ) ,~ % y'(0) = y'(d\kern .5pt ) .~ \eqno (19.31)$$ % The important new feature is that we will take the % potential function $V(x)$ to be periodic, consisting % of a sequence of spikes with spacing $1$, % $$V(x) = 60 (\cos(\pi x))^{16}. \eqno (19.32)$$ % (This is a variant of the {\em Kronig--Penney % model\/} of a 1D crystal.) % To keep the example simple we set $d=8$, so $V$ extends % over just 8 periods. % Here is a plot of the potential, in black, and the % first 33 eigenfunctions, each one raised up a distance equal % to its eigenvalue as in Chapter~6. % d = 8; L = chebop(0,d); L.bc = 'periodic'; V = chebfun('60*cos(pi*x)^16',[0 d],'trig'); L.op = @(x,y) -diff(y,2)/2 + V*y; neigs = 33; [W,D] = eigs(L,neigs); W = simplify(W,1e-4); e = sort(diag(D)); for k = 1:neigs, W(k) = e(k)+3*W(k); end, plot(W,LW,1), hold on xlabel('$x$',FS,10,IN,LT), ylabel('eigenfunctions',FS,9,IN,LT) plot(V,'k',LW,1.4), hold off, ylim([0 100]) title('Fig.~19.7.~~Periodic potential and eigenfunctions',FS,11) %% % \vskip 1.02em %% % The striking thing about this plot is that the % lower eigenfunctions are grouped into bands. At the bottom, hard % to distinguish, are eight eigenfunctions with eigenvalue $\lambda\approx 5$. % Then there are eight eigenfunctions less sharply focussed with % $\lambda \approx 20$ and another eight broadly distributed with % $\lambda \approx 45$. The final 9 eigenfunctions shown correspond to % values $\lambda$ greater than the maximum value of the % potential $V(x)$. Here is a plot. plot(e,'.k',MS,8), ylim([0 100]) title('Fig.~19.8.~~Eigenvalues',FS,11) xlabel('$j$',FS,10), ylabel('eigenvalue $j$',FS,9) %% % \vskip 1.02em %% % % \noindent % The eigenvalues % fall in groups of 8. If we had taken thousand or millions % of periods, these would have fallen % so close together as to produce effectively a continuum % of eigenvalues, more or less like this: % for j = 0:3 plot(j*d+[1 d],e(j*d+[1 d]),'-k',LW,2), hold on end title('Fig.~19.9.~~Schematic for a larger crystal',FS,11) xlabel('$j$',FS,10), ylabel('eigenvalue $j$',FS,9), hold off set(gca,XTL,{' ',' ',' ',' ',' ',' ',' ',' '}) %% % \vskip 1.02em %% % % The implications of this kind of structure % for physics are enormous. A periodic medium may have % {\em forbidden frequencies\/} at which no wave can propagate. % To a condensed matter physicist, these represent an % {\em energy gap\/} or {\em band gap\/} in the spectrum, and the gaps % determine how electrons can propagate. % Roughly speaking, a material is an insulator if it has an energy % gap, a conductor if there is no gap, and a semiconductor if there % is a gap but it is very small. % One place to learn more is in the classic textbook % {\em Introduction to Solid State Physics} % by Kittel, first published in 1966, where it is pointed out % that the difference in % conductivity between an insulator and a conductor (a conductor % of the usual sort, that is, not a superconductor) may % be a factor as high as $10^{32}$. % %% % % Where do band gaps come from? To understand this % we can begin by recalling what the eigenfunctions of the % Schr\"odinger problem would look like if the potential % function $V(x)$ were constant rather than periodic. They would % simply be sines and cosines, $\sin(kx)$ and $\cos(kx)$, % for each value of the wave number~$k$. In a structure with % period $L$, however, interference effects occur when % $k$ is close to an integer multiple of $\pi /L$, corresponding % to a wavelength close to one of the values $2L, 2L/2, 2L/3,\dots .$ % The interfering waves may lead to reflection rather % than propagation, and in X-ray crystallography this is known as % as {\em Bragg reflection.} % %% % % To see the effect concretely, in our computed example, % let us return to Fig.~19.7 and % zoom in on the highest eigenfunction of the first % band, eigenfunction $8$, and the lowest eigenfunction of % the second band, eigenfunction $9$. % plot(W(:,8:9),LW,2), hold on for k = 8:9 w = W{k}; h1 = plot(w,LW,1.6); hold on h2 = plot([0 8],mean(w)*[1 1],'--',LW,0.7); set(h2,CO,get(h1,CO)) end xlabel('$x$',FS,10,IN,LT), ylabel('eigenfunctions',FS,9,IN,LT) plot(V,'k',LW,0.7), hold off, ylim([0 24]) title(['Fig.~19.10.~~Eigenfunctions 8 and 9, ' ... 'straddling the first band gap'],FS,11) %% % \vskip 1.02em %% % % \noindent % Both eigenfunctions have the same wavelength, namely 2, % and if the potential function were constant they would just % be translates of one another corresponding to exactly the same % eigenvalue. However, the periodic potential % has broken the symmetry between the sines and the cosines. % Eigenfunction $8$ is concentrated between the spikes and is little % affected by them, with an eigenvalue not so different from what it % would be if the potential were zero. Eigenfunction $9$, on the other % hand, is concentrated within the spikes, making the potential weigh % heavily so that the eigenvalue is much higher. % %% % % \smallskip % {\sc History.} The study of periodic ODE\kern .5pt s dates to the % 1880s with the work of Hill, Floquet, and Mathieu. % George William Hill (1838--1914), unusually for % mathematicians of that era, was an American, and he was % concerned with astronomical calculations related to the % three-body and four-body problems. % The general theory is due to % Gaston Floquet (1847--1920) in Sur les % \'equations diff\'erentielles lin\'eaires \`a coefficients % p\'eriodiques'' (you don't have to speak French to % translate this one\kern .3pt !), which appeared in {\em Annals of the % \'Ecole Normale Sup\'erieure,} v.~12, pp.~47--88, 1883. % Physicists often give credit to Bloch. % %% % % \smallskip % {\sc Our favorite reference.} % An extensive treatment of % periodic problems can be found in R. Grimshaw, {\em Nonlinear Ordinary % Differential Equations,} CRC Press, 1991. % This book includes a chapter on linear periodic problems, % a chapter on periodic solutions of autonomous nonlinear problems, % and further chapters on more advanced related topics. % %% % % \begin{displaymath} % \framebox[4.7in][c]{\parbox{4.5in}{\vspace{2pt}\sl % {\sc Summary of Chapter 19.} If the coefficients of % an ODE are\/ $T$-periodic for some\/ $T$, then there % may or may not exist\/ $T$-periodic solutions. % If the ODE is linear, there is an associated criticality condition. % If the equation is noncritical, it has a unique solution, whereas % if it is critical, it has no solutions or infinitely many solutions % depending on the inhomogeneous forcing data. % \vspace{2pt}}} % \end{displaymath} % %% % \smallskip\small\parskip=1pt\parindent=0pt % {\em \underline{Exercise $19.1$}. Periodic solutions via shooting.} % Suppose we seek a $2\pi$-periodic solution to % $y' = -y^3 + \sin(t)$. % {\em (a)} Produce a plot showing the solutions emanating % from initial values $y(0) = -1.5,-1.3,\dots, 1.3,1.5$. % Based on this plot, estimate by eye the value of $y(0)$ that % gives a periodic solution. % {\em (\kern .7pt b)} Prove that if $y(0)\in [-1,1]$, then $y(2\pi)\in [-1,1]$. % {\em (c)} Look up the Brouwer fixed point theorem and state it. % Use it to prove that this ODE has a periodic solution. % {\em (d)} Compute this solution with Chebfun using % the \verb|'periodic'| flag. You will need to provide an initial % guess with roughly the form of the true solution. % What is $y(0)\kern .7pt?$ % \par % {\em \underline{Exercise $19.2$}. An unstable variant.} % Now consider the equation of Exercise 19.1 with the other sign, % $w' = w^3 + \sin(t)$. % {\em (a)} Explain why shooting is much more difficult in this % case with reference to a computer plot showing solutions emanating % from several initial values $w(0)$. Instead of the whole interval % $[\kern .3pt 0,2\pi]$ you may use a shorter interval such as % $[\kern .3pt 0,0.5]$ or $[\kern .3pt 0,1]$. % {\em (\kern .7pt b)} On the other hand, the plots suggest an easy way % to reduce this problem to that of Exercise 19.1. Explain this % reduction and use it to derive an analytic formula for $w$ in terms % of the solution $y$ of Exercise 19.1. % \par % {\em Exercise $19.3$. Criticality depends on the interval.} % {\em (a)} Explain how the theorems of this chapter allow for % the possibility that a $2\pi$-periodic linear ODE might have % a unique $2\pi$-periodic solution but not a unique % $4\pi$-periodic solution. % {\em (\kern .7pt b)} Find an example of this behavior. % \par % {\em \underline{Exercise $19.4$}. Two electrons and a nucleus.} % As in Exercise 13.5, consider the idealized problem of two electrons % of mass 1 and charge $-1$ orbiting a nucleus of % charge $+2$ fixed at the origin of the $x$-$y$ plane. Let $z(t)$ % be the position of one electron represented with the % usual complex variable % $z = x+iy$, and suppose the configuration is symmetric % so that the other electron is at $\overline z(t)$. % {\em (a)} Write down the ODE governing the evolution of $z(t)$ % assuming an inverse-square electrostatic force law with constant $1$. % {\em (\kern .7pt b)} Plot the trajectories $z(t)$, $t\in [\kern .3pt 0,40\kern .3pt]$ % corresponding to the initial position $z(0) = ia$ and $z'(0) = 1$ % for $a = 1$ and $2$. % {\em (c)} Find a value of $a$ (to 3 digits of accuracy or more) % that gives a periodic solution. (Use % periodic boundary conditions in Chebfun or not, as you prefer.) % What is the period? % \par % {\em Exercise $19.5$. Nonlinear pendulum.} % Consider $T$-periodic solutions of the nonlinear pendulum % equation $y'' = -\sin(y)$ satisfying $y(0) = 0$. Explain why there % is one such solution for $T\in (0,2\pi]$, three such solutions % for $T\in (2\pi, 4\pi]$, five such solutions % for $T\in (4\pi, 6\pi]$, and so on. % \par % {\em \underline{Exercise $19.6$}. Logistic equation with periodic harvesting.} % As a generalization of Exercise 3.15, consider the equation % $y' = (1-y/Y)y - \sin(t)^2$, where $Y$ is a positive constant. % {\em (a)} Solve the equation with Chebfun with $Y=5$ for $t\in [\kern .3pt 0,15]$ % and make a plot of trajectories % from initial values $y(0) = 0.5,1,\dots, 8$. % {\em (\kern .7pt b)} Use Chebfun to find the oscillatory periodic solution % that the curves are approaching. Use the {\tt mean} command to find % the mean value of this solution. What would the mean be (figure this % part out analytically) if $\sin(t)^2$ were replaced by its average % value $1/2\kern .5pt$? % {\em (c)} Find and plot another periodic solution to this equation. % Again compare the mean to what it would be % if $\sin(t)^2$ were replaced by $1/2$. % \par % {\em \underline{Exercise $19.7$}. Periodic caffeine intake.} % Figure 2.8 showed the caffeine concentration in the bloodstream % for a drinker enjoying three cups of coffee in a 24-hour period. Suppose the % drinker has the same three cups at the same times every day. Compute the % periodic solution, taking the time interval % to be $[-2,22\kern .3pt]$. How does % the maximum caffeine concentration for the periodic solution compare % to that in Fig.~2.8? % \par % {\em \underline{Exercise $19.8$}. Lotka--Volterra system.} % {\em (a)} In Chapter 10 the period $T$ of the Lotka--Volterra system (10.6) % was computed for initial populations $(u_0^{},v_0^{})$ = $(1,1)$ of rabbits % and foxes, respectively. Now let $u_0^{}$ vary from $0.2$ to $4$ and % plot the dependence of $T$ on $u_0^{}$. % {\em (\kern .7pt b)} For $u_0\to 0.2$, determine $T$ analytically. % \par % {\em \underline{Exercise $19.9$}. Tokieda's teacup.} % Place a teacup on a table with the handle positioned in a direction % we shall call north. Practice tapping gently on the rim with a spoon % so that you hear a clean tone. Note that if you tap in the N, E, S, % or W positions, you hear one tone, whereas if you tap in the NE, SE, SW, % or NW positions, you hear a higher tone, often about one semitone higher (a factor % of $2^{1/2}$). (There is a video of Tadashi Tokieda demonstrating % this effect at \verb|https://youtu.be/MfzNJE4CK_s|.) % To explain this, we can imagine that the cup is a ring % oscillating in certain eigenmodes. {\em (a)} First, to model a simple % ring with no handle to break the symmetry, compute the % first five eigenvalues and eigenfunctions % of the differential operator $y''$ with % periodic boundary conditions. Note that eigenvalues $2n$ and $2n+1$ are % equal for each $n\ge 1$ because of rotational symmetry. % {\em (\kern .7pt b)} Now add the handle --- an extra mass at one point % along the ring --- by considering the differential operator $y''/m(x)$, % where $m(x) = 1 + 0.6\exp(-20(x-\pi)^2)$. Again compute the first % five eigenvalues and eigenfunctions. Because of the added mass, all % the frequencies will now be lower. Because of the broken symmetry, % there will no longer be any degeneracies. It is modes 4 and 5 that we % mainly hear when tapping the cup, % corresponding to the ring alternately getting taller/thinner and % shorter/fatter. Measure the ratio of frequency 5 to frequency 4 and % discuss the result. % Plot the eigenfunctions and explain why they are positioned where % they are. {\em (c)} What connection do you see with the Application of % this chapter? (The next exercise is also related.) % \par % {\em \underline{Exercise $19.10$}. Mathieu equation.} % Consider the Mathieu eigenvalue problem % $y'' + b\cos(t)\kern .3pt y = -\lambda y$ % with periodic boundary conditions on $[\kern .3pt 0,2\pi]$. % For the rest of this exercise fix $b = 1$. % {\em (a)} Determine the smallest six eigenvalues $\lambda_1^{} < \cdots % < \lambda_6^{}$ % and plot the associated eigenfunctions. Comment on the relationships % with the case $b=0$. % {\em (\kern .7pt b)} Now solve % $y'' + (a+b\cos(t))y = \exp(\sin(t))$ for $t\in[\kern .3pt 0, % 100\kern .3pt \pi]$ with initial data $y(0)=y'(0) = 0$ % for four choices of $a$: % $(\lambda_2^{}+\lambda_3^{})/2,\dots,$ % $(\lambda_5^{}+\lambda_6^{})/2;$ use {\tt L.maxnorm = 0}. % Plot the results and comment on the % differences between these curves. %