%% 5. Boundary-value problems %% % % \setcounter{page}{48} % %% % % Now it is time to open the big door that will double the % richness of our study of ODE\kern .5pt s. % With second- and higher-order equations, since they need % more than one boundary condition, we have the % prospect of not just IVPs but also {\bf BVPs}, that is, % {\bf boundary-value problems}. % %% % % Boundary-value problems are everywhere in the sciences. % Broadly speaking, % ODE\kern .5pt s whose independent variables represent time % are usually IVPs modeling {\em dynamics.} % ODE\kern .5pt s whose independent variables represent % space are usually BVPs modeling % {\em equilibria.}\footnote{Of course, there are % exceptions. A notable area where BVPs arise with a time % variable is {\em control theory,} in which the engineer seeks to % make a system attain certain function values at certain times. The % Application of Chapter~16 is in this category.} % These differences are fundamental, and % in this book, we will usually change the independent variable from % $t$ to $x$ when working with a BVP.\ \ In % particular, % $y'$ and $y''$ will now denote % $dy/dx$ and $d^{\kern .7pt 2}y/dx^2$ instead of $dy/dt$ and $d^{\kern .7pt 2}y/dt^2$. % As mentioned in a footnote in Chapter 1, we % will also plot solutions of BVPs in blue instead of green. % %% % % As our first example of a BVP, here is the solution of % $$y''=-y,\quad x\in [\kern .3pt 0,60\kern.3pt], % ~y(0)=1,~y(60)=0.\eqno (5.1)$$ % The value of $y$ at $x=0$ is specified, but % not its derivative. A solution is found % to match the value of $y$ specified at the other end of the interval. % ODEformats, L = chebop(0,60); L.op = @(x,y) diff(y,2) + y; L.lbc = 1; L.rbc = 0; y = L\0; plot(y,CO,bvp), axis([0 60 -4 4]) title('Fig.~5.1.~~BVP (5.1) with Dirichlet BCs',FS,11) %% % \vskip 1.02em %% % % Here is the solution to another BVP, in which the condition at the % right involves the derivative $y'$ rather than the value of $y$ itself, % $$y''=-y,\quad x\in [\kern .3pt 0,60\kern.3pt],~y(0)=1,~y'(60)=0.\eqno (5.2)$$ % A condition involving $y$ is called a {\bf Dirichlet} condition, and % a condition involving $y'$ is called a {\bf Neumann} condition. A % condition that mixes % both $y$ and $y'$ is called a {\bf Robin} condition. % (To specify a Robin condition in Chebfun, such as % $y'(60) - y(60) = 1$, for example, we could write % \verb|L.rbc = @(y) diff(y)-y-1|.) % L.rbc = 'neumann'; y = L\0; plot(y,CO,bvp), axis([0 60 -4 4]) title('Fig.~5.2.~~BVP (5.2) with a Neumann BC',FS,11) %% % \vskip 1.02em %% % % The significance of boundary conditions becomes % more striking if we look at equation (4.10) % with $\omega^2 <0$, where one solution % exponentially grows and the other exponentially decays. % To illustrate, we revert to the Dirichlet condition of (5.1) % and take $\varepsilon =0$, $\omega^2=-1$ in (4.10), % $$y''=y,\quad x\in [\kern .3pt 0,60\kern.3pt],~ y(0)=y(60)=1.\eqno (5.3)$$ % Note that the solution shown below is close to zero except near the % endpoints, where it is far from % zero in two regions known as {\bf boundary layers.} % Here we begin to see strongly the difference between IVPs and BVPs. % This curve in principle might have % been the solution of an IVP, if the initial data had been % chosen just right, but such a solution would hardly be likely % to arise in practice. On the contrary, % from the curve it is clear that this % solution is probably being driven by conditions imposed at the % two ends. We shall consider boundary layers in detail % in Chapter~20. % L.op = @(x,y) diff(y,2) - y; L.lbc = 1; L.rbc = 1; y = L\0; plot(y,CO,bvp), axis([0 60 -1.5 1.5]) title('Fig.~5.3.~~Nonoscillatory BVP (5.3)',FS,11) %% % \vskip 1.02em %% % % For the problem above, the characteristic roots (4.11)--(4.12) % are $r_{\pm} = \pm 1$. The solution we have plotted can % be approximated very closely as % $$y(x) \approx e^{-x} + e^{x-60}. \eqno (5.4)$$ % See Exercise 5.1. % %% % Suppose we are given a second-order constant-coefficient scalar linear % homogeneous ODE BVP of the form (4.10). To solve it, following % Theorem 4.1, we would look for a function in the % vector space spanned by $e^{r_+x}$ and $e^{r_-x}$, or % $e^{rx}$ and $xe^{rx}$ in the critically damped case, that % satisfies the boundary conditions. % Just as in the last chapter, a small variation of this approach % applies to an inhomogeneous problem of the form % $$y'' + \varepsilon\kern .3pt y' + \omega^2\kern .3pt y = f(x), \eqno (5.5)$$ % where $f$ is a given forcing function. If we can find % a particular solution $y_{\rm p}^{}$ to (5.5), then the general % solution is obtained by adding elements from this two-dimensional % vector space. When $f$ is simple, a particular solution % can often be obtained by inspection or the method of undetermined % coefficients. %% % % Here is an example. Suppose we put the % linear function $f(x) = (x-20)/20$ on the right-hand side of (5.3), % $$y''=y+(x-20)/20,\quad x\in [\kern .3pt 0,60\kern.3pt],~ y(0)=y(60)=0.\eqno (5.6)$$ % A particular solution to this equation is % $y_{\rm p}^{}(x) = -(x-20)/20$, and from % here one can readily obtain a formula for the solution plotted below % (Exercise 5.2). % L.op = @(x,y) diff(y,2) - y; L.lbc = 0; L.rbc = 0; rhs = chebfun('(x-20)/20',[0 60]); y = L\rhs; plot(y,CO,bvp), axis([0 60 -2 2]) title('Fig.~5.4.~~Inhomogeneous BVP (5.6)',FS,11) %% % \vskip 1.02em %% % For a nonlinear variation, let us % change the term $y$ in (5.6) to $y^3$, % $$y''=y^3+(x-20)/20,\quad x\in % [\kern .3pt 0,60\kern.3pt],~ y(0)=y(60)=0.\eqno (5.7)$$ % The resulting curve is an interesting variation % on Figure 5.4. The zero of $y(x)$ near $x=20$, incidentally, % is not exactly at $x=20$, more like $19.99999999985$. L.op = @(x,y) diff(y,2) - y^3; y = L\rhs; plot(y,CO,bvpnl), axis([0 60 -2 2]) title('Fig.~5.5.~~Nonlinear variation (5.7)',FS,11) %% % \vskip 1.02em %% % % Returning to linear problems, let us explore some of the % variety of boundary conditions that may arise. % Here is an {\em advection-diffusion} equation: % $$0.02y'' + y' + y = 0 , \quad x\in [\kern .3pt 0,1]. % \eqno (5.8)$$ % First let us solve it with two Dirichlet boundary conditions % $y(0) = 0$, $y(1) = 1$. % L = chebop(0,1); L.op = @(x,y) 0.02*diff(y,2) + diff(y) + y; L.lbc = 0; L.rbc = 1; y = L\0; plot(y,CO,bvp) title('Fig.~5.6.~~Solution of (5.8) with two Dirichlet BCs',FS,11) %% % \vskip 1.02em %% % % \noindent % Next we change the right boundary condition to Neumann form, % $y'(1) = 1$. % L.rbc = @(u) diff(u)-1; y = L\0; plot(y,CO,bvp), ylim([-3 0]) title(['Fig.~5.7.~~Solution of (5.8) with Dirichlet ' ... 'and Neumann BCs'],FS,11) %% % \vskip 1.02em %% % % \noindent % Here we couple the two endpoints with the % conditions $y(0) = 1$, $y'(1) = -y'(0)$. % L.lbc = 1; L.rbc = []; L.bc = @(x,y) deriv(y,1) + deriv(y,0); y = L\0; plot(y,CO,bvp) title('Fig.~5.8.~~Solution of (5.8) with coupled BCs',FS,11) %% % \vskip 1.02em %% % % Not all boundary'' conditions must be based at % the boundaries. For example, we might seek a solution % with $y(0) = 1$ and $\int_0^1 y(x) dx = 0$, so that % $y$ has mean zero. A nonstandard condition % like this is sometimes called a {\em side condition.} % L.bc = @(x,y) sum(y); y = L\0; plot(y,CO,bvp) title(['Fig.~5.9.~~Solution of (5.8) with an ' ... 'integral condition'],FS,11) %% % \vskip 1.02em %% % % So far we have considered second-order equations. % With equations of third order % or higher, there is the possibility of mixing boundary and interior % conditions. For example, here is the solution to a third-order % problem with boundary conditions at three distinct % points, marked by dots in the plot for emphasis. % $$y''' + y = 1, \quad x\in [\kern .3pt 0,2], ~ % y(0) = 0, ~ y(1) = y(2) = 1. \eqno (5.9)$$ % L = chebop(0,2); L.op = @(x,y) diff(y,3) + y; L.lbc = 0; L.rbc = 1; L.bc = @(x,y) y(1)-1; y = L\1; plot(y,CO,bvp) hold on, plot(0:2,y(0:2),'.k',MS,14), hold off title(['Fig.~5.10.~~Solution of 3rd-order eq.\ (5.9) ' ... 'with interior condition'],FS,11) %% % \vskip 1.02em %% % % Throughout this discussion we have avoided mentioning % existence and uniqueness. This is % because these matters are not straightforward for % BVPs, even linear ones, because of the possible presence of % {\em eigenvalues.} This will be the subject of the next chapter. % For nonlinear BVPs there are further nonuniqueness % possibilities that are both interesting and scientifically % important, as we will explore in Chapters 16--18. As an example of % nonuniqueness in a nonlinear BVP, consider % the solution that Chebfun finds to the nonlinear equation % $${1\over 2} y'' + y^2 = 1, \quad x\in [-1,1], ~ y(-1) = y(1) = 0. % \eqno (5.10)$$ % N = chebop(-1,1); N.op = @(x,y) .5*diff(y,2) + y^2; N.lbc = 0; N.rbc = 0; y = N\1; plot(y,CO,bvpnl) title(['Fig.~5.11.~~One solution of the ' ... 'nonlinear equation (5.10)'],FS,11) %% % \vskip 1.02em %% % % \noindent % If we ask Chebfun to look for a solution starting from % the initial guess $y(x) = \cos(\pi x/2)$ rather than its % default initial guess $y(x) = 0$, however, a different % and equally valid solution is found, which we superimpose % on the same plot. % N.init = chebfun('cos(pi*x/2)'); y = N\1; hold on, plot(y,'r'), hold off title('Fig.~5.12.~~Another solution of (5.10)',FS,11) %% % \vskip 1.02em %% % % \noindent % Does (5.10) have exactly these two solutions, or % are there more? There is no general theorem that % gives an answer to such questions, but we shall % make some progress on them in Chapter 16. % %% % % \begin{center} % \hrulefill\\[1pt] % {\sc Application: beam theory and the strength of spaghetti}\\[-3pt] % \hrulefill % \end{center} % %% % % In physics and engineering an idealized {\em string\/} can sustain % a tension but has no stiffness --- no resistance % to bending --- and this leads to second-order ODE\kern .5pt s. This is where we would % start, say, in studying the sounds produced by a guitar. An idealized % {\em beam}, on the other hand, has stiffness and resists bending, even % though we still model it as a one-dimensional object. % This leads to fourth-order ODE\kern .5pt s, and the same % equations apply to the girders of the Eiffel tower as % to a piece of dried spaghetti. % %% % In the simplest setting, omitting % dimensional constants and assuming infinitesimal % deflections, the equation takes the form % $$y''''(x) = f(x) , \eqno (5.11)$$ % where $f(x)$ is the vertical force applied at position $x$ along a beam and % $y(x)$ is the vertical deflection at that point. % The derivatives of $y$ can be interpreted as follows: %% % % \par\vskip .7em\par % $y(x)$ and $y'(x)$: position and slope, % \par\vskip .3em\par % $y''(x)$: bending moment, whose square is energy density, % \par\vskip .3em\par % $y'''(x)$: vertical force, % \par\vskip .3em\par % $y''''(x)$: force density. % \par\vskip .7em\par % \noindent % In a time-dependent problem, a nonzero value of $y''''(x)$ would % be associated with vertical acceleration $-y''''(x)$, % but we are considering here the static problem. % Equation (5.11) represents the condition that the force caused % by the bending of the beam is balanced by an external force % $f(x)$, so there is no acceleration. % %% % % For example, suppose a beam of length $2$ and weight density $1$ is % clamped at $x=\pm 1$, % $$y'''' = -1, \quad y(-1) = y'(-1) = y(1) = y'(1) = 0. \eqno (5.12)$$ % An image shows that the beam sags a little in the middle (the % deflection is $1/24$). % We plot $y$ with a thicker line to suggest stiffness. % This and the next few calculations, incidentally, can be % done on paper (Exercise 5.4), though that would change if % complications were introduced such as variable coefficients. % L = chebop(-1,1); L.op = @(x,y) diff(y,4); L.lbc = [0;0]; L.rbc = [0;0]; y = L\-1; d1 = norm(y,inf); plot(y,CO,bvp,LW,3); axis([-1.2 1.2 -.5 .1]) title(['Fig.~5.13.~~Clamped beam (5.12) of length 2: deflection ' ... num2str(d1)],FS,11) %% % \vskip 1.02em %% % % Next, instead of clamping the beam at $x=1$, let us fix the position % there but not the slope. In such a case one has what is called a % {\em natural boundary condition,} % $y''(1) = 0$. The deflection of the beam at the % midpoint doubles, and the maximum deflection, which is now located % to the right of the midpoint, increases by a factor of about 2.079. % (The exact maximal deflection is $(39+ 55\sqrt{33}\kern .7pt )/4096.$) % L.rbc = @(y) [y; diff(y,2)]; y = L\-1; d2 = norm(y,inf); plot(y,CO,bvp,LW,3); axis([-1.2 1.2 -.5 .1]) title(['Fig.~5.14.~~Beam with one freely supported ' ... 'end: deflection ' num2str(d2)],FS,11) %% % \vskip 1.02em %% % % \noindent % Alternatively, we can let both ends be free. Symmetry is restored, % and the maximum deflection is in the middle again, with % value $5/24$. % L.lbc = @(y) [y; diff(y,2)]; y = L\-1; d3 = norm(y,inf); plot(y,CO,bvp,LW,3); axis([-1.2 1.2 -.5 .1]) title(['Fig.~5.15.~~Beam with two freely supported ' ... 'ends: deflection ' num2str(d3)],FS,11) %% % \vskip 1.02em %% % % This brings us to the spaghetti problem. Here in Oxford, we bought a % pack of spaghetti at a local newsagent. % With Scotch tape we taped % some pennies\footnote{Conveniently for our experiment, the British 2p coin % weighs exactly the same as two pennies.} % to the end of a piece, making a cantilever that can be modeled % as a massless beam clamped at $x=0$ with a % weight $w$ applied at $x=d$: % $$% y'''' = 0, \quad y(0) = y'(0) = y''(d) = 0, ~ y'''(d) = w. % \eqno (5.13)$$ % d = 2; L = chebop(0,d); L.op = @(x,y) diff(y,4); L.lbc = [0; 0]; L.rbc = @(y) [diff(y,2); diff(y,3)-0.1]; y = L\0; d4 = norm(y,inf); plot(y,CO,bvp,LW,3) penny = d+1i*(y(end)-.01)+.08*(chebfun('exp(1i*t)',[0, 2*pi])-1i); hold on, fill(real(penny),imag(penny),[.72 .45 .2]) axis equal, xlim([-0.2 d+0.2]), hold off title(['Fig.~5.16.~~Spaghetti cantilever (5.13): deflection ' ... num2str(d4)],FS,11) %% % \vskip 1.02em %% % % \noindent % We could have measured these deflections for various cantilever % lengths $d$ and numbers of pennies, but it seemed more fun to investigate, % for a given number of pennies, % what is the length of cantilever at which the % spaghetti breaks? % With up to four pennies no breaking occurs --- at 26 cm, the % piece of spaghetti is too short. % But for $5,6,\dots, 16$ pennies we found these breaking lengths in cm: % d = [17.1 14.6 12.9 8.8 9.0 7.5 6.7 6.6 5.2 5.3 4.7 4.4]; %% % % \noindent % Here is a plot of the data, showing how the breaking weight decreases % as the length $d$ increases. % xx = linspace(4.7,17.5); plot(xx,80./xx,'-r'), hold on plot(d,5:16,'.k',MS,12), hold off xlabel('length of spaghetti',FS,9) ylabel('number of pennies to break',FS,9) text(10.3,10.4,'inverse',CO,'r',FS,10) text(10.3,9,'linear fit',CO,'r',FS,10) title('Fig.~5.17.~~Breaking a spaghetti cantilever',FS,11) %% % \vskip 1.02em %% % % \noindent % As the red line shows, the relationship appears to be % inverse linear, and we explain this as follows. % A piece of spaghetti breaks % when its curvature $|y''(x)|$ exceeds a critical value for some $x$. % From calculation or analysis of (5.13), one finds that as % a function of $d$, $y$ scales as $wd^{\kern .7pt 3}$, where $w$ again % is the weight. This means that % $y'$ scales as $wd^{\kern .7pt 2}$, and $y''$ scales as $wd$. % If $wd$ is to remain bounded as % $d\to\infty$, we must have $w = O(d^{\kern .7pt -1})$. % %% % % Returning from spaghetti to mathematics for a moment, it is interesting to % note that the solutions of $u''''=0$ are % precisely the polynomials of degree $\le 3$. Thus the % study of idealized beams with no body forces $f(x)$ % is equivalent to the study of % cubic polynomials that satisfy various boundary conditions. % This connection is famous in numerical analysis, where % the classic wooden flexible rulers used for designing % shapes in the aircraft and shipbuilding industries evolved into % the mathematical % {\em splines} used all across computational science and engineering. % %% % % \smallskip % {\sc History.} Beam theory goes back to Euler and % Jacob Bernoulli in the mid-18th century. It was % later that the explicit study of boundary-value problems % became a field of mathematics in its own right with % the investigations of {\em Sturm--Liouville problems} in the % 1830s by Charles Sturm and Joseph Liouville. % French mathematics had a period of extraordinary % productivity following the French Revolution % and Napoleonic wars, with leading figures including Fourier % (born 1768), Poisson (1781), Cauchy (1789), % Sturm (1803), and Liouville (1809). % \smallskip % %% % % {\sc Our favorite reference.} A standard work among numerical % analysts is Ascher, Mattheij, and Russell, % {\em Numerical Solution of Boundary Value Problems for % Ordinary Differential Equations,} SIAM, 1995 (first published in % 1988). Its appealing concrete style is well established by % Section I.2 on p.~7, which begins, % In this section we have collected 22 instances of BVPs which arise in % a variety of application areas.'' % \smallskip % %% % % \begin{displaymath} % \framebox[4.7in][c]{\parbox{4.5in}{\vspace{2pt}\sl % {\sc Summary of Chapter 5.} % When an ODE is of order $d\ge 2$, it will normally have % a $d$-dimensional solution space, making it possible % for different boundary conditions to be specified at different points. % In the simplest case $d=2$, we get a % {\em two-point boundary value problem.} % Problems of elasticity and solid mechanics lead to % fourth-order equations. % Existence and uniqueness of solutions are % not as straightforward for BVPs, even linear ones, as for IVPs. % Nonlinear BVPs may have multiple solutions. % \vspace{2pt}}} % \end{displaymath} % %% % % \small\parindent=0pt\parskip=1pt % {\em Exercise $5.1$. BVPs and IVPs.}\footnote{This exercise is related to % the numerical method for solving ODE BVPs known as {\em shooting}, and % highlights some of the challenges of that method. We shall say more % about shooting in Chapter~16.} % {\em (a)} Determine analytically the exact solution % to (5.3). How much does it differ from the approximation (5.4)? % {\em (\kern .7pt b)} As remarked in the text, though a curve may be % obviously'' the solution to a BVP, in principle it could equally % have been the solution to an IVP.\ \ For this problem (5.3), what % initial condition on $y'(0)$ would be needed to achieve % the same solution? {\em (c)} If that initial condition is % multiplied by $0.99$, is the resulting value % $y(60)$ less than or greater than Avogadro's number? % \par % {\em Exercise $5.2$. Exact solution to $(5.6)$.} {\em (a)} Derive a formula % for the exact solution to (5.6). {\em (\kern .7pt b)} Do the same for the case % where the right-hand side is changed to $(x-20)^2/400$. % \par % {\em \underline{Exercise $5.3$}. Sixth-order BVP.} % The ODE % $y^{(6)} - a y' + y = 0$, $x\in [-4,4]$, where $a\in [\kern .3pt 0,2\kern .3pt]$ is % a parameter, has boundary conditions $y(-4) = y'(-4) = y''(-4)=0$, % $y(4) = 1$, $y'(4)=y''(4) = 0$. If $y(0) = 0$, what is $a$? % \par % {\em {Exercise $5.4$}. Exact beam formulas.} Find exact formulas % for the solutions presented in Figures {\em (a)} 5.13, % {\em (\kern .7pt b)} 5.15, and {\em (c)} 5.14. % The reason for this ordering is that the solutions of % Figures 5.13 and 5.15 are simpler because of symmetry. % \par % {\em {Exercise $5.5$}. Exact spaghetti formula.} Find an % exact formula for the solution plotted in Figure 5.16. % \par % {\em {Exercise $5.6$}. Cubic and quartic dependence.} % It is readily confirmed that on % an interval of length $d$, the deflections of Figures 5.13--5.15 % scale as $O(d^{\kern .7pt 4})$, whereas in Figure 5.16, as noted % in the text, it is $O(d^{\kern .7pt 3})$. % Thus for large enough $d$ the clamped beam will sag more than % the cantilever. Yet the clamped beam is % supported at {\em both\/} ends, so surely % it should bend less. Resolve this paradox. % \par % {\em \underline{Exercise $5.7$}. Blasius equation.} % A classic nonlinear BVP due to Blasius is % $2\kern .3pt y''' + y\kern .3pt y'' = 0$, % $y(0) = y'(0) = 0$, $y'(\infty) = 1$. % Solve this with Chebfun, replacing $\infty$ by 10. Plot the % result and report the values of $y''(0)$ and $y(10) - 10$, both of % which are good approximations to what one % would get on $[\kern .3pt 0,\infty)$. % (For more information see J. P. Boyd, The Blasius function % in the complex plane,'' {\em Experimental Mathematics} 8 (1999), % pp.~381--394.) % \par % {\em {Exercise $5.8$}. Eliminating inhomogeneous boundary conditions.} % Suppose we have a linear inhomogeneous BVP with inhomogeneous % boundary data. To be specific, let us say it is a % second-order problem $Ly = f$, $y(0) = \alpha,$ % $y(1) = \beta$. % {\em (a)} Let $y_1^{}$ and $y_2^{}$ be two linearly independent % solutions to the homogeneous problem and let $y_{\rm left}^{}$, % $y_{\rm inner}^{}$, and $y_{\rm right}^{}$ be particular solutions % with $f=\beta=0$, $\alpha=\beta=0$, and $f=\alpha= 0$, % respectively. Write expressions for a particular solution of % the BVP and for the general solution. % {\em (\kern .7pt b)} A different approach is as follows. Let % $Y$ be {\em any\/} smooth function that satisfies $Y(0) = \alpha$ and % $Y(1) = \beta$; $Y$ does not have to satisfy the differential equation. % Show how the general solution of the BVP can be derived using $Y$. % \par % {\em \underline{Exercise $5.9$}. Water droplet.} The height of % the surface of a water droplet satisfies $y'' = (y-1)(1+(y')^2)^{1.5}$ % with $y(\pm 1) = 0$. What is the height at the midpoint $x=0$? % \par % {\em \underline{Exercise $5.10$}. Nonlinear problems % and Newton iteration.} Chebfun solves linear BVPs by solving % linear systems of equations obtained by discretizing in $x$, as % described in Appendix~A.\ \ For % nonlinear BVPs, there is an additional Newton iteration % involved, discussed in Appendix~A and also briefly in Chapter~16. % {\em (a)} Rerun the problem of equation (5.3) and Figure 5.3 using % {\tt tic} and {\tt toc} to measure how long the computation takes. % (As with all timing experiments, do this several times to make % sure the numbers have settled down.) Plot the solution $y$ on % a semilog scale and also report the length of the chebfun % as in Exercise~1.6. % {\em (\kern .7pt b)} Now do the same calculations again but with the right-hand % side of (5.3) changed from $y$ to $y^{1.01}$, % and superimpose the new solution curve % $y$ on the same plot. How is the computing time affected? How is % the length of the chebfun affected? % {\em (c)} Show the outputs that result when {\tt L} is executed % without a semicolon for these two computations. Note the indication % that in part {\em (\kern .7pt b)}, Chebfun has determined that the operator % is nonlinear and must be treated by different numerical methods. %