\documentclass[twoside]{article} \usepackage{amsfonts, amsmath, graphicx} \pagestyle{myheadings} \markboth{\hfil An adaptive numerical method for the wave equation \hfil EJDE/Conf/10} {EJDE/Conf/10 \hfil A. S. Ackleh, K. Deng, \& J. Derouen \hfil} \begin{document} \setcounter{page}{23} \title{\vspace{-1in}\parbox{\linewidth}{\footnotesize\noindent Fifth Mississippi State Conference on Differential Equations and Computational Simulations, \newline Electronic Journal of Differential Equations, Conference 10, 2003, pp 23--31. \newline http://ejde.math.swt.edu or http://ejde.math.unt.edu \newline ftp ejde.math.swt.edu (login: ftp)} \vspace{\bigskipamount} \\ % An adaptive numerical method for the wave equation with a nonlinear boundary condition % \thanks{ {\em Mathematics Subject Classifications:} 35B40, 35L05, 35L20, 65M25, 65N50. \hfil\break\indent {\em Key words:} Time-adaptive numerical method, blow-up time, blow-up rate. \hfil\break\indent \copyright 2003 Southwest Texas State University. \hfil\break\indent Published February 28, 2003. } } \date{} \author{Azmy S. Ackleh, Keng Deng, \& Joel Derouen} \maketitle \begin{abstract} We develop an efficient numerical method for studying the existence and non-existence of global solutions to the initial-boundary value problem \begin{gather*} u_{tt}=u_{xx}\quad 00,\\ -u_{x}(0,t)=h(u(0,t)) \quad t>0,\\ u(x,0)=f(x),\quad u_{t}(x,0)=g(x) \quad 00,\\ -u_{x}(0,t)=h(u(0,t)) \quad t>0,\\ u(x,0)=f(x),\quad u_{t}(x,0)=g(x) \quad 00$. Moreover, if$h(u)$is Lipschitz continuous, then the solution is unique. \end{theorem} \begin{theorem} Suppose that$|h(u)|\leq\rho(|u|)$with$\rho(r)>0$continuous, nondecreasing on$[0,\infty)$, and such that $\int^{\infty}\frac{dr}{\rho(r)}=\infty,$ then all mild solutions of (\ref{1}) are global. \end{theorem} \begin{theorem} Suppose that$f(t)+\int_{0}^{t}g(s)ds\geq0\ (\not \equiv0)$on$[0,\infty)$and that$h(u)\geq\sigma(|u|)$with$\sigma(r)>0$continuous, nondecreasing on$[0,\infty)$, and such that $\int^{\infty}\frac{dr}{\sigma(r)}<\infty,$ then every mild solution of (\ref{1}) blows up in finite time. \end{theorem} \begin{theorem} Suppose that$\int_{0}^{\infty}f(t)dt+\int_{0}^{\infty}\int_{0}^{t}g(s)dsdt>0$and$h(u)\geq c|u|^{p}$($p>1$,$c>0$), then the mild solution of (\ref{1}) blows up in finite time. \end{theorem} In \cite{AD}, we point out that the blow-up occurs on the boundary$x=0$only. Moreover, using asymptotic techniques for integral equations \cite{BH} we establish the following blow-up rates: Letting$T_{b}$be the blow-up time, \begin{itemize} \item If$h(u)\sim u^{p}$, then$u(0,t)\sim\big(\frac{1}{p-1}\big) ^{\frac{1}{p-1}}(T_{b}-t)^{-\frac{1}{p-1}}$as$t\to T_{b}$; \item If$h(u)\sim e^{u}$, then$u(0,t)\sim\log\big( \frac{1}{T_{b}-t}\big)$as$t\to T_{b}$. \end{itemize} The goal of this paper is to develop a numerical method for solving (\ref{1}). In Section 2 we discuss the numerical approximation while in Section 3, we present numerical examples. In Section 4, we conclude with some remarks. \section{Time-Adaptive Method} We begin this section by integrating (\ref{1}) along characteristics to obtain the following integral representation of solutions: For$t\leq x$, $$u(x,t)=\frac{1}{2}[f(x+t)+f(x-t)]+\frac{1}{2}\int_{x-t}^{x+t}g(s)ds, \label{21}%$$ and for$t>x, \begin{aligned} u(x,t)=&\frac{1}{2}[f(t+x)+f(t-x)]+\frac{1}{2}\Big[ \int _{0}^{t+x}g(s)ds+\int_{0}^{t-x}g(s)ds\Big] \\ &+\int_{0}^{t-x}h(u(0,\tau))d\tau. \end{aligned} \label{22a} A solution to the integral equations (\ref{21})-(\ref{22a}) defines a mild solution to the problem (\ref{1}). Furthermore, if the initial dataf$and$g$are smooth and satisfy some compatibility conditions, then one can argue that a solution of (\ref{21})-(\ref{22a}) is also a strong solution of (\ref{1}). Our numerical method will focus on the approximation of (\ref{21}% )-(\ref{22a}) rather than (\ref{1}). This provides an efficient scheme which does not require a rather strong regularity assumption on the initial data. Substituting$x=0$in (\ref{22a}), we get the Volterra integral equation $$u(0,t)=f(t)+\int_{0}^{t}g(s)ds+\int_{0}^{t}h(u(0,\tau))d\tau. \label{23}$$ Since blow-up occurs only on the boundary$x=0$, a special attention will be devoted to the development of an approximation of$u(0,t)$particularly near the blow-up time$T_{b}$. Once this is achieved, the approximations of the blow-up time$T_{b}$and$u(0,t)$are used to compute$u(x,t)$from the equations (\ref{21})-(\ref{22a}). To this end, differentiating (\ref{23}) we get the following differential equation for$u(0,t)$: $\frac{du(0,t)}{dt}=\frac{df(t)}{dt}+g(t)+h(u(0,t)).$ Let$\Delta t>0$be sufficiently small. Using Taylor approximation (formally) we observe that $u(0,t+\Delta t)-u(0,t)=\Delta t\frac{du(0,t)}{dt}+\frac{d^{2}u(0,\xi)}{dt^{2}% }\Delta t^{2},\quad \xi\in(t,t+\Delta t).$ A key idea in our scheme is to adapt the time step in order to keep the quantity$|u(0,t+\Delta t)-u(0,t)|\sim|\Delta t\frac{du(0,t)}{dt}|$bounded by a fixed constant. Since$h(u)\to\infty$as$u\to\infty$and blow-up occurs at$T_{b}$we see that$ \frac{du(0,t)}{dt}\to\infty$, as$t\to T_{b}$. In particular, as$t\to T_{b}$the size of the time step must approach zero if the magnitude of$\Delta t\frac{du(0,t)}{dt}$is to remain bounded by a fixed constant. This forces the numerical approximation not to go beyond the blow-up time. Making use of this fact we now present a time-adaptive algorithm for computing$u(0,t)$and the blow-up time$T_{b}$. Let$\Delta t_{\min}$and$\Delta t_{\max}$be fixed numbers with$0<\Delta t_{\min}<\Delta t_{\max}<\infty$. Let$u_{0}^{i}$be the approximation of$u(0,t_{i})$with$t_{0}=0$and$\Delta t_{i}=t_{i}-t_{i-1}\in\left[ \Delta t_{\min},\Delta t_{\max}\right]$. Denote by $(u_{t})_{0}^{i}=\frac{u_{0}^{i}-u_{0}^{i-1}}{\Delta t_{i}}$ the difference approximations of$u_{t}(0,t_{i})$. Guess an initial time step$\Delta t_{1}$and fix a scaling factor$\alpha>1$. Choose constants$d_{l}$and$d_{u}$such that$d_{l}d_{u}$then the approximated quantity$|u_{0}% ^{i+1}-u_{0}^{i}|>d_{u}$. In this case the time step is decreased by a factor of$1/\alpha$and the solution is recomputed at the new time step$(1/\alpha)\Delta t_{i}$. The second case is that if the current time step$\Delta t_{i}<\Delta t_{\max}$,$|u_{0}^{i+1}-u_{0}^{i}|\leq d_{l}$and$|u_{0}^{i}-u_{0}^{i-1}|\leq d_{l}$, then this indicates that the time steps used for the last two iterations are very conservative. Hence, the scheme increases this time step to$\min(\alpha\Delta t_{i},\Delta t_{\max})$in order to save computation time. \ It is easy to see that near the blow-up time, the time step$\Delta t_{i}$will decrease until it reaches$\Delta t_{\min}$. When this happens the computation stops, and the current time is an approximation of the blow-up time$T_{b}$. We remark that the accuracy of the approximations of$T_{b}$depends on the choice of$\Delta t_{\min}$. To compute$u_{0}^{i}$we combine the Runge-Kutta numerical method (see for example, \cite{FB}) with the above time-adaptive algorithm: Let$u_{0}% ^{0}=f\left( 0\right) and% \begin{align*} k_{1} & =\Delta t_{i+1}y\big( t_{i},u_{0}^{i}\big) \\ k_{2} & =\Delta t_{i+1}y\big( t_{i}+\frac{\Delta t_{i+1}}{2},u_{0}% ^{i}+\frac{1}{2}k_{1}\big) \\ k_{3} & =\Delta t_{i+1}y\big( t_{i}+\frac{\Delta t_{i+1}}{2},u_{0}% ^{i}+\frac{1}{2}k_{2}\big) \\ k_{4} & =\Delta t_{i+1}y\big( t_{i+1},u_{0}^{i}+k_{3}\big) ,\, \end{align*} wherei=0,1,2,\dots$, and$\Delta t_{i+1}$is determined by the time-adaptive method developed above. Compute$u_{0}^{i+1}$as follows: $u_{0}^{i+1}=u_{0}^{i}+\frac{1}{6}\left( k_{1}+2k_{2}+2k_{3}+k_{4}\right).$ Now, to approximate the solution of (\ref{21})-(\ref{22a}) we choose$x_{\max }>0$and divide the interval$[0,x_{\max}]$into uniform mesh$x_{j}$with$\Delta x=x_{j}-x_{j-1}$,$j=0,1,\dots,m$. Denote by$S^{n}(a,b,I)$the Simpson's numerical method for integrating a function$I(t)$on the interval$(a,b)$with$n$subdivisions, and let$P^{h}(t)$be the cubic interpolant of the function$h(u(0,t))$at the mesh points$t_{i}$. Then we let$u_{j}^{i}$be the approximation of$u(x_{j},t_{i})$and compute$u_{j}^{i}$as follows: For$t_{i}\leq x_{j}$, $u_{j}^{i}=\frac{1}{2}[f(x_{j}+t_{i})+f(x_{j}-t_{i})]+\frac{1}{2}S^{n}% (x_{j}-t_{i},x_{j}+t_{i},g),$ and for$t_{i}>x_{j}, \begin{align*} u_{j}^{i} =&\frac{1}{2}[f(t_{i}+x_{j})+f(t_{i}-x_{j})]\\ & +\frac{1}{2}\left[ S^{n}(0,t_{i}+x_{j},g)+S^{n}(0,t_{i}-x_{j},g)\right] +S^{n}(0,t_{i}-x_{j},P^{h}). \end{align*} In the next section we present numerical results which indicate the accuracy of such an adaptive numerical scheme in computing bothu(x,t)$and the blow-up time$T_{b}$. \section{Numerical Results} The numerical method developed in the previous section is now used to corroborate and complement theoretical results in our earlier paper \cite{AD}. For the rest of this section let$\Delta t_{\max}=10^{-3}$,$\Delta t_{\min }=10^{-7}$,$\alpha=2$,$d_u = 1$,$d_l=0.1,n=10,x_{\max}=5$, and$m=200$. In the first example we present the accuracy of our method. To this end, we choose$f=0$,$g=1$and$h(u)=u^{2}$. It is not difficult to show that$u(0,t)=\tan t,$and hence blow-up occurs at$t=\pi/2$. In Figure 1 we show the relative error$\dfrac{|u_{0}^{i}-\tan t_{i}|}{\Delta t_{i}}$. The computed blow-up time$T_{b}=1.5704$. \begin{figure} [ptb] \begin{center} \includegraphics[width=0.7\textwidth]{fig1.eps} \caption{The relative error between the computed function$u(0,t)$and the exact solution$\tan t$.} \end{center} \end{figure} In the second example we let$f(x)=-(x-2)^{2}+4$,$g(x)=0$and$h(u)=u^{3}$. Notice that this choice of initial data does not satisfy the assumptions of Theorems 1.3-1.4 in Section 1. However, the numerical results presented in Figures 2-3 indicate that blow-up occurs for this choice of functions with an approximated blow-up time$T_{b}=0.5118$. \begin{figure}[ptb] \begin{center} \includegraphics[width=0.7\textwidth]{fig2.eps}\caption{The computed function$u(0,t)$for the data$f(x)=-(x-2)^{2}+4$,$g(x)=0$and$h(u)=u^{3}$.} \end{center} \end{figure} \begin{figure} [ptbptb] \begin{center} \includegraphics[width=0.7\textwidth]{fig3.eps} \caption{The solution$u(x,t)$for the data$f(x)=-(x-2)^{2}+4$,$g(x)=0$and$h(u)=u^{3}$.} \end{center} \end{figure} In our third numerical experiment we examine whether blow-up occurs for nonlinearities such as$h(u)=(1+u)[\log(1+u)]^{p}$with initial data that do not satisfy the assumptions of Theorem 1.4. In Figure 4 we present the numerical results of$u(0,t)$for the case$p=6$,$f(x)=3e^{-x}\cos(20x)-0.1$and$g(x)=0$, and in Figure 5 we display the 3-D graph of the function$u(x,t)$. We remark that the blow-up time is$T_{b}=0.22296$. \begin{figure} [ptb] \begin{center} \includegraphics[width=0.7\textwidth]{fig4.eps} \caption{The computed function$u(0,t)$for the data$f(x)=3e^{-x}% \cos(20x)-0.1$and$g(x)=0$and$h(u)=(1+u)[\log(1+u)]^{6}$.}% \end{center} \end{figure} \begin{figure} [ptb] \begin{center} \includegraphics[width=0.7\textwidth]{fig5.eps} \caption{The computed solution$u(x,t)$for the data$f(x)=3e^{-x}\cos(20x)-0.1$and$g(x)=0$and$h(u)=(1+u)[\log(1+u)]^{6}$.}% \end{center} \end{figure} Using our numerical scheme, we have successfully verified the blow-up rates given in Section 1 for the functions$e^{u}$and$u^{p}\,(p>1)$. We now use this method to examine the blow-up rate for the function$h(u)=(1+u)[\log (1+u)]^{p}$. Before presenting the numerical results we formally derive such a rate. Near the blow-up time the values$\frac{df(t)}{dt}$and$g(t)$are negligible when compared to$u(0,t)$, and hence $\frac{du(0,t)}{dt}\sim(1+u(0,t))\left[ \log(1+u(0,t))\right] ^{p}.$ Integrating the above we find $\int_{u(0,t)}^{\infty}\frac{du}{(1+u)\left[ \log(1+u)\right] ^{p}% }\sim\int_{t}^{T_{b}}dt.$ Solving for$u$we get $$u(0,t)\sim e^{(\frac{1}{(p-1)(T_{b}-t)})^{\frac{1}{p-1}}}-1.\label{31}%$$ In Table 1 we give numerical results that verify such a blow-up rate. For this computational purpose we use the following equivalent form of (\ref{31}) $\frac{1}{p-1}=(T_{b}-t)[\log(1+u(0,t))]^{p-1}.$ \begin{table}[ht] \caption{The blow-up rate for the function$h(u)=(1+u)(\log (1+u))^{p}$.} \begin{center} \begin{tabular} [c]{|l|l|l|l|l|}\hline \multicolumn{1}{|c|}{$p$} & 4 & 6 & 8 & 10\\\hline Conjectured:$\frac{1}{p-1}$&$0.3333$&$0.2$&$0.1429$&$0.1111$\\\hline Approximation &$0.3205$&$0.1973$&$0.1411$&$0.1106$\\\hline \end{tabular} \end{center} \end{table} \section{Concluding Remarks} The objective of this paper is to develop a numerical approximation for studying the existence and non-existence of global solutions to the wave equation with a nonlinear boundary condition. Our numerical results indicate that such a scheme is very accurate and efficient for computing the blow-up time, the blow-up rate, and the solution. These results also open up several theoretical questions: 1) How much can the conditions on the initial data$f$and$g$be relaxed for blow-up to occur? 2) Can one improve Theorem 1.4 for weaker nonlinearties such as$h(u)=(1+u)[\log(1+u)]^{p}$($p>1)\$? 3) Can one prove the blow-up rate given by (\ref{31}) for such nonlinearities? Our future research efforts will focus on such questions as well as the application of time-adaptive methods to a system of wave equations coupled in the boundary conditions discussed in \cite{AD2}. Finally, it is worth mentioning that one can also devise a numerical method by directly approximating the Volterra integral equation (\ref{23}) using a combination of the time-adaptive method presented here and numerical quadrature methods for Volterra integral equations \cite{BM}. \begin{thebibliography}{0} \frenchspacing \bibitem{AD} A. S. Ackleh and K. Deng, Existence and nonexistence of global solutions of the wave equation with a nonlinear boundary condition,\emph{\ Quarterly of Applied Mathematics,} \textbf{59 }(2001), 153-158. \bibitem{AD2} A. S. Ackleh and K. Deng, Global existence and blow-up for a system of wave equations coupled in the boundary conditions, \emph{Dynamics of Continuous, Discrete and Impulsive Systems,} {\bf 8} (2001), 415-424. \bibitem{BM} C. T. H. Baker and G. F. Miller, Treatment of Integral Equations by Numerical Methods, Acedemic Press, New York, 1982. \bibitem{BH} N. Bleistein and R. A. Handelsman, Asymptotic Expansion of Integrals, Holt, Rinehart and Winston, New York, 1975. \bibitem{FB} J. D. Faires and R. L. Burden, Numerical Methods, PWS-Publishing Company, Boston, 1993. \end{thebibliography} \noindent\textsc{Azmy S. Ackleh} (e-mail: ackleh@louisiana.edu)\\ \textsc{Keng Deng} (e-mail: deng@louisiana.edu)\\ \textsc{Joel Derouen} (e-mail: jbd8438@louisiana.edu)\\[2pt] Department of Mathematics\\ University of Louisiana at Lafayette\\ Lafayette, Louisiana 70504, USA. \end{document}