This article studies a posteriori error analysis of fully discrete finite element approximations for semilinear parabolic optimal control problems. Based on elliptic reconstruction approach introduced earlier by Makridakis and Nochetto [25], a residual based a posteriori error estimators for the state, co-state and control variables are derived. The space discretization of the state and co-state variables is done by using the piecewise linear and continuous finite elements, whereas the piecewise constant functions are employed for the control variable. The temporal discretization is based on the backward Euler method. We derive a posteriori error estimates for the state, co-state and control variables in the L^∞(0,T;L²(Ω))-norm. Finally, a numerical experiment is performed to illustrate the performance of the derived estimators.

We consider numerical methods to solve the Allen-Cahn equation using the secondorder Crank-Nicolson scheme in time and the second-order central difference approach in space. The existence of the finite difference solution is proved with the help of Browder fixed point theorem. The difference scheme is showed to be unconditionally convergent in L_∞ norm by constructing an auxiliary Lipschitz continuous function. Based on this result, it is demonstrated that the difference scheme preserves the maximum principle without any restrictions on spatial step size and temporal step size. The numerical experiments also verify the reliability of the method.

In this paper, we study a nonlinear first-order singularly perturbed Volterra integrodifferential equation with delay. This equation is discretized by the backward Euler for differential part and the composite numerical quadrature formula for integral part for which both an a priori and an a posteriori error analysis in the maximum norm are derived. Based on the a priori error bound and mesh equidistribution principle, we prove that there exists a mesh gives optimal first order convergence which is robust with respect to the perturbation parameter. The a posteriori error bound is used to choose a suitable monitor function and design a corresponding adaptive grid generation algorithm. Furthermore, we extend our presented adaptive grid algorithm to a class of second-order nonlinear singularly perturbed delay differential equations. Numerical results are provided to demonstrate the effectiveness of our presented monitor function. Meanwhile, it is shown that the standard arc-length monitor function is unsuitable for this type of singularly perturbed delay differential equations with a turning point.

The present paper conjectures a topological condition which classifies a T-spline into standard, semi-standard and non-standard. We also provide the basic framework to prove the conjecture on the classification of semi-standard T-splines and give the proof for the semi-standard of bi-degree (1, d) and (d, 1) T-splines.

This paper focuses on a fast and high-order finite difference method for two-dimensional space-fractional complex Ginzburg-Landau equations. We firstly establish a three-level finite difference scheme for the time variable followed by the linearized technique of the nonlinear term. Then the fourth-order compact finite difference method is employed to discretize the spatial variables. Hence the accuracy of the discretization is $\mathcal{O}$(τ² + $h_1^4$ + $h_2^4$) in L₂-norm, where τ is the temporal step-size, both h₁ and h₂ denote spatial mesh sizes in x- and y- directions, respectively. The rigorous theoretical analysis, including the uniqueness, the almost unconditional stability, and the convergence, is studied via the energy argument. Practically, the discretized system holds the block Toeplitz structure. Therefore, the coefficient Toeplitz-like matrix only requires $\mathcal{O}$(M₁M₂) memory storage, and the matrix-vector multiplication can be carried out in $\mathcal{O}$(M₁M₂(log M₁ + log M₂)) computational complexity by the fast Fourier transformation, where M₁ and M₂ denote the numbers of the spatial grids in two different directions. In order to solve the resulting Toeplitz-like system quickly, an efficient preconditioner with the Krylov subspace method is proposed to speed up the iteration rate. Numerical results are given to demonstrate the well performance of the proposed method.

In this paper, the numerical methods for semi-linear stochastic delay integro-differential equations are studied. The uniqueness, existence and stability of analytic solutions of semi-linear stochastic delay integro-differential equations are studied and some suitable conditions for the mean-square stability of the analytic solutions are also obtained. Then the numerical approximation of exponential Euler method for semi-linear stochastic delay integro-differential equations is constructed and the convergence and the stability of the numerical method are studied. It is proved that the exponential Euler method is convergent with strong order 1/2 and can keep the mean-square exponential stability of the analytical solutions under some restrictions on the step size. In addition, numerical experiments are presented to confirm the theoretical results.

In this paper, we present a unified finite volume method preserving discrete maximum principle (DMP) for the conjugate heat transfer problems with general interface conditions. We prove the existence of the numerical solution and the DMP-preserving property. Numerical experiments show that the nonlinear iteration numbers of the scheme in [24] increase rapidly when the interfacial coefficients decrease to zero. In contrast, the nonlinear iteration numbers of the unified scheme do not increase when the interfacial coefficients decrease to zero, which reveals that the unified scheme is more robust than the scheme in [24]. The accuracy and DMP-preserving property of the scheme are also verified in the numerical experiments.

In this article, we consider the primal-dual path-following method and the trust-region updating strategy for the standard linear programming problem. For the rank-deficient problem with the small noisy data, we also give the preprocessing method based on the QR decomposition with column pivoting. Then, we prove the global convergence of the new method when the initial point is strictly primal-dual feasible. Finally, for some rankdeficient problems with or without the small noisy data from the NETLIB collection, we compare it with other two popular interior-point methods, i.e. the subroutine pathfollow.m and the built-in subroutine linprog.m of the MATLAB environment. Numerical results show that the new method is more robust than the other two methods for the rank-deficient problem with the small noise data.

This paper deals with discontinuous dual reciprocity boundary element method for solving an inverse source problem. The aim of this work is to determine the source term in elliptic equations for nonhomogenous anisotropic media, where some additional boundary measurements are required. An equivalent formulation to the primary inverse problem is established based on the minimization of a functional cost, where a regularization term is employed to eliminate the oscillations of the noisy data. Moreover, an efficient algorithm is presented and tested for some numerical examples.

A nonlinear fully implicit finite difference scheme with second-order time evolution for nonlinear diffusion problem is studied. The scheme is constructed with two-layer coupled discretization (TLCD) at each time step. It does not stir numerical oscillation, while permits large time step length, and produces more accurate numerical solutions than the other two well-known second-order time evolution nonlinear schemes, the Crank-Nicolson (CN) scheme and the backward difference formula second-order (BDF2) scheme. By developing a new reasoning technique, we overcome the difficulties caused by the coupled nonlinear discrete diffusion operators at different time layers, and prove rigorously the TLCD scheme is uniquely solvable, unconditionally stable, and has second-order convergence in both space and time. Numerical tests verify the theoretical results, and illustrate its superiority over the CN and BDF2 schemes.

In this paper, we derive a residual based a posteriori error estimator for a modified weak Galerkin formulation of second order elliptic problems. We prove that the error estimator used for interior penalty discontinuous Galerkin methods still gives both upper and lower bounds for the modified weak Galerkin method, though they have essentially different bilinear forms. More precisely, we prove its reliability and efficiency for the actual error measured in the standard DG norm. We further provide an improved a priori error estimate under minimal regularity assumptions on the exact solution. Numerical results are presented to verify the theoretical analysis.

In this paper, we consider the knot placement problem in B-spline curve approximation. A novel two-stage framework is proposed for addressing this problem. In the first step, the l_∞,1-norm model is introduced for the sparse selection of candidate knots from an initial knot vector. By this step, the knot number is determined. In the second step, knot positions are formulated into a nonlinear optimization problem and optimized by a global optimization algorithm — the differential evolution algorithm (DE). The candidate knots selected in the first step are served for initial values of the DE algorithm. Since the candidate knots provide a good guess of knot positions, the DE algorithm can quickly converge. One advantage of the proposed algorithm is that the knot number and knot positions are determined automatically. Compared with the current existing algorithms, the proposed algorithm finds approximations with smaller fitting error when the knot number is fixed in advance. Furthermore, the proposed algorithm is robust to noisy data and can handle with few data points. We illustrate with some examples and applications.

We present a stochastic trust-region model-based framework in which its radius is related to the probabilistic models. Especially, we propose a specific algorithm termed STRME, in which the trust-region radius depends linearly on the gradient used to define the latest model. The complexity results of the STRME method in nonconvex, convex and strongly convex settings are presented, which match those of the existing algorithms based on probabilistic properties. In addition, several numerical experiments are carried out to reveal the benefits of the proposed methods compared to the existing stochastic trust-region methods and other relevant stochastic gradient methods.

This paper presents a second-order direct arbitrary Lagrangian Eulerian (ALE) method for compressible flow in two-dimensional cylindrical geometry. This algorithm has half-face fluxes and a nodal velocity solver, which can ensure the compatibility between edge fluxes and the nodal flow intrinsically. In two-dimensional cylindrical geometry, the control volume scheme and the area-weighted scheme are used respectively, which are distinguished by the discretizations for the source term in the momentum equation. The two-dimensional second-order extensions of these schemes are constructed by employing the monotone upwind scheme of conservation law (MUSCL) on unstructured meshes. Numerical results are provided to assess the robustness and accuracy of these new schemes.

In this paper, ETD3-Padé and ETD4-Padé Galerkin finite element methods are proposed and analyzed for nonlinear delayed convection-diffusion-reaction equations with Dirichlet boundary conditions. An ETD-based RK is used for time integration of the corresponding equation. To overcome a well-known difficulty of numerical instability associated with the computation of the exponential operator, the Padé approach is used for such an exponential operator approximation, which in turn leads to the corresponding ETD-Padé schemes. An unconditional L² numerical stability is proved for the proposed numerical schemes, under a global Lipshitz continuity assumption. In addition, optimal rate error estimates are provided, which gives the convergence order of ^O(k³ + h^r) (ETD3- Padé) or O(k⁴ + h^r) (ETD4-Padé) in the L² norm, respectively. Numerical experiments are presented to demonstrate the robustness of the proposed numerical schemes.

We study spatially semidiscrete and fully discrete two-scale composite finite element method for approximations of the nonlinear parabolic equations with homogeneous Dirichlet boundary conditions in a convex polygonal domain in the plane. This new class of finite elements, which is called composite finite elements, was first introduced by Hackbusch and Sauter [Numer. Math., 75 (1997), pp. 447-472] for the approximation of partial differential equations on domains with complicated geometry. The aim of this paper is to introduce an efficient numerical method which gives a lower dimensional approach for solving partial differential equations by domain discretization method. The composite finite element method introduces two-scale grid for discretization of the domain, the coarse-scale and the fine-scale grid with the degrees of freedom lies on the coarse-scale grid only. While the fine-scale grid is used to resolve the Dirichlet boundary condition, the dimension of the finite element space depends only on the coarse-scale grid. As a consequence, the resulting linear system will have a fewer number of unknowns. A continuous, piecewise linear composite finite element space is employed for the space discretization whereas the time discretization is based on both the backward Euler and the Crank-Nicolson methods. We have derived the error estimates in the L^∞(L²)-norm for both semidiscrete and fully discrete schemes. Moreover, numerical simulations show that the proposed method is an efficient method to provide a good approximate solution.

In this paper, we present a two-grid discretization scheme for semilinear parabolic integro-differential equations by H¹-Galerkin mixed finite element methods. We use the lowest order Raviart-Thomas mixed finite elements and continuous linear finite element for spatial discretization, and backward Euler scheme for temporal discretization. Firstly, a priori error estimates and some superclose properties are derived. Secondly, a two-grid scheme is presented and its convergence is discussed. In the proposed two-grid scheme, the solution of the nonlinear system on a fine grid is reduced to the solution of the nonlinear system on a much coarser grid and the solution of two symmetric and positive definite linear algebraic equations on the fine grid and the resulting solution still maintains optimal accuracy. Finally, a numerical experiment is implemented to verify theoretical results of the proposed scheme. The theoretical and numerical results show that the two-grid method achieves the same convergence property as the one-grid method with the choice h = H ².

Implicit-explicit Runge-Kutta-Rosenbrock methods are proposed to solve nonlinear stiff ordinary differential equations by combining linearly implicit Rosenbrock methods with explicit Runge-Kutta methods. First, the general order conditions up to order 3 are obtained. Then, for the nonlinear stiff initial-value problems satisfying the one-sided Lipschitz condition and a class of singularly perturbed initial-value problems, the corresponding errors of the implicit-explicit methods are analysed. At last, some numerical examples are given to verify the validity of the obtained theoretical results and the effectiveness of the methods.

We propose a θ-L approach for solving a sharp-interface model about simulating solidstate dewetting of thin films with isotropic/weakly anisotropic surface energies. The sharpinterface model is governed by surface diffusion and contact line migration. For solving the model, traditional numerical methods usually suffer from the severe stability constraint and/or the mesh distribution trouble. In the θ-L approach, we introduce a useful tangential velocity along the evolving interface and utilize a new set of variables (i.e., the tangential angle θ and the total length L of the interface curve), so that it not only could reduce the stiffness resulted from the surface tension, but also could ensure the mesh equidistribution property during the evolution. Furthermore, it can achieve second-order accuracy when implemented by a semi-implicit linear finite element method. Numerical results are reported to demonstrate that the proposed θ-L approach is efficient and accurate.

In this paper, matrix representations of the best spline quasi-interpolating operator over triangular sub-domains in S₂¹(Δ_mn⁽²⁾), and coefficients of splines in terms of B-net are calculated firstly. Moreover, by means of coefficients in terms of B-net, computation of bivariate numerical cubature over triangular sub-domains with respect to variables x and y is transferred into summation of coefficients of splines in terms of B-net. Thus concise bivariate cubature formulas are constructed over rectangular sub-domain. Furthermore, by means of module of continuity and max-norms, error estimates for cubature formulas are derived over both sub-domains and the domain.

A combined scheme of the improved two-grid technique with the block-centered finite difference method is constructed and analyzed to solve the nonlinear time-fractional parabolic equation. This method is considered where the nonlinear problem is solved only on a coarse grid of size H and two linear problems based on the coarse-grid solutions and one Newton iteration is considered on a fine grid of size h. We provide the rigorous error estimate, which demonstrates that our scheme converges with order O(Δt^2-α+h²+H⁴) on non-uniform rectangular grid. This result indicates that the improved two-grid method can obtain asymptotically optimal approximation as long as the mesh sizes satisfy h=O(H²). Finally, numerical tests confirm the theoretical results of the presented method.

This is one of our series works on numerical methods for mean-field forward backward stochastic differential equations (MFBSDEs). In this work, we propose an explicit multistep scheme for MFBSDEs which is easy to implement, and is of high order rate of convergence. Rigorous error estimates of the proposed multistep scheme are presented. Numerical experiments are carried out to show the efficiency and accuracy of the proposed scheme.

We consider optimal shape design in Stokes flow using H¹ shape gradient flows based on the distributed Eulerian derivatives. MINI element is used for discretizations of Stokes equation and Galerkin finite element is used for discretizations of distributed and boundary H¹ shape gradient flows. Convergence analysis with a priori error estimates is provided under general and different regularity assumptions. We investigate the performances of shape gradient descent algorithms for energy dissipation minimization and obstacle flow. Numerical comparisons in 2D and 3D show that the distributed H¹ shape gradient flow is more accurate than the popular boundary type. The corresponding distributed shape gradient algorithm is more effective.

We analyze the convergence of the weighted nonlocal Laplacian (WNLL) on the high dimensional randomly distributed point cloud. Our analysis reveals the importance of the scaling weight, μ~|P|/|S|with|P|and|S|being the number of entire and labeled data, respectively, in WNLL. The established result gives a theoretical foundation of the WNLL for high dimensional data interpolation.

Precision matrix estimation is an important problem in statistical data analysis. This paper proposes a sparse precision matrix estimation approach, based on CLIME estimator and an efficient algorithm GISS^ρ that was originally proposed for l₁ sparse signal recovery in compressed sensing. The asymptotic convergence rate for sparse precision matrix estimation is analyzed with respect to the new stopping criteria of the proposed GISS^ρ algorithm. Finally, numerical comparison of GISS^ρ with other sparse recovery algorithms, such as ADMM and HTP in three settings of precision matrix estimation is provided and the numerical results show the advantages of the proposed algorithm.

This note introduces a method for sampling Ising models with mixed boundary conditions. As an application of annealed importance sampling and the Swendsen-Wang algorithm, the method adopts a sequence of intermediate distributions that keeps the temperature fixed but turns on the boundary condition gradually. The numerical results show that the variance of the sample weights is relatively small.

We prove a theorem concerning the approximation of generalized bandlimited multivariate functions by deep ReLU networks for which the curse of the dimensionality is overcome. Our theorem is based on a result by Maurey and on the ability of deep ReLU networks to approximate Chebyshev polynomials and analytic functions efficiently.

In this paper, we consider a modified alternating positive semidefinite splitting preconditioner for solving the saddle point problems arising from the finite element discretization of the hybrid formulation of the time-harmonic eddy current model. The eigenvalue distribution and an upper bound of the degree of the minimal polynomial of the preconditioned matrix are studied for both simple and general topology. Numerical results demonstrate the effectiveness of the proposed preconditioner when it is used to accelerate the convergence rate of Krylov subspace methods such as GMRES.

Schwarz method is put forward to solve second order backward stochastic differential equations (2BSDEs) in this work. We will analyze uniqueness, convergence, stability and optimality of the proposed method. Moreover, several simulation results are presented to demonstrate the effectiveness; several applications of the 2BSDEs are investigated. It is concluded from these results that the proposed the method is powerful to calculate the 2BSDEs listing from the financial engineering.

A common way to handle the Tikhonov regularization method for the first kind Fredholm integral equations, is first to discretize and then to work with the final linear system. This unavoidably inflicts discretization errors which may lead to disastrous results, especially when a quadrature rule is used. We propose to regularize directly the integral equation resulting in a continuous Tikhonov problem. The Tikhonov problem is reduced to a simple least squares problem by applying the Golub-Kahan bidiagonalization (GKB) directly to the integral operator. The regularization parameter and the iteration index are determined by the discrepancy principle approach. Moreover, we study the discrete version of the proposed method resulted from numerical evaluating the needed integrals. Focusing on the nodal values of the solution results in a weighted version of GKB-Tikhonov method for linear systems arisen from the Nyström discretization. Finally, we use numerical experiments on a few test problems to illustrate the performance of our algorithms.

This paper concerns the reconstruction of a scalar coefficient of a second-order elliptic equation in divergence form posed on a bounded domain from internal data. This problem finds applications in multi-wave imaging, greedy methods to approximate parameterdependent elliptic problems, and image treatment with partial differential equations. We first show that the inverse problem for smooth coefficients can be rewritten as a linear transport equation. Assuming that the coefficient is known near the boundary, we study the well-posedness of associated transport equation as well as its numerical resolution using discontinuous Galerkin method. We propose a regularized transport equation that allow us to derive rigorous convergence rates of the numerical method in terms of the order of the polynomial approximation as well as the regularization parameter. We finally provide numerical examples for the inversion assuming a lower regularity of the coefficient, and using synthetic data.

Iterative ILU factorizations are constructed, analyzed and applied as preconditioners to solve both linear systems and eigenproblems. The computational kernels of these novel Iterative ILU factorizations are sparse matrix-matrix multiplications, which are easy and efficient to implement on both serial and parallel computer architectures and can take full advantage of existing matrix-matrix multiplication codes. We also introduce level-based and threshold-based algorithms in order to enhance the accuracy of the proposed Iterative ILU factorizations. The results of several numerical experiments illustrate the efficiency of the proposed preconditioners to solve both linear systems and eigenvalue problems.

Shifted symmetric higher-order power method (SS-HOPM) has been proved effective in solving the nonlinear eigenvalue problem oriented from the Bose-Einstein Condensation (BEC-like NEP for short) both theoretically and numerically. However, the convergence of the sequence generated by SS-HOPM is based on the assumption that the real eigenpairs of BEC-like NEP are finite. In this paper, we will establish the point-wise convergence via Lojasiewicz inequality by introducing a new related sequence.

We propose a novel algorithm, based on physics-informed neural networks (PINNs) to efficiently approximate solutions of nonlinear dispersive PDEs such as the KdV-Kawahara, Camassa-Holm and Benjamin-Ono equations. The stability of solutions of these dispersive PDEs is leveraged to prove rigorous bounds on the resulting error. We present several numerical experiments to demonstrate that PINNs can approximate solutions of these dispersive PDEs very accurately.

In this paper, the unconditional error estimates are presented for the time-dependent Navier-Stokes equations by the bilinear-constant scheme. The corresponding optimal error estimates for the velocity and the pressure are derived unconditionally, while the previous works require certain time-step restrictions. The analysis is based on an iterated timediscrete system, with which the error function is split into a temporal error and a spatial error. The τ-independent (τ is the time stepsize) error estimate between the numerical solution and the solution of the time-discrete system is proven by a rigorous analysis, which implies that the numerical solution in L^∞-norm is bounded. Thus optimal error estimates can be obtained in a traditional way. Numerical results are provided to confirm the theoretical analysis.

Flexible GMRES (FGMRES) is a variant of preconditioned GMRES, which changes preconditioners at every Arnoldi step. GMRES often has to be restarted in order to save storage and reduce orthogonalization cost in the Arnoldi process. Like restarted GMRES, FGMRES may also have to be restarted for the same reason. A major disadvantage of restarting is the loss of convergence speed. In this paper, we present a heavy ball flexible GMRES method, aiming to recoup some of the loss in convergence speed in the restarted flexible GMRES while keep the benefit of limiting memory usage and controlling orthogonalization cost. Numerical tests often demonstrate superior performance of the proposed heavy ball FGMRES to the restarted FGMRES.

By employing EQ₁^rot nonconforming finite element, the numerical approximation is presented for multi-term time-fractional mixed sub-diffusion and diffusion-wave equation on anisotropic meshes. Comparing with the multi-term time-fractional sub-diffusion equation or diffusion-wave equation, the mixed case contains a special time-space coupled derivative, which leads to many difficulties in numerical analysis. Firstly, a fully discrete scheme is established by using nonconforming finite element method (FEM) in spatial direction and L1 approximation coupled with Crank-Nicolson (L1-CN) scheme in temporal direction. Furthermore, the fully discrete scheme is proved to be unconditional stable. Besides, convergence and superclose results are derived by using the properties of EQ₁^rot nonconforming finite element. What's more, the global superconvergence is obtained via the interpolation postprocessing technique. Finally, several numerical results are provided to demonstrate the theoretical analysis on anisotropic meshes.

In this paper, we systematically construct two classes of structure-preserving schemes with arbitrary order of accuracy for canonical Hamiltonian systems. The one class is the symplectic scheme, which contains two new families of parameterized symplectic schemes that are derived by basing on the generating function method and the symmetric composition method, respectively. Each member in these schemes is symplectic for any fixed parameter. A more general form of generating functions is introduced, which generalizes the three classical generating functions that are widely used to construct symplectic algorithms. The other class is a novel family of energy and quadratic invariants preserving schemes, which is devised by adjusting the parameter in parameterized symplectic schemes to guarantee energy conservation at each time step. The existence of the solutions of these schemes is verified. Numerical experiments demonstrate the theoretical analysis and conservation of the proposed schemes.

The purpose of this paper is to discuss representations of high order C⁰ finite element spaces on simplicial meshes in any dimension. When computing with high order piecewise polynomials the conditioning of the basis is likely to be important. The main result of this paper is a construction of representations by frames such that the associated L² condition number is bounded independently of the polynomial degree. To our knowledge, such a representation has not been presented earlier. The main tools we will use for the construction is the bubble transform, introduced previously in[1], and properties of Jacobi polynomials on simplexes in higher dimensions. We also include a brief discussion of preconditioned iterative methods for the finite element systems in the setting of representations by frames.

This paper develops and analyzes a fully discrete finite element method for a class of semilinear stochastic partial differential equations (SPDEs) with multiplicative noise. The nonlinearity in the diffusion term of the SPDEs is assumed to be globally Lipschitz and the nonlinearity in the drift term is only assumed to satisfy a one-sided Lipschitz condition. These assumptions are the same ones as the cases where numerical methods for general nonlinear stochastic ordinary differential equations (SODEs) under “minimum assumptions” were studied. As a result, the semilinear SPDEs considered in this paper are a direct generalization of these nonlinear SODEs. There are several difficulties which need to be overcome for this generalization. First, obviously the spatial discretization, which does not appear in the SODE case, adds an extra layer of difficulty. It turns out a spatial discretization must be designed to guarantee certain properties for the numerical scheme and its stiffness matrix. In this paper we use a finite element interpolation technique to discretize the nonlinear drift term. Second, in order to prove the strong convergence of the proposed fully discrete finite element method, stability estimates for higher order moments of the H¹-seminorm of the numerical solution must be established, which are difficult and delicate. A judicious combination of the properties of the drift and diffusion terms and some nontrivial techniques is used in this paper to achieve the goal. Finally, stability estimates for the second and higher order moments of the L²-norm of the numerical solution are also difficult to obtain due to the fact that the mass matrix may not be diagonally dominant. This is done by utilizing the interpolation theory and the higher moment estimates for the H¹-seminorm of the numerical solution. After overcoming these difficulties, it is proved that the proposed fully discrete finite element method is convergent in strong norms with nearly optimal rates of convergence. Numerical experiment results are also presented to validate the theoretical results and to demonstrate the efficiency of the proposed numerical method.