\documentclass[reqno]{amsart} \usepackage{hyperref} \usepackage{graphicx} \AtBeginDocument{{\noindent\small Eighth Mississippi State - UAB Conference on Differential Equations and Computational Simulations. {\em Electronic Journal of Differential Equations}, Conf. 19 (2010), pp. 75--83.\newline ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu \newline ftp ejde.math.txstate.edu} \thanks{\copyright 2010 Texas State University - San Marcos.} \vspace{9mm}} \begin{document} \setcounter{page}{75} \title[\hfilneg EJDE-2010/Conf/19/\hfil Numerical solution] {Numerical solution for nonlocal Sobolev-type differential equations} \author[S. A. Dubey\hfil EJDE/Conf/19 \hfilneg] {Shruti A. Dubey} \address{Shruti A. Dubey \newline Laboratoire d'Analyse et d'Architecture des Systemes du CNRS (LAAS-CNRS),\newline Toulouse, France} \email{shrurkd@gmail.com} \thanks{Published September 25, 2010.} \subjclass[2000]{35L65, 65M06, 65M12} \keywords{Sobolev-type differential equations; Laplace transform; \hfill\break\indent numerical inversion; nonlocal conditions} \begin{abstract} We present a numerical approximate solution to Sobolev-type differential equation subject to nonlocal initial boundary conditions. A Laplace transform method is described for the solution of considered equation. Following Laplace transform of the original problem, an appropriate method of solving differential equations is used to solve the resultant time-independent modified equation and solution is inverted numerically back into the time domain. Numerical results are provided to show the accuracy of the proposed method. \end{abstract} \maketitle \numberwithin{equation}{section} \newtheorem{theorem}{Theorem}[section] \newtheorem{lemma}[theorem]{Lemma} \newtheorem{example}[theorem]{Example} \allowdisplaybreaks \section{Introduction} This work is concerned with the Sobolev-type partial differential equation \begin{equation} \frac{\partial}{\partial t}\Big(u(x,t)-\frac{\partial^2}{\partial x^2}u(x,t)\Big)-\frac{\partial ^2}{\partial x^2}u(x,t)=f(x,t),\quad x\in [0,1],\;t>0,\label{eq1} \end{equation} subject to the conditions \begin{equation} \begin{gathered} u(x,0)=u_0(x),\quad x\in [0,1], \\ u(0,t)=\int^1_0 a(x)u(x,t)dx + p(t),\quad t>0,\\ u(1,t)=\int^1_0 b(x)u(x,t)dx + q(t),\quad t>0. \end{gathered} \label{eq2} \end{equation} where, $f$, $u_0$, $a,\;b,\;p,\;q$ are prescribed continuous functions and $u(x,t)$ is an unknown function which is a solution of \eqref{eq1} and satisfies conditions \eqref{eq2} at the same time. Sobolev-type equation appears in a variety of physical problems such as flow of fluid through fissured rocks, thermodynamics and propagation of long waves of small amplitude. There is an extensive literature in which Sobolev type of equations are investigated, in the abstract framework, see for instance \cite{ab,bk,kt,lr,sh2}. Some other models of nonlocal boundary conditions are numerically solved by Dehghan \cite{d1,d2}. The first results on Sobolev-type equation were obtained by Hilbert space methods \cite{ss}. Subsequently several authors studied and discussed same type of problem subject to classical as well as nonlocal conditions. Very strong and complete results are known concerning existence,uniqueness and properties of solutions. Brill \cite{b} and Showalter \cite{sh1} established the existence of solutions of semilinear evolution equations of Sobolev type in Banach space. Balachandran and Park \cite{bp} investigated integrodifferential equation of Sobolev type with nonlocal condition and proved the existence of mild and strong solutions using semigroup theory and Schauder fixed point theorem. Recently, the existence of solution to semilinear Sobolev type equation with integral conditions is studied by Bouziani and Merazga \cite{bm} using Rothe time discretization method. So far not much seems to be done for obtaining an explicit solution of Sobolev type partial differential equations, however the solvability of these equations have been theoretically studied in terms of the existence and uniqueness of a solution. The purpose of present article is to give a method of solution to Equation \eqref{eq1}, \eqref{eq2} using Laplace transform technique. In recent years, Laplace transform method has been used to approximate the solution of different class of linear partial differential equations \cite{a,bs,kw,ms}. Suying et al \cite{smzw}, established a numerical method based on Laplace transform for solving initial problem of nonlinear dynamic differential equations. The main difficulty in using Laplace transform method consists in finding its inverse. Numerical inversion methods are then used to overcome this difficulty. There are many numerical techniques available in literature to invert Laplace transform. In this paper we focus exclusively on the Stehfest inversion algorithm \cite{s} in order to efficiently and accurately invert the Laplace transform (which cannot be done analytically). The plan of the paper is as follows. In section 2, we develop a method of solution and find out a solution of the problem in Laplace domain. Obtained solution is then converted in to real domain by the means of numerical inversion algorithm. Some examples are given in Section 3 to demonstrate the competence of the method, followed by conclusion in Section 4. \section{Method based on Laplace transform} \subsection{Laplace transform technique} Laplace transform is widely used in the area of engineering technology and mathematical science. There are many problems whose solution may be found in terms of a Laplace transform. In fact, it is an efficient method for solving various differential equations. Laplace transform of a function $w(x,t)$ with respect to $t$ can be expressed as $$ W(x;s)=\mathcal{L}\{w(x,t);t\to s\} =\int_0^{\infty} w(x,t)e^{-st}dt, $$ where, $s$ is known as Laplace variable. A capital letter $W$ represents Laplace transform of a function $w$, i.e., $W$ is a function in Laplace domain. Starting with this approach and taking Laplace transform on both sides of \eqref{eq1}, \eqref{eq2} with respect to $t$, we get \begin{gather*} \frac{d^2}{dx^2}U(x;s)-\frac{s}{(1+s)}U(x;s)=-\frac{1}{(1+s)} \Big[F(x;s)+U(x;0)-\frac{d^2}{dx^2}U(x;0)\Big], \\ U(x;0)=u_0(x), \\ U(0;s)=\int^1_0a(x)U(x;s)dx+P(s), \\ U(1;s)=\int^1_0b(x)U(x;s)dx+Q(s). \end{gather*} Using the initial condition into the differential equation, we obtain \begin{equation} \begin{gathered} \frac{d^2}{dx^2}U(x;s)-\frac{s}{(1+s)}U(x;s)=-\frac{1}{(1+s)} \Big[F(x;s)+u_0(x)-\frac{d^2}{dx^2}u_0(x)\Big], \\ U(0;s)=\int^1_0a(x)U(x;s)dx+P(s), \\ U(1;s)=\int^1_0b(x)U(x;s)dx+Q(s), \end{gathered} \label{eq3} \end{equation} where, $U(x;s)=\mathcal{L}\{u(x,t);t\to s\}$, $F(x;s)=\mathcal{L}\{f(x,t);t\to s\}$, $P(s)=\mathcal{L}\{p(t);t\to s\}$ and $Q(s)=\mathcal{L}\{q(t);t\to s\}$. Thus, considered equation is reduced in boundary value problem governed by a second order inhomogeneous ordinary differential equation. On solving it, we obtain a general solution of \eqref{eq3} as \begin{equation} \begin{aligned} U(x;s)&=\frac{-1}{(1+s)}\sqrt{\frac{(s+1)}{s}}\int _0^x[F(\tau;s)+u_0(\tau)-\frac{d^2}{d\tau^2}u_0(\tau)]\\ &\quad\times \sinh\Big(\sqrt{\frac{s}{s+1}}(x-\tau)\Big)d\tau +C_1(s)e^{-\sqrt{\frac{s}{s+1}}x} +C_2(s)e^{\sqrt{\frac{s}{s+1}}x}, \end{aligned}\label{eq4} \end{equation} where, $C_1$ and $C_2$ are arbitrary functions of $s$. Substituting \eqref{eq4} in to the boundary conditions, we have \begin{align*} &C_1(s)\Big[1-\int_0^1a(x)e^{-\sqrt{\frac{s}{s+1}}x}dx\Big] +C_2(s)\Big[1-\int_0^1a(x)e^{\sqrt{\frac{s}{s+1}}x}dx\Big]\\ &= \frac{-1}{(s+1)}\sqrt{\frac{s+1}{s}}\int_0^1 \Big[[F(\tau;s)+u_0(\tau)-\frac{d^2}{d\tau^2}u_0(\tau)] \\ &\quad\times \int_\tau^1a(x)\sinh\Big(\sqrt{\frac{s}{s+1}}(x-\tau) \Big)dx\Big]d\tau+P(s), \end{align*} \begin{align*} &C_1(s)\Big[e^{-\sqrt{\frac{s}{s+1}}}-\int_0^1b(x) e^{-\sqrt{\frac{s}{s+1}}x}dx\Big] +C_2(s)\Big[e^{\sqrt{\frac{s}{s+1}}} -\int_0^1b(x)e^{\sqrt{\frac{s}{s+1}}x}dx\Big]\\ &=\frac{-1}{(s+1)}\sqrt{\frac{s+1}{s}}\int_0^1 \Big[[F(\tau;s)+u_0(\tau)-\frac{d^2}{d\tau^2}u_0(\tau)]\\ &\quad\times \int_\tau^1b(x)\sinh \Big(\sqrt{\frac{s}{s+1}}(x-\tau)\Big)dx\Big]d\tau+Q(s)\\ &\quad+\frac{1}{(s+1)}\sqrt{\frac{s+1}{s}} \int_0^1[F(\tau;s)+u_0(\tau)-\frac{d^2}{d\tau^2}u_0(\tau)]\sinh \Big(\sqrt{\frac{s}{s+1}}(1-\tau)\Big)d\tau. \end{align*} It gives \begin{equation} \begin{pmatrix} C_1(s) \\ C_2(s) \end{pmatrix} =\begin{pmatrix} a_{11}(s)\quad a_{12}(s)\\ a_{21}(s)\quad a_{22}(s) \end{pmatrix}^{-1} \begin{pmatrix} b_1(s)\\b_2(s) \end{pmatrix},\label{eq5} \end{equation} where, \begin{equation} \begin{gathered} \label{eq6} a_{11}(s)= 1-\int_0^1a(x)e^{-\sqrt{\frac{s}{s+1}}x}dx, \quad a_{12}(s)=1-\int_0^1a(x)e^{\sqrt{\frac{s}{s+1}}x}dx, \\ a_{21}(s)=e^{-\sqrt{\frac{s}{s+1}}} -\int_0^1b(x)e^{-\sqrt{\frac{s}{s+1}}x}dx, \quad a_{22}(s)=e^{\sqrt{\frac{s}{s+1}}} -\int_0^1b(x)e^{\sqrt{\frac{s}{s+1}}x}dx,\\ \begin{aligned} b_1(s)&=-\frac{1}{(s+1)}\sqrt{\frac{s+1}{s}}\int_0^1 \Big[[F(\tau;s)+u_0(\tau)\\ &\quad -\frac{d^2}{d\tau^2}u_0(\tau) \int_\tau^1a(x)\sinh\Big(\sqrt{\frac{s}{s+1}}(x-\tau)\Big)dx\Big] d\tau+P(s), \end{aligned} \\ \begin{aligned} b_2(s)&=\frac{-1}{(s+1)}\sqrt{\frac{s+1}{s}}\int_0^1 \Big[[F(\tau;s)+u_0(\tau)-\frac{d^2}{d\tau^2}u_0(\tau)] \\ &\quad \times \int_\tau^1b(x) \sinh\Big(\sqrt{\frac{s}{s+1}}(x-\tau)\Big)dx\Big]d\tau+Q(s) +\frac{1}{(s+1)}\sqrt{\frac{s+1}{s}}\\ &\quad \times \int_0^1[F(\tau;s) +u_0(\tau)-\frac{d^2}{d\tau^2}u_0(\tau)] \sinh\Big(\sqrt{\frac{s}{s+1}}(1-\tau)\Big)d\tau . \end{aligned} \end{gathered} \end{equation} Thus, to find out the solution in Laplace domain one has to evaluate all the integrals appear in \eqref{eq4} and \eqref{eq6}. This can be done for known functions $F,P,Q,a,b,u_0$, however, in many cases, the resulting function is not easy to integrate exactly. Therefore, there is need for numerical approximation of the integrals. A well known Gaussian Quadrature formula exists for computing integrals numerically (see \cite{as}). Using this formula we have the following approximations of the above integrals: \begin{gather*} \int_0^1\begin{pmatrix}a(x)\\b(x) \end{pmatrix} e^{\pm \sqrt{\frac{s}{s+1}}x}\,dx \simeq \frac{1}{2}\sum_{i=1}^{n}w_i \begin{pmatrix} a(\frac{1}{2}(x_i+1))\\ b(\frac{1}{2}(x_i+1)) \end{pmatrix} )e^{\pm \sqrt{\frac{s}{s+1}}(\frac{1}{2}(x_i+1))},\\ \begin{aligned} &\int_0^x\Big[F(\tau;s)+u_0(\tau)-\frac{d^2}{d\tau^2}u_0(\tau)\Big] \sinh\Big(\sqrt{\frac{s}{s+1}}(x-\tau)\Big)d\tau \\ &\simeq \frac{x}{2}\sum_{i=1}^nw_i\Big[F\left(\frac{x}{2} (x_i+1);s\right)+u_0\left(\frac{x}{2}(x_i+1)\right) -\frac{d^2}{d\tilde{\tau}^2}u_0\left(\frac{x}{2}(x_i+1)\right)\Big]\\ &\quad \times \sinh\Big(\sqrt{\frac{s}{s+1}} \left(x-\frac{x}{2}(x_i+1)\right)\Big), \end{aligned} \end{gather*} where $\tilde{\tau}=\frac{x}{2}(x_i+1)$. \begin{align*} &\int_0^1\Big[\Big(F(\tau;s)+u_0(\tau)-\frac{d^2}{d\tau^2}u_0(\tau) \Big)\int_\tau^1 \begin{pmatrix}a(x)\\b(x) \end{pmatrix} \sinh\Big(\sqrt{\frac{s}{s+1}}(x-\tau)\Big)dx\Big]d\tau \\ &\simeq \frac{1}{2}\sum_{i=1}^nw_i \Big[F\left(\frac{1}{2}(x_i+1);s\right)+u_0\left(\frac{1}{2}(x_i+1) \right)-\frac{d^2}{d\tilde{\tau}^2}u_0\left(\frac{1}{2}(x_i+1)\right) \Big]\\ &\quad \times\left(\frac{1-\frac{1}{2}(x_i+1)}{2}\right) \sum_{j=1}^{n}w_j \begin{pmatrix} a\big(\frac{1-\frac{1}{2}(x_i+1)}{2}x_j +\frac{1+\frac{1}{2}(x_i+1)}{2}\big)\\ b\big(\frac{1-\frac{1}{2}(x_i+1)}{2}x_j +\frac{1+\frac{1}{2}(x_i+1)}{2}\big)\\ \end{pmatrix}\\ &\quad\times \sinh\Big(\sqrt{\frac{s}{s+1}}\Big( \big(\frac{1-\frac{1}{2}(x_i+1)}{2}\big)x_j +\frac{1+\frac{1}{2}(x_i+1)}{2}-\frac{1}{2}(x_i+1)\Big)\Big), \end{align*}\\ where, $\tilde{\tau}=\frac{1}{2}(x_i+1)$. $x_i$ and $w_i$ are defined by $$ x_i: i^{th} \text{ zero of } P_n(x),\quad w_i=2/(1-x_i^2)[P'_n(x_i)]^2 $$ and known as abscissas and weights respectively. Their tabulated values can be found in \cite{as} for different values of $n$. \subsection{Numerical inversion of Laplace transform} Now, we have a solution in Laplace transform domain as given in \eqref{eq4}. So we expect to obtain a solution of original problem by means of inverting the Laplace transform. Simple transforms can often be inverted using readily available table. More complex functions can be analytically inverted through the complex inversion formula $$ g(t)=\frac{1}{2\pi j}\int_{c-j\infty}^{c+j\infty} e^{st}G(s)ds, $$ where, $c$ is a positive real number such that all the poles of the function $G(s)$ lie at the left of the line $\operatorname{Re}(s)=c$. Sometimes analytical inversion of a Laplace domain solution is difficult to obtain, therefore, a numerical inversion method must be used. A variety of different methods for numerically inverting the Laplace transform are available that can be employed. There exists no universal method but different types of methods work well for different classes of functions. A nice comparison of four frequently used numerical Laplace inversion algorithms is given by Hasan et al. \cite{hd}. We use the Stehfest algorithm \cite{s} in the work as it is easy to implement and leads to result of sufficient accuracy throughout the time range. This numerical technique was first introduced by Graver \cite{g} and its algorithm then offered by Stehfest. Stehfest's algorithm approximates the time domain solution as $$ u(x,t)\approx \frac{ \ln 2}{t}\sum_{k=1}^{2m}\beta_kU\Big(x; \frac{\ln 2}{t}k\Big), $$ where, $m$ is the positive integer and $$ \beta_k=(-1)^{m+k}\sum_{l=\left[\frac{k+1}{2}\right]}^{\min(k,m)} \frac{l^{m}(2l)!}{(m-l)!\;l!\;(l-1)!(k-l)!(2l-k)!}. $$ Here $[r]$ denotes the integer part of $r$. The parameter $m$ is a free parameter that should be optimized by trial and error. It was seen that with increasing $m$ accuracy of result increases up to a point and then owing to the rounding errors it decreases \cite{s}. Thus, for choosing optimum $m$, it is beneficial to apply an algorithm repeatedly for different values of $m$ and study its effect on the solution. The other way to choose optimal value of $m$ could be, to apply the Stehfest's algorithm for inverting the Laplace transform of some elementry functions which are known. \noindent\textbf{Remarks:} * Stehfest's method gives accurate results for many problems including diffusion problem, fractional functions in the Laplace domain. However, it fails to predict $e^t$ type functions or those with oscillatory behavior such as sine and wave functions (see \cite{hd}). * Note that more than one numerical inversion algorithm can also be performed to check the accuracy of the result. \section{Numerical results} This section illustrates some of the numerical examples that we carried out in order to test the reliability of the method proposed. We consider the problem \eqref{eq1}, \eqref{eq2} with different choice of functions as follows: \begin{example} \label{exa1} \rm Here, we take \begin{gather*} f(x,t)=\Big(x(x-1)-\frac{1}{7}\Big)e^{-t},\quad u_0(x)=x(1-x)+\frac{1}{7},\\ a(x)=b(x)=6/13,\quad p(t)=q(t)=0. \end{gather*} In this case exact solution is $$ u(x,t)=\Big(x(1-x)+\frac{1}{7}\Big)e^{-t}. $$ We use developed method for finding the numerical solution. It is found that the best choice of the parameter $m$ for this and subsequent example is 5. Thus we take m=5, n=8 and follow exactly the same steps described in the previous section for the method of solution. Obtained numerical result is compared with exact solution in Fig. \ref{fig1} and Fig. \ref{fig2}. \end{example} Fig. \ref{fig1} provides the contrast of the numerical result and the exact solution for $t=0.5 $ and $x\in [0,1]$. In Fig. \ref{fig2}, the comparison between these solutions are given for a fixed value of $x=0.2 $ and for $t\in[0.1,1]$. Results show that in both cases, approximate solution agrees well with the exact solution. \begin{figure}[ht] \begin{center} \includegraphics[width=0.7\textwidth]{fig1} % graph1_x.eps \end{center} \caption{Comparison between numerical and exact solutions for $t=0.5$, $x\in[0,1]$} \label{fig1} \end{figure} \begin{figure}[ht] \begin{center} \includegraphics[width=0.7\textwidth]{fig2} % graph1_t.eps \end{center} \caption{Comparison between numerical and exact solutions for $x=0.2$, $t\in[0.1, 1]$} \label{fig2} \end{figure} \begin{example} \label{exa2} \rm Here, we take \begin{gather*} f(x,t)=\frac{-2(x^2+t+1)}{(t+3)^3},\quad u_0(x)=x^2/9, \quad a(x)=b(x)=x,\\ p(t)=-\frac{1}{4(t+3)^2},\quad q(t)=\frac{3}{4(t+3)^2}. \end{gather*} Its exact solution is $$ \frac{x^2}{(t+3)^2}. $$ \end{example} To apply the numerical method proposed, we need the Laplace transform of the above mentioned functions $f(x,t)$, $p(t)$, $q(t)$ with respect to $t$ which can be given by \begin{gather*} F(x;s)=\int_0^\infty \frac{-2s(sx^2+\xi+s)}{(\xi+3s)^3}e^{-\xi}d\xi,\\ P(s)=\int_0^\infty\frac{-s}{4(\xi+3s)^2}e^{-\xi}d\xi,\\ Q(s)=\int_0^\infty\frac{3s}{4(\xi+3s)^2}e^{-\xi}d\xi. \end{gather*} Integrals involved here are evaluated with the help of numerical approximation formula \cite[(25.45)]{as}. Proceeding in the similar manner as of the previous example and choosing $m=5$ and $n=8$, a comparison of numerical and exact solutions is done. In Fig. \ref{fig3} graphs of the approximate and exact solution are given for $t=0.4$ and $x\in [0,1]$. For $x=0.6$ and $t\in[0.1, 1]$ results are drawn in Fig. \ref{fig4}. Presented graphs clearly show that the approximate and exact solutions are almost superposed. \noindent\textbf{Remarks:} * The numerical solution matches with exact solution up to at least 4 significant places of decimal. * The simulation results show that our proposed method achieves good performance. \begin{figure}[ht] \begin{center} \includegraphics[width=0.7\textwidth]{fig3} % graph2_x.eps \end{center} \caption{Comparison between numerical and exact solutions for $t=0.4$, $x\in[0,1]$} \label{fig3} \end{figure} \begin{figure}[ht] \begin{center} \includegraphics[width=0.7\textwidth]{fig4} % graph2_t.eps \end{center} \caption{Comparison between numerical and exact solutions for $x=0.6$, $t\in[0.1,1]$} \label{fig4} \end{figure} \subsection{Conclusion} A method of solution based on Laplace transform to the considered nonlocal problem is described. The benefit of the presented method is that it gives explicitly a numerical approximate solution of the problem, however, several theoretical approaches to existence and uniqueness of solution to Sobolev-type equations are summarized in literature since last few decades. Numerical results show that the time domain solution evaluated by presented method is comparable with the exact solution. \begin{thebibliography}{99} \bibitem{as} Abramowitz M., Stegun I. A.; \emph{Handbook of Mathematical functions}, Dover, New York, 1972. \bibitem{ab} Agarwal S., Bahuguna D.; Existence of solutions to Sobolev-type partial neutral differential equations, {\em J. Appl. Math. Stoc. Anal.} 2006 (2006) DOI 10.1155/JAMSA/2006/16308. \bibitem{a} Ang W. T.; A method of solution for the one-dimensional heat equation subject to nonlocal conditions, {\em Southeast Asian Bulletin of Mathematics} 26 (2002) 185-191. \bibitem{bs} Babolian E., Shamloo A. S.; Numerical solution of Volterra integral and integro-differential equations of convolution type by using operational matrices of piecewise constant orthogonal functions, {\em J. Comp. Appl. Math.} 214 (2008), 495--508. \bibitem{bk} Balachandran Krishnan, Karunanithi Subbarayan; Regularity of solutions of Sobolev-type semilinear integrodifferential equations in Banach spaces, {\em Electron. J. Diff. Eq.} Vol. 2003, No. 114 (2003), 1-8. \bibitem{bp} Balachandran K., Park J. Y.; Sobolev type integrodifferential equation with nonlocal condition in Banach spaces, {\em Taiwanese Journal of Mathematics} 7 (2003) 155-163. \bibitem{bm} Bouziani Abdelfatah, Merazga Nabil; Solution to a semilinear pseudoparabolic problem with integral conditions, {\em Elec. J. Diff. Eq.} Vol. 2006 No. 115 (2006) 1-18. \bibitem{b} Brill H.; A semilinear Sobolev evolution equation in Banach space, {\em J. Diff. Eq.} 24 (1977) 412-425. \bibitem{d1} Dehghan M.; On the solution of an initial boundary value problem that combines Neumann and integral condition for the wave equation, {\em Numerical Methods for Partial Differantial Equations}, 21(1) (2005), 24--40. \bibitem{d2} Dehghan M.; The one-dimensional heat equation subject to a boundary integral specification, {\em Chaos, Solitons Fractals}, 32(2) (2007), 661--675. \bibitem{g} Gaver D. P.; Observing stochastic processes and approximate transform inversion, {\em Oper. Res} 14 (1966) 444-459. \bibitem{hd} Hassanzadeh Hassan; Pooladi-Darvish Mehran, Comparison of different numerical Laplace inversion methods for engineering applications, {\em Appl. Math. Comp.} 189 (2007) 1966-1981. \bibitem{kt} Kojo Tomomi, Tsutsumi Masayoshi; Initial value problem for abstract differential equations of Sobolev type, {\em 6th Int. Coll. on differential equations} (1996) 135-141. \bibitem{kw} Kong Weibin, Wu Xionghua; The Laplace transform and polynomial Trefftz method for a class of time dependent PDEs, {\em Appl. Math. Mod.} in press. \bibitem{lr} Lightbourne III J. H., Rankin III S. M.; A partial functional differential equation of Sobolev type, {\em J. Math. Anal. Appl.} 93 (1983) 328-337. \bibitem{ms} Maleknejad K., Shamloo A. S.; Numerical solution of singular Volterra integral equations system of convolution type by using operational matrices, {\em Appl. Math. Comp.} 195 (2008), 500--505. \bibitem{sh1} Showalter R. E.; Existence and representation theorem for a semilinear Sobolev equation in Banach space, {\em SIAM J. Math. Anal.} 3 (1972) 527-543. \bibitem{sh2} Showalter R. E.; A nonlinear parabolic Sobolev equation, {\em J. Math. Anal. Appl.} 50 (1975) 183-190. \bibitem{ss} Sobolev S.; Some new problems in mathematical physics, {\em Izv. Akad. Nauk SSSR Ser. Mat.} 18 (1954) 3-50. \bibitem{s} Sthefest H.; Numerical inversion of the Laplace transform, {\em Comm. ACM} 13(1) (1970) 47-49. \bibitem{smzw} Suying Zhang, Minzhen Zhang, Zichen Deng, Wencheng Li; Solution of nonlinear dynamic differential equations based on numerical Laplace transform inversion, {\em Appl.Math. Comp.} 189 (2007) 79-86. \end{thebibliography} \end{document}