1.
Introduction
Nonlinear evolution equations (NEEs) are mathematical equations that describe the behavior of complex physical systems. These equations are often used in various fields of science, including material and biological sciences, physics, mathematics, and mechanics [1,2,3,4,5]. One of the most commonly used NEE is the nonlinear Schrödinger equation, which is used in quantum mechanics to describe the behavior of wave packets [6,7]. Another example is the Klein-Gordon equation, which is used in particle physics to describe the behavior of spinless particles. The Cahn-Hilliard equation is another important equation used in materials science to describe the phase separation behaviour of binary mixtures. This equation has applications in the fields of metallurgy, polymer science, and biomaterials. The Navier-Stokes equation is a well-known NEE that describes the behavior of fluids, and is used in fluid mechanics to study the flow of fluids in various applications, including aerospace engineering, civil engineering, and oceanography [8,9,10,11].
NEEs are important in many fields of science because they can provide a mathematical model to describe the behavior of complex physical systems. By studying these equations and developing new numerical techniques to solve them, scientists can gain a better understanding of the physical phenomena they describe, and develop new technologies and materials that are based on this knowledge. To study the behavior of the physical problem (see [12,13,14,15]), either analytical or numerical solutions are needed. However, finding an analytical solution for NNEs is not always simple. As a result, developing numerical methods for such problems is useful. In this article, we investigate the nonlinear Klein-Gordon and sinh-Gordon equations, which are given as follows:
where u=u(x,t), α is constant, and ψ is some nonlinear function. For
the Eq (1.1) becomes the sinh-Gordon equation [16]. When
the Eq (1.1) becomes the Klein-Gordon equation, which is given as follows:
with the following initial and boundary conditions:
The Klein-Gordon equation has numerous applications in different areas of physics and engineering and plays a crucial role in our understanding of the behavior of particles and fields in the universe. The Klein-Gordon equation is a relativistic wave equation that describes the behavior of particles with a zero spin, such as the Higgs boson particle. It has several applications in different fields of physics; for example, the Klein-Gordon equation is fundamental in the quantum field theory. It describes the behavior of scalar fields, which are fields that have a single value at each point in space-time [17,18]. Additionally, it can be used to describe the behavior of phonons, which are the quantized vibrations of a crystal lattice, and is used to describe the behavior of the inflaton field, which is thought to be responsible for the rapid expansion of the universe during the Big Bang [19]. Due to the wide range of applications, numerous researchers have investigated the solution of this equation utilizing various methodologies, including collocation points with the thin-plate splines radial basis function [20], the sine-cosine and tanh methods for travelling wave solutions [19], and the generalized auxiliary equation methods in [21]. Pseduospectral formulation with a perfectly matched layer has been suggested for the non-relativistic Klein-Gordon equation in [22]. Three different numerical methods were introduced to solve the rotating Klein-Gordon equation from relativistic regimes to non-relativistic regimes in polar and Cartesian coordinates [23]. Various numerical methods, including finite difference time domain methods, the time-splitting method, an exponential wave integrator, a limit integrator, a multiscale time integrator, the two-scale formulation method, and an iterative exponential integrator, were compared in terms of their accuracy while solving the Klein-Gordon equation (see [24]). The authors used a finite element technique based on cubic B-splines, the decomposition method [25], the fourth-order implicit Runge-Kutta-Nystrom, the fourth-order compact finite difference method [26], the boundary integral equation method [27], and the cubic B-spline collocation approach [28] to solve the generalized form of the nonlinear Klein-Gordon equation. For the solution of this model, a variety of numerical techniques have been reported in [29,30]. The two-dimensional sinh-Gordon equation is an essential partial differential equation (PDE) model with applications in a variety of fields.
The following is a concise review of the sinh-Gordon equation. In [31], the (G'/G)-expansion scheme was proposed for the generalized form of the double sinh-Gordon equations. For the analytical solutions to the sinh-Gordon equation, the tanh method was used in [32]. The element-free Galerkin technique was suggested in [33] to solve the generalized sinh-Gordon equation. The radial basis function was utilized by the authors in [34] to solve the sinh-Gordon equation. In [35], the author employed an implicit Lie-group iterative technique for the solution of the sinh-Gordon equation, whereas the Pseudo spectral and Kansa's methods based on radial basis functions were used in [36] to solve the two-dimensional sinh-Gordon equation. In [37], two effective approaches, namely the generalized leap-frog and the method of lines, were utilized for the sinh-Gordon equation, whereas the polynomial differential quadrature method was examined in [38].
In this work, we utilize the Fibonacci polynomials combined with Störmer's method to compute the numerical solutions of the one-dimensional Klein-Gordon and the two-dimensional sinh-Gordon equations. One of the key advantages of the suggested method is the ease implementation in higher-order derivatives using the relation of the Lucas and Fibonacci polynomials. Furthermore, the proposed approach improves the accuracy even for a small number of collocation points, thus lowering the computing costs. These polynomials have a wide range of applications in differential equations; for instance, the relationship between Chebyshev and Lucas polynomials was discussed in [39], and accurate solutions for the boundary value problems were obtained. Higher-order differential equations were solved using the Lucas polynomial technique in [40], while the Volterra-Fredholm integral differential equations were solved using a Fibonacci polynomial approach in [41]. Delay difference equations were solved using a hybrid Taylor-Lucas polynomial approach [42]. The Fibonacci wavelet approach was implemented for non-linear reaction diffusion equations of a Ficher-type, non-linear Hunter-Saxton and time fractional telegraph equation in [43,44]. The author and his co-worker used the Gegenbauer wavelet method for time fractional PDEs in [45,46]. For the first time, the authors proposed a hybrid Lucas and Fibonacci polynomial scheme to solve the time-dependent PDEs in [16,47]. More publications using Lucas polynomials and finite-differences to find efficient numerical solutions for different types of PDE models can be found in [48,49,50,51].
The rest of the paper is organized as follows: in Section 2, fundamental concepts and definitions are discussed; in Section 3, the proposed methodology of the scheme is described; the method is verified with the help of numerical experiment in Section 4, and finally, the concluding remarks are drawn at the end.
2.
Fundamental concepts
In this section, we will go over several key definitions and properties related to the Fibonacci and Lucas polynomials.
2.1. Fibonacci polynomials
The Fibonacci polynomial is an extension of Fibonacci numbers described by the linear recurrence relation, which is shown as follows [52]:
with the initial values
Their explicit form is as follows:
where
and ⌊k⌋ denotes the integer floor function. Equation (2.1) generates a sequence for x=1 of the Fibonacci numbers.
2.2. Lucas polynomials
The Lucas polynomials can be defined as follows [52]:
with
Their explicit form is
when x=1. Equation (2.2) generates a sequence of the Lucas numbers.
2.3. Function approximation
Let u(x) be a continuous function that may be approximated in terms of the Lucas series in the following manner:
where L(x) are the Lucas polynomials and λk are the unknown coefficients. Utilizing the finite terms of the Lucas series, the mth-order derivative of a function u(x) can be approximated as follows:
Utilizing the differentiation matrix D and Fibonacci polynomials, L(m)k(x) can be approximated [52] as follows:
where
and d is a N×N matrix
3.
Proposed methodology
To describe the proposed hybrid scheme, we consider Eq (1.1)
with the following conditions:
First, re-arrange Eq (3.1)
Now, using Störmer's method for the time discretization of Eq (3.4), we have the following:
where
dt is time stepsize for t∈[0,T].
Then, let us discretize the spacial domain [a,b] into the M+1 number of the nodes xi,i=0,1,2,…,M, where x0,xM are the boundary points while the remaining (i.e., xi,i=1,2,…,M−1) are interior points of the domain, which are computed as follows:
where
is spatial step size. Here, we approximate the unknown function u(x) by Lucas polynomials at the nth time level, which is denoted as un(x):
where C is the vector of unknown coefficients dependent on time, and L is the (M+1)×1 vector of Lucas polynomials, such that
Now, for the collocation points xi,
Let
then, in the matrix notation,
where
and
Similarly, the kth order derivative of unknown functions are approximated by the kth order derivatives of Lucas polynomials
In the matrix form,
where
and
Using the above approximation in Eq (3.5), we have the following:
where
and
For n = 1, u1 and u0 can be obtained from the initial conditions. To satisfy the boundary conditions we collocated the boundary conditions that corresponded to each of the following:
That is, for
we then update the solution vector at boundary points to have
Remark 3.1. The stability condition for the second order Störmer method is as follows:
where α is the coefficient of the Laplace operator. Then the method is then said to be stable, by the rule of thumb, if all the eigenvalues, λ, of the discretized operator,
scaled by dt2 lies within the stability region of the time-stepping method [53], i.e.,
4.
Numerical results and discussion
This section demonstrates the numerical solution of one dimension Klien-Gordon and two dimension sinh-Gordon equations utilizing the proposed scheme. All the computations are performed using MATLAB (R2012a) on Dell PC Laptop with an Intel(R) Core(TM)i5-2450M CPU 2.50 GHz 2.50 GHz 8 GB RAM. To check the accuracy and convergence, we used the following error norms:
where E is the exact solution.
4.1. (1+1)-dimensional case
Test Problem 4.1. Consider Eq (1.2) with α=β=−1, γ=0, g(x,t)=0; we have
with the following initial and boundary conditions:
and the actual solution is as follows:
Using the proposed methodology for this particular equation the iterative scheme (3.14) becomes the following:
The above equation is used to generate numerical results, which are provided in Table 1 in terms of the L∞ and L2 error norms. These findings are obtained using nodal points M=10 and various values of the time step size dt while preserving the final time T=0.1, and for various values of T and fixing dt=0.0001. The accuracy increases as the value of dt decreases, as shown in this table. Based on the comparison with [29], it is obvious that the suggested scheme is more accurate than those reported in the recent literature. The numerical results in terms of L∞, L2, and Lrms, as well as central processing unit (CPU) times in seconds, are computed using the suggested approach for various values of M in Table 2. The table illustrates that as the number of collocation points rises, the current approach provides an improved accuracy. Furthermore, the proposed technique may be seen to be quite efficient, whereas the comparison between the exact and approximate solutions are visualized in Figure 1.
Test Problem 4.2. Now assume
and
in Eq (1.2) with the following conditions:
The exact solution of Eq (1.2) for the above parameter is given in [29], which is
In Test Problem 4.2, the numerical results are computed utilizing the proposed technique. These results are shown in Table 3 using nodal points
and various values of dt while preserving the final time
and for various values of T and fixing
In this situation, as indicated in the table, the method's accuracy improves as the value of dt declines. Furthermore, in comparison with [29], we came to know that the current method is more accurate than other articles published in the recent literature. Similarly to the previous test problem, the numerical results are computed in the form of L∞, L2, and Lrms, as well as the CPU times in seconds, using the suggested approach for various values of M in Table 4. As the number of collocation points increases, the current approach efficiently improves the accuracy, as seen in the table whereas profiles of the exact and approximate solutions are given in Figure 2.
Test Problem 4.3. Consider Eq (1.2) with the exact solution u(x,t)=x3t3, when
and
with associated conditions
Tables 5 and 6 show the numerical results for Test Problem 4.3 for various values of dt and M, respectively. We compared our results with the methods in [29] for
in Table 5. We can see from this table that the method's accuracy improves as the number of iterations increases, with a very good accuracy of 10−8 for
In comparison to the method in [29], the proposed method produces more accurate results. In Table 6, the results are computed in terms of L∞, L2, and Lrms for
and time T=1, as well as the CPU time in seconds. The suggested procedure is accurate and efficient, as shown in this table. The profile of the exact and approximate solutions are shown in Figures 3 and 4, respectively.
Test Problem 4.4. Consider
and
in Eq (1.2) with the following conditions:
The exact solution of Eq (1.2) for the above parameter given as follows:
The solution is computed for the fixed point T=0.1 and M=20, while decreasing the time step. The results are reported in Table 7, which shows that the accuracy improves as the number of iteration increases. In Table 8, the results are recorded for various numbers of nodal points M, and notice that as number of collocation point increases, the numerical convergence accuracy improves. Along with, maximum error, the CPU time in seconds is also reported in the tables. The solution profiles of the exact and approximate solutions are illustrated in Figure 5, which shows that both the solutions are well matched.
Test Problem 4.5. (1+2)-dimensional case: Finally, consider the two dimensional sinh-Gordon equation
with initial conditions
and boundary conditions
The exact solution
is considered, as given in [47]. Additionally the source term is as follows:
The results of the suggested method for Test Problem 4.5 are shown in Table 9 for various values of dt with M=4 and T=1. This table shows that the method delivers a good accuracy for very coarse nodes, and that the method's accuracy improves to some extent as the number of iterations increases. Additionally in this table, we also compared the computed results to the approach in [47]. The results of the proposed approach in contrast to the method presented in [47] are shown in Table 10 for various values of M for the same test problem. This comparison reveals that the proposed approach is more accurate than [47]. The approximate solution of the suggested approach is compared to the exact solution in Figure 6. Although the absolute errors are indicated in Figure 7, the numerical solution is in a good agreement with the exact solution, as demonstrated in Figure 6.
5.
Conclusions
A numerical method based on Fibonacci polynomials combined with Störmer's method was developed to solve the nonlinear Klein/Sinh-Gordon equation. On a variety of test cases, numerical evaluation exhibited an improved performance of the proposed scheme against the available results in the literature. As evidenced by the comparison, the proposed technique is simple to execute and faster in convergence. Furthermore, it was observed that by reducing the spatial and time steps size, the errors were reduced. With slight modifications, the proposed method can be applied to a variety of PDEs. This modification is mainly dictated by the problem at hand with the associated initial and boundary conditions that are dictated by the particular problem of interest. As such, we are interested in applying the proposed method to time-dependent Neumann and/or mixed boundary value problems with possible discontinuous initial condition(s). The limitations, as is the case in spectral methods, was the high-condition numbers of the collocated/system matrices. However, this can be alleviated by the state-of-the-art pre-conditioned algorithms developed for spectral methods. In the future, the proposed method can be expanded for applications to complex fractal-fractional problems.
Conflict of interest
The authors declare that there are no conflicts of interest in this paper.
DownLoad:
