# Solving Parabolic and Hyperbolic Equations with Variable Coefficients Using Space-Time Localized Radial Basis Function Collocation Method

Authors: Mohammed Hamaidi; Ahmed Naji; Ahmed Taik
Journal: Modelling and Simulation in Engineering
 (2021)

---

## Abstract

In this paper, we investigate the numerical approximation solution of parabolic and hyperbolic equations with variable coefficients and different boundary conditions using the space-time localized collocation method based on the radial basis function. The method is based on transforming the originald-dimensional problem in space into d+1-dimensional one in the space-time domain by combining the d-dimensional vector space variable and 1-dimensional time variable in one d+1-dimensional variable vector. The advantages of such formulation are (i) time discretization as implicit, explicit, θ-method, method-of-line approach, and others are not applied; (ii) the time stability analysis is not discussed; and (iii) recomputation of the resulting matrix at each time level as done for other methods for solving partial differential equations (PDEs) with variable coefficients is avoided and the matrix is computed once. Two different formulations of the d-dimensional problem as a d+1-dimensional space-time one are discussed based on the type of PDEs considered. The localized radial basis function meshless method is applied to seek for the numerical solution. Different examples in two and three-dimensional space are solved to show the accuracy of such method. Different types of boundary conditions, Neumann and Dirichlet, are also considered for parabolic and hyperbolic equations to show the sensibility of the method in respect to boundary conditions. A comparison to the fourth-order Runge-Kutta method is also investigated.


---

## Body

## 1. Introduction

Second-order parabolic and hyperbolic partial differential equations with variable coefficients are one of the most important problems in Mathematical and engineering fields. They attract the attention of a large number of scientists. The general heat transfer problem is considered one of the more interesting examples. Such initial boundary value problems (parabolic or hyperbolic) can be solved by coupling time-stepping algorithms with different numerical methods such as finite element, finite volume, boundary elements methods, meshless methods, fundamental solutions, and spectral and wavelet methods. In [1], Parzlivand and Shahrezaee presented a numerical technique based on a combination of collocation method and radial basis functions to solve parabolic partial differential equations with time-dependent coefficients. Dehghan and Shokri [2] have employed the thin-plate splines RBFs in the collocation RBF method to solve the two-space dimensional linear hyperbolic equation with variable coefficients, subject to appropriate initial and Dirichlet boundary conditions. A new version of the method of approximate particular solutions using radial basis functions has been proposed in [3] by Jiang et al. for solving elliptic partial differential equations with variable coefficients. In this paper [3], a comparison to Kansa’s method and the method of fundamental solutions has also been investigated. Some other investigation of meshless methods to PDEs with variable coefficients is given in [4–6]. We can also mention the works of Gu et al. [7, 8], where they used local versions of meshless methods leading to sparse matrices.In most published works, these methods are based on first discretizing the time variable by applying any time time-stepping algorithms as implicit, explicit, Runge-Kutta, or others; and seeking the approximate solution at each instantt in a space domain problem.The space-time methods have seen some investigation these last decades. Among them, we can mention the space-time finite element method developed by Tayfun et al. in [9] for the computation of fluid-structure interaction problems. Klaij et al. have also developed the space-time discontinuous Galerkin finite element method for solving compressible Navier-Stokes equations in [10] and advection-diffusion problems in [11]. The technique has also been applied to shallow water flows by Ambati et al. in [12].Although few works were published concerning the space-time meshless method for solving PDEs with independent variable coefficients [13–18], up to date and to our best knowledge, there is no investigation on the application of space-time meshless method for solving PDEs with either variable or time-dependent variable coefficients. In this paper, we investigate the application of space-time localized meshless collocation radial basis functions developed in [13] to a general second-order parabolic and hyperbolic problems, with variable coefficients.The method is based on firstly transforming the consideredd-dimensional evolutionary problem in space as a d+1-dimensional space-time one. The formulation starts by combining the d-dimensional vector space variable and 1-dimensional time variable in one d+1-dimensional variable vector and defining the space-time domain and its boundary. Then, the boundary conditions are sited on the global space-time boundary domain for PDEs with the first derivative with respect to time and on just one part of the boundary for PDEs with the second derivative with respect to time. The method does not need first to discretize all derivatives with respect to time and solve the problem in the space domain at each time level, as it is usually the case with many numerical methods. The developed formulation leads always to a square algebraic resulting system and benefit from the following advantages:
(i)
The time stability analysis is not discussed as it is the case for other time-stepping schemes as implicit, explicit,θ-method, etc.(ii)
Reducing the computational time as there is no need to recompute the matrix for the resulting algebraic system at each time level, unlike the case for others time integration methods used to solve PDEs with time-dependent coefficients(iii)
Solving hyperbolic partial differential equations with variable coefficients as inverse problem, since just the hyperbolic equation with constant coefficients is less discussed in the literature and needs more sophisticated time integration techniqueThe paper is organized as follows. In Section2, we introduce the formulation of the parabolic and hyperbolic problem as space-time problem and the space-time localized RBF method implementation. Section 3 is devoted to the discussion of results obtained by solving different parabolic and hyperbolic examples in two- and three-dimensions in regular and irregular domains. A comparison of the given technique to the fourth-order Runge-Kutta method is also given in Section 4. We conclude in Section 5.


## 2. Space-Time Localized RBF Method Formulation

To recall the space-time localized RBF method defined in [13], let Ω⊂ℝd be a bounded domain with a sufficiently regular boundary ∂Ω and consider the following time-dependent boundary value problem
(1)∂∂tux,t+Px,tux,t=fx,t∀x∈Ω,∀t∈0,T,(2)Bux,t=gx,t∀x∈∂Ω,∀t∈0,T,(3)ux,0=u0x∀x∈Ω,or
(4)∂2∂t2ux,t+Px,tux,t=fx,t∀x∈Ω,∀t∈0,T,(5)Bux,t=gx,t∀x∈∂Ω,∀t∈0,T,(6)ux,0=u0x∀x∈Ω,(7)∂∂tux,0=u1x∀x∈Ω,where Px,t is a differential operator of second order with variable coefficients of the form:
(8)Px,tu=∑i=1i=d∑j=1j=dax,t∂2u∂xi∂xj+bx,tu,and B is a linear boundary operator depending on the kind of problem treated. The given functions f, g, u0, and u1 are assumed to be sufficiently regular.The formulation of the evolutionary problems given by (1-2-3) and (4-5-6-7) as a space-time one starts by combining the space variablex and the time variable t in one d+1-dimensional vector. The constructed variable vector belongs to the space-time domain ΩT=Ω×0,T represented by Figure 1 for some simple cases. The boundary of the new formulated domain ΩT is given by ∂ΩT=Ω×t=T,∂Ω×0,T,Ω×t=0.Figure 1
Space-time domain concept.As discussed in [13], the formulation of the d-dimensional problem as a d+1-dimensional space-time one depends on the type of equations considered. For the problem given by Eqs. (1)–(3), Eq. (1) given by ∂/∂tu+Px,tu=fx,t,∀x,t∈Ω×0,T is considered as a domain equation in ΩT and Eqs. (2) and (3) as boundary conditions on ∂Ω×0,T and Ω×t=0, respectively. To complete the set of boundary conditions on the ∂ΩT, we define Eq. (1) as a boundary condition on Ω×t=T. The formulation is completed, and the new form of the system is
(9)∂∂tux,t+Px,tux,t=fx,t∀x,t∈ΩT,(10)∂∂tux,t+Px,tux,t=fx,t∀x,t∈Ω×t=T,(11)Bux,t=gx,t∀x∈∂Ω,×0,T,(12)ux,t=u0x∀x∈Ω×t=0.For the case of the hyperbolic equation, we do not need any boundary condition on the part on the part of the space-time boundary characterized byΩ×t=T. The problem is solved as an ill-posed one, and the algebraic system is square (see [1]).Then, the system has the new form:(13)∂2∂t2ux,t+Px,tux,t=fx,t∀x,t∈ΩT,(14)∂∂tux,t=u1x∀x∈Ω×t=0,(15)Bux,t=gx,t∀x∈∂Ω,×0,T,(16)ux,t=u0x∀x∈Ω×t=0.Figure2 summarizes the defined space-time problems for the two considered cases.Figure 2
Domain of space-time formulation.In general, to construct the approximate solution of this news system, we select a number of distinct centers pointsx^j=xj,tjj=1Ni and x^j=xj,tjj=Ni+1N in the domain ΩT and on its boundary ∂ΩT, respectively. The finite-dimensional space VΦ that contains the RBF interpolation functions is defined by
(17)VΦ=spanΦ∥.−x^1∥,Φ∥.−x^2∥,⋯,Φ∥.−x^N∥+ℙmd+1,where m is the order of the used conditionally positive defined function Φ. The approximate function is given by
(18)sVx^=sVx,t=∑j=1NαjΦjx^.We can mention that the space-time domainΩT=Ω×0,T satisfies the interior cone condition since the domain Ω and the interval 0,T are also satisfying the interior cone condition. So, the lemma of error interpolation cited in [19] is satisfied for the case of the space-time domain (see also [13]).The general form of the system of Eqs. (9)–(12) or (13)–(16) can be written under the following form:
(19)Lx^ux^=fx^∀x^∈ΩT,Bx^ux^=hx^∀x^∈∂ΩT,where Lx^=∂/∂t+Px,t on ΩT, Bx^=Lx^B1 and h=fgu0t on ∂ΩT for the parabolic case and Lx^=∂2/∂t2+Px,t on ΩT, Bx^=∂/∂tB1 and h=u1gu0t on ∂ΩT for the hyperbolic case.Following Li et al. [20] and by Yao et al. [21], a localized influence domain ΩTs for each point x^s∈ΩT is selected. It contains a number ns of neighboring points x^ksk=1ns∈ΩTs to the selected center x^s. (see Figure 3). Then, the proximate function in the subdomain ΩTs is given by the expansion
(20)ux^s≃u^x^s=∑k=1nsαkΦ∥x^s−x^ks∥,where αkk=1ns are the unknown coefficients, · is the Euclidean norm, and Φ is the chosen RBF.Figure 3
Schematic showing five-node, nine-node, and thirteen-node local domainsΩs for 2D (left) and nine-node local domain for 3D (right).Applying Eq. (20) to x^ksk=1ns⊂ΩTs set points, we get ns×ns linear equations given by
(21)u^s=Φsαs,where Φs=Φ∥x^is−x^js∥1≤i,j≤ns.Then, the problem of seeking the expansion coefficientsαkk=1ns is transformed into a determination of the values of the solution u^ at each center points x^ksk=1ns⊂ΩTs by using the equation
(22)αs=Φs−1u^s.Forx^s∈ΩTs, we apply the differential operator Lx^ to Eq. (20) to obtain the following equation
(23)Lx^u^x^s=∑k=1nsαkLx^Φ∥x^s−x^ks∥=Ysu^s,where u^s=ux^1s,ux^2s,⋯,ux^nss and Ys=Lx^ΦsΦs−1. The vector u^=ux^1,ux^2,⋯,ux^N is incorporated in the system of (23) by adding zeros at the proper locations based on the mapping of u^s to u^, and considering the Y1×N as the global expansion of Y1×nss. The global system of Eq. (23) is then written under the form
(24)fx^s=Lx^u^x^s=Yx^su^.In the same way, selecting a centerx^s on the boundary ∂ΩT, we have
(25)Bx^ux^s=∑k=1nsαkBx^Φ∥x^s−x^ks∥=Θsu^,where Θs=Bx^ΦsΦs−1. In a global form, the system is then written as
(26)hx^s=Bu^x^s=Θx^su^,where Θ is the expansion of Θs by adding zeros.By collocating at all centers pointsx^jj=1N using Eqs. (24) and (26), we get the following sparse linear system of equations
(27)YΘu^=fh,where
(28)Y=Yx^1....Yx^NiandΘ=Θx^Ni+1....Θx^N.The approximate solution at the interpolation pointsV^x^jj=1N can be obtained by solving the above sparse linear system of equations.


## 3. Numerical Simulations

The investigated test examples in this section are two- and three-dimensional parabolic and hyperbolic equations with variable coefficients given as follows:(29)∂u∂t+ax,tΔu+bx,tu=finΩ×0,T,(30)∂2u∂t2+ax,tΔu+bx,tu=finΩ×0,T.Knowing the shape of the space-time domainΩT, the distributed nodes can be done uniformly or randomly without making any distinction between space and time variables. Through these numerical simulations, ns denotes the number of neighboring points in an influence space-time domain ΩT. To balance between the approximation quality of PDEs solution and sparsity of the algebraic matrix, a reasonable value of ns is chosen to be at least 2dimΩT+1 when dealing with Dirichlet boundary condition. More points can be added in the case of the Neumann condition depending the on problem treated (see [13] for more details). Seeking optimal values of shape parameter is still an open subject in the area of numerical approximation of PDEs. Many techniques have been proposed in the literature without any theoretical analysis. The N=Nd+Nb×Nt is the total number of nodes used, where Nb and Nd are the number of boundary and interior nodes of space domain Ω, respectively, and Nt is the number of nodes on the time axe which can be either uniform or random.In all proposed simulations, Neumann and Dirichlet boundary conditions are considered for parabolic and hyperbolic equations to show the sensibility of the method in respect to different boundary conditions on the space domain. Although any RBF can be used to test the technique, the multiquadric radial basis function (MQ) is the function used for all given simulations. The MQ-RBF is defined byΦr=1+ε2r2, where ε is the shape parameter and r is the distance between two nodes.The validity and the accuracy of the presented technique are demonstrated by considering the following discrete error-norms; the maximum absolute error (MAE), the root mean squared error (RMSE), and theLer1 relative errors are defined by
(31)MAE=max1≤j≤Nu^j−uj,RMSE=1N∑j=1Nu^j−uj2,Ler1=∑j=1Nu^j−uj∑j=1Nuj,where uj and u^j are the exact and the approximate solutions at xj,tj, respectively.


### 3.1. One-Dimensional Test and Comparison with Runge-Kutta Method

To show the efficiency of the space-time method in word of stability and CPU time, we compare it to the 4th order Runge-Kutta method. The investigated example is the one-dimensional convection-diffusion PDE with variable coefficients defined as follows:(32)∂u∂t−x+t∂2u∂x2+2x+3t∂u∂x=0.The initial and the boundary conditions are defined according to the analytical solution given byux,t=exp2x−t2, and the space domain is Ω=0,1.As the coefficients are depending ont, the matrix construction is called 4 times at each time step using the 4th order Runge-Kutta method, which enlarges the CPU time. According to Table 1, the CPU computation time of the Runge-Kutta method is much larger than that of the space-time technique. We can also mention that time stability is not guaranteed in the Runge-Kutta method or any other time-stepping algorithms. Besides the implementation simplicity of the demonstrated technique, results show that our method is faster and more accurate than the used Runge-Kutta method.Table 1
Comparison of numerical results obtained using space-time method (top) and RK4 method (bottom).
T0.10.20.511Space-timeNt2020202040CPU time0.072310.067170.066570.068440.18470MAE8.80E-041.59E-033.22E-031.02E-038.63E-04RMSE5.96E-041.12E-032.29E-037.31E-046.16E-04RK4 methoddt0.0010.0010.0010.0010.0005CPU time1.687752.924397.5114813.7163827.19279MAE3.95E-034.63E-035.72E-03—2.72E-03RMSE1.42E-031.83E-032.50E-03—6.92E-04Table1 shows that the space-time method is stable and accurate compared to the RK4 method. The space-time method converges for different values of the final time T, keeping the same value of Nt. From the same Table 1, it can be remarked that the RK4 method diverges for the simulation with T=1 and dt=0.001, which means that a smaller value of dt is needed to get the convergence. Regarding to the Table 1 and incase the results of the space-time method is discussed according to dt, we can remark that dt is ranged from 0.05 to 0.025 and results still more accurate. All results are obtained using the same number of nodes on the space domain Nx=20 and the final time T ranges from 0.1 to 1. The shape parameter of the MQ-RBF is chosen to be ε=0.1.


### 3.2. Two-Dimensional Parabolic Equations

For examples treated herein, we consider the considered as computation space domainΩ, a two-planar domain given by Ω=x,y∈ℝ2:x,y≤1. Two different distributions of nodes, uniform and circular, illustrated in Figure 4, are considered. For uniform distribution of nodes, we take Nd=326 and Nb=50 and Nd=320 and Nb=50 for circular one.Example 1.
As a first simulation, we consider the parabolic example with the analytical solution given by(33)ux,y,t=sinx+yexp−t,x,t∈Ω×0,1where ax,y,t=1/1+x2+y2 and bx,y,t=x+2y+2.Figure 4
Uniform and circular distributions of nodes for the 2D space domain and the 3D space-time domain.The problem was solved using different values of bothNt and shape parameter ε. Considering Dirichlet boundary condition on ∂Ω, Tables 2 and 3 illustrate results for Nt=10 and Nt=20, respectively. They also show the influence of shape parameter ε on the approximate solution. Using ε=0.25 and ε=0.43 for uniform and circular distributions of nodes, respectively, Figure 5 shows the absolute error of the approximate solution at t=1 on the computation space-domain Ω for Dirichlet boundary condition on ∂Ω×0,T. Using uniform distribution of nodes, the obtained maximum absolute error and the root mean square error are MAE=2.71×10−5 and RMSE=1.09×10−5, respectively. For the circular distributed nodes, we obtained MAE=5.47×10−5 and RMSE=2.37×10−5. As the ΩT⊂ℝ3 and the boundary condition are Dirichlet type, the used value of the number of neighboring points is ns=7.Table 2
Results for different values ofε in Example 1 with Nt=10. Uniformly distributed nodes (left) and circular distributed nodes (right), with Dirichlet conditions.
εMAERMSELer1MAERMSELer10.11.61E-045.60E-051.14E-043.70E-041.21E-042.50E-040.21.44E-045.03E-051.02E-043.38E-041.10E-042.28E-040.31.15E-044.07E-058.46E-052.84E-049.25E-051.91E-040.47.25E-052.78E-055.91E-052.07E-046.72E-051.39E-040.59.80E-051.91E-052.80E-051.11E-043.64E-057.69E-050.61.57E-043.65E-056.99E-051.45E-043.37E-056.10E-050.72.33E-047.19E-051.49E-042.42E-048.56E-051.76E-04Table 3
Results for different values ofε in Example 1 with Nt=20. Uniform spaced nodes (left) and circular nodes (right), using Dirichlet conditions.
εMAERMSELer1MAERMSELer10.17.80E-052.78E-056.04E-053.08E-049.59E-052.03E-040.26.70E-051.93E-054.17E-052.68E-048.22E-051.73E-040.36.37E-057.79E-061.33E-052.04E-045.94E-051.23E-040.45.91E-052.18E-054.94E-051.25E-042.91E-055.80E-050.51.44E-045.19E-051.15E-041.03E-043.10E-056.78E-050.62.58E-049.19E-052.02E-042.28E-048.58E-051.89E-040.74.04E-041.43E-043.12E-044.29E-041.59E-043.51E-04Figure 5
Absolute errors att=1 with Nt=10. Uniform (a) and circular (b). Example 1 using Dirichlet conditions.
(a)(b)To demonstrate the influence of the boundary condition on the approximate solution and then on the technique, mixed Neumann-Dirichlet boundary condition is considered on∂Ω. Dirichlet condition was used on the part of the space domain boundary defined by Γ1=x,y∈∂Ω,y>0 and Neumann condition on the other part Γ2=∂Ω−Γ1. For this case, the number of neighboring points used is ns=9, if the collocation point considered belongs to the boundary of Neumann condition Γ2 and ns=7 for the all other centers. Tables 4 and 5 present the results obtained for Nt=10 and 20 for different distribution of nodes.Table 4
Results for different values ofε in Example 1 with Nt=10. Uniform spaced nodes (left) and circular nodes (right), using Dirichlet and Neumann conditions.
εMAERMSELer1MAERMSELer10.14.41E-041.13E-042.28E-043.89E-041.58E-043.63E-040.24.18E-041.04E-042.09E-043.79E-041.46E-043.33E-040.33.79E-048.98E-051.77E-043.67E-041.26E-042.83E-040.43.22E-046.90E-051.32E-043.84E-049.96E-052.15E-040.52.45E-044.52E-057.97E-054.05E-047.71E-051.48E-040.61.61E-044.20E-058.33E-054.26E-048.79E-051.68E-040.72.73E-048.41E-051.82E-044.61E-041.47E-043.09E-04Table 5
Results for different values ofε in Example 1 with Nt=20. Uniform spaced nodes (left) and circular nodes (right), using Dirichlet and Neumann conditions.
εMAERMSELer1MAERMSELer10.13.25E-046.79E-051.34E-043.90E-041.19E-042.76E-040.22.90E-045.54E-051.06E-043.79E-041.04E-042.36E-040.32.32E-043.69E-056.42E-053.61E-048.18E-051.74E-040.41.83E-043.11E-056.57E-053.36E-046.73E-051.35E-040.51.67E-046.63E-051.57E-043.61E-049.30E-051.86E-040.63.10E-041.25E-042.93E-045.38E-041.61E-043.46E-040.75.26E-042.02E-044.67E-047.90E-042.58E-045.69E-04Figure6 shows the maximum and root mean errors of the approximate solution at t=1 for mixed boundary condition and different distributed nodes. Their values are MAE=8.55×10−5 and RMSE=4.20×10−5 for uniformly distributed collocation points using ε=0.47. For the circular distributed nodes, we obtain MAE=1.78×10−4 and RMSE=4.91×10−5 for ε=0.43. We can conclude that accurate results have been obtained for this first example.Example 2.
The second simulated case is the problem with the analytical solution represented by(34)ux,y,t=12−texp−x2+y22−t,Ω×0,1.where ax,y,t=1+x22−t and bx,y,t=xyt+1.Figure 6
Absolute Errors att=1 with Nt=10. Uniform (a) and circular (b). Example 1 using Dirichlet and Neumann conditions.
(a)(b)The computational domain and boundaries are the same as for Example1. The simulation of this second example shows some difficulties compared to Example 1 because of the pseudosingularity that present at points with t=2. Keeping the same neighboring points value ns=7, Tables 6 and 7 give obtained results for different value of Nt considering Dirichlet boundary condition.Table 6
Results for different values ofε in Example 2 with Nt=10. Uniform spaced nodes (left) and circular nodes (right), with Dirichlet conditions.
εMAERMSELer1MAERMSELer10.11.71E-032.99E-043.80E-047.25E-041.04E-041.14E-040.21.55E-032.52E-043.06E-045.67E-041.21E-041.41E-040.31.29E-031.78E-041.83E-048.37E-041.88E-042.55E-040.49.04E-041.12E-041.01E-041.23E-033.11E-044.54E-040.59.00E-041.94E-042.82E-041.76E-034.88E-047.31E-040.61.58E-033.95E-045.92E-042.46E-037.23E-041.09E-030.72.43E-036.64E-041.00E-033.34E-031.02E-031.55E-03Table 7
Results for different values ofε in Example 2 with Nt=20. Uniform spaced nodes (left) and circular nodes (right), with Dirichlet conditions.
εMAERMSELer1MAERMSELer10.11.55E-033.46E-044.76E-045.66E-049.00E-051.15E-040.21.39E-032.98E-044.04E-043.97E-046.58E-058.62E-050.31.11E-032.17E-042.81E-044.02E-041.01E-041.47E-040.47.08E-041.04E-041.12E-046.63E-042.17E-043.43E-040.51.76E-048.89E-051.55E-041.11E-033.88E-046.14E-040.67.13E-042.87E-044.71E-041.74E-036.14E-049.67E-040.71.52E-035.49E-048.73E-042.56E-039.01E-041.42E-03From Figure7, we can remark that at t=1, we obtain MAE=2.51×10−4 with ε=0.53 for the uniform distributed nodes and MAE=4.13×10−4 with ε=0.26 for the circular distributed nodes.Figure 7
Absolute errors att=1 with Nt=10. Uniform (a) and circular (b). Example 2 with Dirichlet conditions.
(a)(b)Following the same step of simulation done before, Example2 was then solved by considering mixed Neumann-Dirichlet boundary condition: Neumann condition on Γ2 and Dirichlet conditions on Γ1. On the boundary with the Neumann condition, the number of selected neighboring point is set to be ns=9. Tables 8 and 9 give the found results found.Example 3.
The investigated problem is the widely used convection-diffusion flow problem of a Gaussian pulse given by the equation:(35)∂u∂tx,y,t+vx,y,t∇ux,y,t−νΔux,y,t=fx,y,t∀x,y,t∈Ω×0,T,where vx,y,t=−4y,4xT is the velocity vector, and ν is a constant. Initial and boundary conditions are taken according to the analytical solution:
(36)ux,y,t=σ2σ2+4νtexp−x−0.52+y−0.52σ2+4νt.Table 8
Results for different values ofε in Example 2 with Nt=10. Uniform spaced nodes (left) and circular nodes (right), with Neumann conditions.
εMAERMSELer1MAERMSELer10.11.72E-033.29E-044.53E-041.19E-032.07E-042.92E-040.21.51E-032.71E-043.59E-041.43E-031.91E-042.33E-040.31.42E-032.39E-042.95E-041.85E-033.09E-043.91E-040.41.85E-033.67E-044.88E-042.46E-035.87E-048.68E-040.52.44E-036.59E-049.97E-043.29E-039.97E-041.55E-030.63.23E-031.08E-031.71E-034.37E-031.54E-032.45E-030.74.23E-031.63E-032.63E-035.76E-032.24E-033.58E-03Table 9
Results for different values ofε in Example 2 with Nt=20. Uniform spaced nodes (left) and circular nodes (right), with Neumann conditions.
εMAERMSELer1MAERMSELer10.11.61E-033.93E-045.79E-048.84E-042.48E-043.81E-040.21.39E-033.21E-044.64E-041.11E-031.69E-042.43E-040.31.08E-032.30E-043.28E-041.51E-031.89E-042.24E-040.41.47E-032.61E-043.25E-042.09E-034.39E-046.63E-040.52.03E-035.11E-047.80E-042.88E-038.30E-041.33E-030.62.76E-039.03E-041.47E-033.90E-031.35E-032.18E-030.73.70E-031.41E-032.32E-035.20E-031.99E-033.21E-03In our simulation, we takeσ2=0.001 and ν=0.01. The computational space-time domain ΩT is 0,1×0,1×0,π/2. In Figure 8, we show the numerical solution at final time t=T face to the analytical solution, the total number of nodes is N=303. Table 10 displays the absolute maximum error obtained for different values of N; it can be remarked that the rate of convergence is near quadric. In all simulations for Example 3, the shape parameter is chosen to be ε=5.Figure 8
Exact (a) and numerical (b) solutions att=π/2 for the Example 3 with N=303 and ε=5.
(a)(b)Table 10
Maximum absolute error and rate of convergence att=T for different values of N in Example 3 using ε=5 and ΩT=0,1×0,1×0,π/2.
NMAEhROC1536.70E-040.0714—2033.74E-040.05261.912532.35E-040.04171.993031.62E-040.03451.96


### 3.3. Two-Dimensional Hyperbolic Equations

Although it has been shown that the technique gives good results solving the parabolic equation with different boundary conditions, the hyperbolic equation is still an important and interesting equation to be tested using the described methodology. In all simulations, the problem is treated as an ill-posed one as it has mentioned before and the algebraic matrix obtained is square.Example 4.
Considering a two-dimensional wave problem defined by the following equation on the domainΩ×0,T with Ω=0,12(37)∂2u∂t2x,y,t−aΔux,y,t=fx,y,t∀x,y,t∈Ω×0,T,with two initial conditions ux,0=u0x and ∂u/∂tx,y,0=hx,y and Dirichlet boundary condition on ∂Ω.Takingax,y,t=xye−2t, the functions u0, h, and f are chosen according to the analytical solution given by ux,t=sinxsinysint.Applying the “ill-posed-problem” technique [22], we obtain accurate results by setting the number of neighboring points to be ns=9 in the entire domain. Table 11 shows results obtained of different errors, taking ε=0.1 and using the space-time domain 0,13. Figure 9 shows that the MAE error decreases when increasing the number of nodes.Table 11
Errors for different values ofNx, Ny, and Nt of Example 4 using ε=0.1 and ΩT=0,13.
NxNyNtMAERMSELer15551.80E-035.13E-042.69E-031010104.59E-041.26E-048.08E-041515151.91E-045.23E-053.53E-042020201.06E-042.86E-051.97E-042525256.68E-051.81E-051.25E-04Figure 9
Maximum absolute errors vs.Nx with Nx=Ny=Nt of Example 4 using ε=0.1.As in the example of the one-dimensional problem, the convergence rate of the technique using different uniformly distributed sources points andε=0.1 is given in Table 12 showing that the order of convergence is quadratic.Table 12
Rate of convergence for wave Example4 using ε=0.1.
hMAEROCRMSEROCLer1ROC0.250001.70E-03—4.61E-04—2.44E-03—0.111114.46E-041.64601.22E-041.63637.86E-041.39700.071431.78E-042.07964.92E-052.05963.34E-041.93650.052631.02E-041.81902.80E-051.85161.94E-041.78330.041676.43E-051.98351.78E-051.93831.25E-041.8771Att=1, the maximum absolute error is MAE=6.68×10−5, and the root mean square error is RMSE=1.81×10−5; these values are obtained with ε=0.1 and N=253. We remark that the maximum error is obtained at t=T (see Figure 10).Figure 10
Absolute errors withNt=25 and ε=0.1 at t=1, for Example 4.


### 3.4. Three-Dimensional Parabolic and Hyperbolic Examples

In this section, we show the robustness and the accuracy of the proposed method implemented in high dimension problems. The traditional three-dimensional problems are transformed into 4D problems.Example 5.
The problem investigated is the 3D heat transfer problem given by the equation:(38)∂u∂tx,y,z,t−ax,y,z,tΔux,y,z,t=fx,y,z,t.Withax,y,z,t=xyzet, the initial and boundary conditions are chosen in such a way that the analytical solution is given by the function ux,y,z,t=e−tsinxsinysinz. The space-time domain is chosen Ω=0,14 and the right-hand side is fx,y,z,t=sinxsinysinze−t−1+3xyzet.Takingns=9, Table 13 shows the errors found in two simulations. The total number of nodes is N=Nx×Nt×Ny×Nz=54 and N=94.Table 13
Errors for Example5 using two different values of N.
NεMAERMSELer1540.94.88E-49.28E-55.45E-4940.82.23E-43.62E-52.75E-4We can remark that good accuracy is obtained even with weak number of nodes. We note that similar accuracy is found using different values of the shape parameterε in the interval 0,1.5Example 6.
The considered hyperbolic 3D problem is given by(39)∂2u∂t2x,y,z,t−ax,y,z,y,tΔux,y,z,t=fx,y,z,t.We useax,y,z,t=xyze−t and ε=0.12. The computational space-time domain is Ω=0,14. The boundary and initial conditions and the function f are chosen according to the solution ux,y,z,t=sinxsinysinzsint. In this example, the number of neighboring points is set to ns=11 on the entire domain. The total number of nodes is N=Nx×Nt×Ny×Nz=94. The errors obtained for this example are MAE=1.54×10−3, RMSE=9.23×10−5, and Ler1=6.00×10−4. The same ill-posed technique is used with the squared algebraic linear system.


## 3.1. One-Dimensional Test and Comparison with Runge-Kutta Method

To show the efficiency of the space-time method in word of stability and CPU time, we compare it to the 4th order Runge-Kutta method. The investigated example is the one-dimensional convection-diffusion PDE with variable coefficients defined as follows:(32)∂u∂t−x+t∂2u∂x2+2x+3t∂u∂x=0.The initial and the boundary conditions are defined according to the analytical solution given byux,t=exp2x−t2, and the space domain is Ω=0,1.As the coefficients are depending ont, the matrix construction is called 4 times at each time step using the 4th order Runge-Kutta method, which enlarges the CPU time. According to Table 1, the CPU computation time of the Runge-Kutta method is much larger than that of the space-time technique. We can also mention that time stability is not guaranteed in the Runge-Kutta method or any other time-stepping algorithms. Besides the implementation simplicity of the demonstrated technique, results show that our method is faster and more accurate than the used Runge-Kutta method.Table 1
Comparison of numerical results obtained using space-time method (top) and RK4 method (bottom).
T0.10.20.511Space-timeNt2020202040CPU time0.072310.067170.066570.068440.18470MAE8.80E-041.59E-033.22E-031.02E-038.63E-04RMSE5.96E-041.12E-032.29E-037.31E-046.16E-04RK4 methoddt0.0010.0010.0010.0010.0005CPU time1.687752.924397.5114813.7163827.19279MAE3.95E-034.63E-035.72E-03—2.72E-03RMSE1.42E-031.83E-032.50E-03—6.92E-04Table1 shows that the space-time method is stable and accurate compared to the RK4 method. The space-time method converges for different values of the final time T, keeping the same value of Nt. From the same Table 1, it can be remarked that the RK4 method diverges for the simulation with T=1 and dt=0.001, which means that a smaller value of dt is needed to get the convergence. Regarding to the Table 1 and incase the results of the space-time method is discussed according to dt, we can remark that dt is ranged from 0.05 to 0.025 and results still more accurate. All results are obtained using the same number of nodes on the space domain Nx=20 and the final time T ranges from 0.1 to 1. The shape parameter of the MQ-RBF is chosen to be ε=0.1.


## 3.2. Two-Dimensional Parabolic Equations

For examples treated herein, we consider the considered as computation space domainΩ, a two-planar domain given by Ω=x,y∈ℝ2:x,y≤1. Two different distributions of nodes, uniform and circular, illustrated in Figure 4, are considered. For uniform distribution of nodes, we take Nd=326 and Nb=50 and Nd=320 and Nb=50 for circular one.Example 1.
As a first simulation, we consider the parabolic example with the analytical solution given by(33)ux,y,t=sinx+yexp−t,x,t∈Ω×0,1where ax,y,t=1/1+x2+y2 and bx,y,t=x+2y+2.Figure 4
Uniform and circular distributions of nodes for the 2D space domain and the 3D space-time domain.The problem was solved using different values of bothNt and shape parameter ε. Considering Dirichlet boundary condition on ∂Ω, Tables 2 and 3 illustrate results for Nt=10 and Nt=20, respectively. They also show the influence of shape parameter ε on the approximate solution. Using ε=0.25 and ε=0.43 for uniform and circular distributions of nodes, respectively, Figure 5 shows the absolute error of the approximate solution at t=1 on the computation space-domain Ω for Dirichlet boundary condition on ∂Ω×0,T. Using uniform distribution of nodes, the obtained maximum absolute error and the root mean square error are MAE=2.71×10−5 and RMSE=1.09×10−5, respectively. For the circular distributed nodes, we obtained MAE=5.47×10−5 and RMSE=2.37×10−5. As the ΩT⊂ℝ3 and the boundary condition are Dirichlet type, the used value of the number of neighboring points is ns=7.Table 2
Results for different values ofε in Example 1 with Nt=10. Uniformly distributed nodes (left) and circular distributed nodes (right), with Dirichlet conditions.
εMAERMSELer1MAERMSELer10.11.61E-045.60E-051.14E-043.70E-041.21E-042.50E-040.21.44E-045.03E-051.02E-043.38E-041.10E-042.28E-040.31.15E-044.07E-058.46E-052.84E-049.25E-051.91E-040.47.25E-052.78E-055.91E-052.07E-046.72E-051.39E-040.59.80E-051.91E-052.80E-051.11E-043.64E-057.69E-050.61.57E-043.65E-056.99E-051.45E-043.37E-056.10E-050.72.33E-047.19E-051.49E-042.42E-048.56E-051.76E-04Table 3
Results for different values ofε in Example 1 with Nt=20. Uniform spaced nodes (left) and circular nodes (right), using Dirichlet conditions.
εMAERMSELer1MAERMSELer10.17.80E-052.78E-056.04E-053.08E-049.59E-052.03E-040.26.70E-051.93E-054.17E-052.68E-048.22E-051.73E-040.36.37E-057.79E-061.33E-052.04E-045.94E-051.23E-040.45.91E-052.18E-054.94E-051.25E-042.91E-055.80E-050.51.44E-045.19E-051.15E-041.03E-043.10E-056.78E-050.62.58E-049.19E-052.02E-042.28E-048.58E-051.89E-040.74.04E-041.43E-043.12E-044.29E-041.59E-043.51E-04Figure 5
Absolute errors att=1 with Nt=10. Uniform (a) and circular (b). Example 1 using Dirichlet conditions.
(a)(b)To demonstrate the influence of the boundary condition on the approximate solution and then on the technique, mixed Neumann-Dirichlet boundary condition is considered on∂Ω. Dirichlet condition was used on the part of the space domain boundary defined by Γ1=x,y∈∂Ω,y>0 and Neumann condition on the other part Γ2=∂Ω−Γ1. For this case, the number of neighboring points used is ns=9, if the collocation point considered belongs to the boundary of Neumann condition Γ2 and ns=7 for the all other centers. Tables 4 and 5 present the results obtained for Nt=10 and 20 for different distribution of nodes.Table 4
Results for different values ofε in Example 1 with Nt=10. Uniform spaced nodes (left) and circular nodes (right), using Dirichlet and Neumann conditions.
εMAERMSELer1MAERMSELer10.14.41E-041.13E-042.28E-043.89E-041.58E-043.63E-040.24.18E-041.04E-042.09E-043.79E-041.46E-043.33E-040.33.79E-048.98E-051.77E-043.67E-041.26E-042.83E-040.43.22E-046.90E-051.32E-043.84E-049.96E-052.15E-040.52.45E-044.52E-057.97E-054.05E-047.71E-051.48E-040.61.61E-044.20E-058.33E-054.26E-048.79E-051.68E-040.72.73E-048.41E-051.82E-044.61E-041.47E-043.09E-04Table 5
Results for different values ofε in Example 1 with Nt=20. Uniform spaced nodes (left) and circular nodes (right), using Dirichlet and Neumann conditions.
εMAERMSELer1MAERMSELer10.13.25E-046.79E-051.34E-043.90E-041.19E-042.76E-040.22.90E-045.54E-051.06E-043.79E-041.04E-042.36E-040.32.32E-043.69E-056.42E-053.61E-048.18E-051.74E-040.41.83E-043.11E-056.57E-053.36E-046.73E-051.35E-040.51.67E-046.63E-051.57E-043.61E-049.30E-051.86E-040.63.10E-041.25E-042.93E-045.38E-041.61E-043.46E-040.75.26E-042.02E-044.67E-047.90E-042.58E-045.69E-04Figure6 shows the maximum and root mean errors of the approximate solution at t=1 for mixed boundary condition and different distributed nodes. Their values are MAE=8.55×10−5 and RMSE=4.20×10−5 for uniformly distributed collocation points using ε=0.47. For the circular distributed nodes, we obtain MAE=1.78×10−4 and RMSE=4.91×10−5 for ε=0.43. We can conclude that accurate results have been obtained for this first example.Example 2.
The second simulated case is the problem with the analytical solution represented by(34)ux,y,t=12−texp−x2+y22−t,Ω×0,1.where ax,y,t=1+x22−t and bx,y,t=xyt+1.Figure 6
Absolute Errors att=1 with Nt=10. Uniform (a) and circular (b). Example 1 using Dirichlet and Neumann conditions.
(a)(b)The computational domain and boundaries are the same as for Example1. The simulation of this second example shows some difficulties compared to Example 1 because of the pseudosingularity that present at points with t=2. Keeping the same neighboring points value ns=7, Tables 6 and 7 give obtained results for different value of Nt considering Dirichlet boundary condition.Table 6
Results for different values ofε in Example 2 with Nt=10. Uniform spaced nodes (left) and circular nodes (right), with Dirichlet conditions.
εMAERMSELer1MAERMSELer10.11.71E-032.99E-043.80E-047.25E-041.04E-041.14E-040.21.55E-032.52E-043.06E-045.67E-041.21E-041.41E-040.31.29E-031.78E-041.83E-048.37E-041.88E-042.55E-040.49.04E-041.12E-041.01E-041.23E-033.11E-044.54E-040.59.00E-041.94E-042.82E-041.76E-034.88E-047.31E-040.61.58E-033.95E-045.92E-042.46E-037.23E-041.09E-030.72.43E-036.64E-041.00E-033.34E-031.02E-031.55E-03Table 7
Results for different values ofε in Example 2 with Nt=20. Uniform spaced nodes (left) and circular nodes (right), with Dirichlet conditions.
εMAERMSELer1MAERMSELer10.11.55E-033.46E-044.76E-045.66E-049.00E-051.15E-040.21.39E-032.98E-044.04E-043.97E-046.58E-058.62E-050.31.11E-032.17E-042.81E-044.02E-041.01E-041.47E-040.47.08E-041.04E-041.12E-046.63E-042.17E-043.43E-040.51.76E-048.89E-051.55E-041.11E-033.88E-046.14E-040.67.13E-042.87E-044.71E-041.74E-036.14E-049.67E-040.71.52E-035.49E-048.73E-042.56E-039.01E-041.42E-03From Figure7, we can remark that at t=1, we obtain MAE=2.51×10−4 with ε=0.53 for the uniform distributed nodes and MAE=4.13×10−4 with ε=0.26 for the circular distributed nodes.Figure 7
Absolute errors att=1 with Nt=10. Uniform (a) and circular (b). Example 2 with Dirichlet conditions.
(a)(b)Following the same step of simulation done before, Example2 was then solved by considering mixed Neumann-Dirichlet boundary condition: Neumann condition on Γ2 and Dirichlet conditions on Γ1. On the boundary with the Neumann condition, the number of selected neighboring point is set to be ns=9. Tables 8 and 9 give the found results found.Example 3.
The investigated problem is the widely used convection-diffusion flow problem of a Gaussian pulse given by the equation:(35)∂u∂tx,y,t+vx,y,t∇ux,y,t−νΔux,y,t=fx,y,t∀x,y,t∈Ω×0,T,where vx,y,t=−4y,4xT is the velocity vector, and ν is a constant. Initial and boundary conditions are taken according to the analytical solution:
(36)ux,y,t=σ2σ2+4νtexp−x−0.52+y−0.52σ2+4νt.Table 8
Results for different values ofε in Example 2 with Nt=10. Uniform spaced nodes (left) and circular nodes (right), with Neumann conditions.
εMAERMSELer1MAERMSELer10.11.72E-033.29E-044.53E-041.19E-032.07E-042.92E-040.21.51E-032.71E-043.59E-041.43E-031.91E-042.33E-040.31.42E-032.39E-042.95E-041.85E-033.09E-043.91E-040.41.85E-033.67E-044.88E-042.46E-035.87E-048.68E-040.52.44E-036.59E-049.97E-043.29E-039.97E-041.55E-030.63.23E-031.08E-031.71E-034.37E-031.54E-032.45E-030.74.23E-031.63E-032.63E-035.76E-032.24E-033.58E-03Table 9
Results for different values ofε in Example 2 with Nt=20. Uniform spaced nodes (left) and circular nodes (right), with Neumann conditions.
εMAERMSELer1MAERMSELer10.11.61E-033.93E-045.79E-048.84E-042.48E-043.81E-040.21.39E-033.21E-044.64E-041.11E-031.69E-042.43E-040.31.08E-032.30E-043.28E-041.51E-031.89E-042.24E-040.41.47E-032.61E-043.25E-042.09E-034.39E-046.63E-040.52.03E-035.11E-047.80E-042.88E-038.30E-041.33E-030.62.76E-039.03E-041.47E-033.90E-031.35E-032.18E-030.73.70E-031.41E-032.32E-035.20E-031.99E-033.21E-03In our simulation, we takeσ2=0.001 and ν=0.01. The computational space-time domain ΩT is 0,1×0,1×0,π/2. In Figure 8, we show the numerical solution at final time t=T face to the analytical solution, the total number of nodes is N=303. Table 10 displays the absolute maximum error obtained for different values of N; it can be remarked that the rate of convergence is near quadric. In all simulations for Example 3, the shape parameter is chosen to be ε=5.Figure 8
Exact (a) and numerical (b) solutions att=π/2 for the Example 3 with N=303 and ε=5.
(a)(b)Table 10
Maximum absolute error and rate of convergence att=T for different values of N in Example 3 using ε=5 and ΩT=0,1×0,1×0,π/2.
NMAEhROC1536.70E-040.0714—2033.74E-040.05261.912532.35E-040.04171.993031.62E-040.03451.96


## 3.3. Two-Dimensional Hyperbolic Equations

Although it has been shown that the technique gives good results solving the parabolic equation with different boundary conditions, the hyperbolic equation is still an important and interesting equation to be tested using the described methodology. In all simulations, the problem is treated as an ill-posed one as it has mentioned before and the algebraic matrix obtained is square.Example 4.
Considering a two-dimensional wave problem defined by the following equation on the domainΩ×0,T with Ω=0,12(37)∂2u∂t2x,y,t−aΔux,y,t=fx,y,t∀x,y,t∈Ω×0,T,with two initial conditions ux,0=u0x and ∂u/∂tx,y,0=hx,y and Dirichlet boundary condition on ∂Ω.Takingax,y,t=xye−2t, the functions u0, h, and f are chosen according to the analytical solution given by ux,t=sinxsinysint.Applying the “ill-posed-problem” technique [22], we obtain accurate results by setting the number of neighboring points to be ns=9 in the entire domain. Table 11 shows results obtained of different errors, taking ε=0.1 and using the space-time domain 0,13. Figure 9 shows that the MAE error decreases when increasing the number of nodes.Table 11
Errors for different values ofNx, Ny, and Nt of Example 4 using ε=0.1 and ΩT=0,13.
NxNyNtMAERMSELer15551.80E-035.13E-042.69E-031010104.59E-041.26E-048.08E-041515151.91E-045.23E-053.53E-042020201.06E-042.86E-051.97E-042525256.68E-051.81E-051.25E-04Figure 9
Maximum absolute errors vs.Nx with Nx=Ny=Nt of Example 4 using ε=0.1.As in the example of the one-dimensional problem, the convergence rate of the technique using different uniformly distributed sources points andε=0.1 is given in Table 12 showing that the order of convergence is quadratic.Table 12
Rate of convergence for wave Example4 using ε=0.1.
hMAEROCRMSEROCLer1ROC0.250001.70E-03—4.61E-04—2.44E-03—0.111114.46E-041.64601.22E-041.63637.86E-041.39700.071431.78E-042.07964.92E-052.05963.34E-041.93650.052631.02E-041.81902.80E-051.85161.94E-041.78330.041676.43E-051.98351.78E-051.93831.25E-041.8771Att=1, the maximum absolute error is MAE=6.68×10−5, and the root mean square error is RMSE=1.81×10−5; these values are obtained with ε=0.1 and N=253. We remark that the maximum error is obtained at t=T (see Figure 10).Figure 10
Absolute errors withNt=25 and ε=0.1 at t=1, for Example 4.


## 3.4. Three-Dimensional Parabolic and Hyperbolic Examples

In this section, we show the robustness and the accuracy of the proposed method implemented in high dimension problems. The traditional three-dimensional problems are transformed into 4D problems.Example 5.
The problem investigated is the 3D heat transfer problem given by the equation:(38)∂u∂tx,y,z,t−ax,y,z,tΔux,y,z,t=fx,y,z,t.Withax,y,z,t=xyzet, the initial and boundary conditions are chosen in such a way that the analytical solution is given by the function ux,y,z,t=e−tsinxsinysinz. The space-time domain is chosen Ω=0,14 and the right-hand side is fx,y,z,t=sinxsinysinze−t−1+3xyzet.Takingns=9, Table 13 shows the errors found in two simulations. The total number of nodes is N=Nx×Nt×Ny×Nz=54 and N=94.Table 13
Errors for Example5 using two different values of N.
NεMAERMSELer1540.94.88E-49.28E-55.45E-4940.82.23E-43.62E-52.75E-4We can remark that good accuracy is obtained even with weak number of nodes. We note that similar accuracy is found using different values of the shape parameterε in the interval 0,1.5Example 6.
The considered hyperbolic 3D problem is given by(39)∂2u∂t2x,y,z,t−ax,y,z,y,tΔux,y,z,t=fx,y,z,t.We useax,y,z,t=xyze−t and ε=0.12. The computational space-time domain is Ω=0,14. The boundary and initial conditions and the function f are chosen according to the solution ux,y,z,t=sinxsinysinzsint. In this example, the number of neighboring points is set to ns=11 on the entire domain. The total number of nodes is N=Nx×Nt×Ny×Nz=94. The errors obtained for this example are MAE=1.54×10−3, RMSE=9.23×10−5, and Ler1=6.00×10−4. The same ill-posed technique is used with the squared algebraic linear system.


## 4. Conclusion

The local space-time RBF collocation method developed in [13] is extended to solve parabolic and hyperbolic equations with variable coefficients in the space-time domain. Following the same technique used in [13], the parabolic problem is solved using the governing equation as boundary condition on the boundary characterized by the final time T, and the hyperbolic equation is solved as an ill-posed problem with incomplete data boundary condition on Ω×t=T. The first one-dimensional tested numerical simulation concerns the comparison of the 4th order Runge-Kutta method to the presented one. It has been shown that although the Runge-Kutta method did not need the inversion of the matrix and it is considered as an explicit technique, it needs more CPU time than the proposed one.The numerical results obtained for other two- and three-dimensional parabolic and hyperbolic problems show the performance of the technique and its accuracy for solving different examples. It has been demonstrated that the technique is simple, straightforward. The formulation of the hyperbolic equation as an ill-posed problem and also its application to three-dimensional problems and quasisingular problems with variables coefficients are the most important ideas of the paper. We can also remark that the main advantages of the technique which are(1)
The time stability analysis is not discussed as it is the case for other time-stepping algorithms as implicit, explicit,θ-method, etc.(2)
Reducing the computational time as there is no need to recompute the matrix for the resulting algebraic system at each time level, unlike the case for other time integration methods used to solve PDEs with time-dependent coefficients and(3)
Solving hyperbolic partial differential equations with variable coefficients as an inverse problem, since just the hyperbolic equation with constant coefficients is less discussed and needs more sophisticated time integration techniqueare more representative for partial differential equation with variable coefficients than with constant one. Further work will focus on the application of the developed technique to nonlinear equations. Stability analysis of the technique should also be investigated in further work.


## 5. Further Work

Further work will focus on applying the method to more real PDEs like the equation of the telegraph and the transient heat conduction problem in an anisotropic medium, all this PDEs will be investigated with variable coefficients.



---

*Source: 6688806-2021-02-08.xml*