# Frozen Jacobian Multistep Iterative Method for Solving Nonlinear IVPs and BVPs

Authors: Fayyaz Ahmad; Shafiq Ur Rehman; Malik Zaka Ullah; Hani Moaiteq Aljahdali; Shahid Ahmad; Ali Saleh Alshomrani; Juan A. Carrasco; Shamshad Ahmad; Sivanandam Sivasankaran
Journal: Complexity
 (2017)

---

## Abstract

In this paper, we present and illustrate a frozen Jacobian multistep iterative method to solve systems of nonlinear equations associated with initial value problems (IVPs) and boundary value problems (BVPs). We have used Jacobi-Gauss-Lobatto collocation (J-GL-C) methods to discretize the IVPs and BVPs. Frozen Jacobian multistep iterative methods are computationally very efficient. They require only one inversion of the Jacobian in the form of LU-factorization. The LU factors can then be used repeatedly in the multistep part to solve other linear systems. The convergence order of the proposed iterative method is5m-11, where m is the number of steps. The validity, accuracy, and efficiency of our proposed frozen Jacobian multistep iterative method is illustrated by solving fifteen IVPs and BVPs. It has been observed that, in all the test problems, with one exception in this paper, a single application of the proposed method is enough to obtain highly accurate numerical solutions. In addition, we present a comprehensive comparison of J-GL-C methods on a collection of test problems.


---

## Body

## 1. Introduction

Iterative methods for solving systems of nonlinear equations associated with IVPs and BVPs are important because, in general, it is hard to find a closed form solution. Generally, nonlinear IVPs and BVPs are solved in two steps. For the first step, we discretize the nonlinear problem to obtain a system of nonlinear equations. In the second step, we use some iterative method to solve the system of nonlinear equations.Pseudospectral collocation methods offer excellent accuracy, and we use them for the discretization of IVPs and BVPs. These methods convert the targeted partial or ordinary differential equations into a system of algebraic equations. Recently, Doha et al. [1] used a J-GL-C method for the discretization of a nonlinear 1+1 Schrödinger equation to obtain a system of ordinary differential equations and solved it by using an implicit Runge-Kutta method. Highly accurate results for four different kinds of nonlinear 1+1 Schrödinger equations were obtained. In [2], nonlinear reaction-diffusion equations were solved by using a J-GL-C method. In [3], a J-GL-C method was also used to solve nonlinear complex generalized Zakharov systems of equations. Many authors (see, e.g., [4–7]) have used pseudospectral collocation method to solve nonlinear IVPs and BVPs. J-GL-C methods are attractive because they depend on two parameters and by changing the values of these parameters we can obtain different pseudospectral collocation methods, for example, the Legendre, Chebyshev, and Gegenbauer pseudospectral collocation methods. J-GL-C methods are based on Jacobi polynomials. The Jacobi polynomials are the eigenfunctions of the singular Sturm-Liouville problem [8],(1)1-x2ϑ′′x+θ-ϕ+θ+ϕ+2xϑ′x+nn+θ+ϕ+1ϑx=0,and can be generated by using the following recurrence relation: (2)Jk+1θ,ϕx=akθ,ϕ-bkθ,ϕJkθ,ϕx,-ckθ,ϕJk-1θ,ϕx,k≥1,J0θ,ϕx=1,J1θ,ϕx=12θ+ϕ+2x+12ϕ-θ,where (3)akθ,ϕ=2k+θ+ϕ+12k+θ+ϕ+22k+1k+θ+ϕ+1,bkθ,ϕ=θ2-ϕ22k+θ+ϕ+12k+1k+θ+ϕ-12k+θ+ϕ,ckθ,ϕ=k+θk+ϕ2k+θ+ϕ+2k+1k+θ+ϕ-12k+θ+ϕ.The pth derivative of the kth-degree Jacobi polynomial Jk(θ,ϕ) is given by (4)DpJkθ,ϕx=Γj+θ+ϕ+p+12pΓj+θ+ϕ+1Jk-pθ+p,ϕ+px.It should be noticed that Jacobi polynomials are orthogonal over the domain [-1,1] with the weight function (1+x)θ(1-x)ϕ. We are interested in approximating the differential operators by using J-GL-C methods. The easiest way to work with differential operators of different orders is to construct the differentiation matrix; for details on the construction of numerical differentiation matrix on Jacobi polynomials over the domain [-1,1], see [9]. Suppose Q is the Jacobi-Gauss-Lobatto differentiation matrix for the first-order derivative operator over the domain [-1,1]. Then, a derivative of order p can be approximated over the interval [a,b] as (5)dpdxp≈2b-aQp.After having a setup for the discretization of IVPs and BVPs by using J-GL-C methods, we are looking for an efficient iterative method for the solution of the associated system of nonlinear equations. The first iterative method that comes into mind to solve a system of nonlinear equationsF(y)=0 is the Newton-Raphson (NR) method [10, 11], which can be written as(6)NR=y0=initialguessyn+1=yn-F′yn-1Fyn,where det⁡F′(yn)≠0. The multistep frozen Jacobian version of the Newton-Raphson method (MNR) can be written as (7)MNR=Number  of  steps=m≥1Convergence  order=m+1Function  evaluations=mJacobian  evaluations=1Number  of  LU-factorization=1Number  of  lower  and  upper  triangular  systems=mBasemethod⟶Number  of  stepNumber  of  stepNumber  of  stepMulti-step  part⟶Number  of  stepy0=initialguessF′y0ϕ1=Fy0y1=y0-ϕ1forjsfors=1,m-1forsF′y0ϕs+1=Fysforsys+1=ys-ϕs+1endy0=ym.Frozen Jacobian methods have the advantage of computing LU decomposition of the Jacobian, because this makes the solution of linear systems, having the Jacobian as their matrix, inexpensive. However, the convergence order of MNR as a function of the number of function evaluations is relatively low, and many researchers [12–23] have obtained frozen Jacobian multistep iterative methods with higher convergence orders for the same number of function evaluations. For example, HJ, FTUC, and MSF are presented in [22], [17], and [18], respectively. These methods can be describe as(8)HJ=Numberofsteps=m≥2Convergenceorder=2mFunctionevaluations=m-1Jacobianevaluations=2LU-factorization=1Matrix-vectormultiplications=mVector-vectormultiplications=2mNumberofsolutionsofloweranduppertriangularsystemsofequations=2m-1Number  of  stepBasemethod⟶Number  of  stepNumber  of  stepNumber  of  stepNumber  of  stepNumber  of  stepNumber  of  stepMulti-step  part⟶Number  of  stepy0=initialguessF′y0ϕ1=Fy0y1=y0-23ϕ1F′y0ϕ2=F′y1ϕ1F′y0ϕ3=F′y1ϕ2y2=y0-238ϕ1+3ϕ2-98ϕ3fors=3,mforsF′y0ϕ4=Fys+1forsF′y0ϕ5=F′y1ϕ4forsys=ys-1-52ϕ4+32ϕ5endy0=ym,FTUC=Number  of  steps=m≥3Convergence-order=3m-4Function  evaluations=m-1Jacobian  evaluations=2LU-factorization=1Matrix-vector  multiplications=m-1Vector-vector  multiplications=m+1Number  of  solutions  of  lower  and  upper  triangular  systems  of  equations=2m-2Number  of  stepBasemethod⟶Number  of  stepNumber  of  stepNumber  of  stepNumber  of  stepMulti-step  part⟶Number  of  stepy0=initialguessF′y0ϕ1=Fy0y1=y0-ϕ1F′y0ϕ2=Fy1y2=y1-3ϕ2F′y0ϕ3=F′y2ϕ2F′y0ϕ4=F′y2ϕ3y3=y1-74ϕ2+12ϕ3+14ϕ4fors=4,mforsF′y0ϕ5=FysforsF′y0ϕ6=F′y2ϕ5forsys=ys-1-2ϕ5+ϕ6endy0=ym,MSF=Number  of  steps=mConvergence-order=3mFunction  evaluations=mJacobian  evaluations=2Second-order  Fréchet  derivative=1LU-factorization=1Matrix-vector  multiplications=2m-2Vector-vector  multiplications=m+2Number  of  solutions  of  lower  and  upper  triangular  systems  of  equations=3m-1Number  of  stepBasemethod⟶Number  of  stepNumber  of  stepNumber  of  stepNumber  of  stepNumber  of  stepMulti-step  part⟶Number  of  stepy0=initialguessF′y0ϕ1=Fy0F′y0ϕ2=F′′y0ϕ12y1=y0-ϕ1-12ϕ2forjsssssfors=2,mforsF′y0ϕ3=FysforsF′y0ϕ4=F′y1ϕ3forsF′y0ϕ5=F′y1ϕ4forsys=ys-1-3ϕ3+3ϕ4-ϕ5endy0=ym.In this paper, we propose another frozen Jacobian multistep iterative method for solving systems of nonlinear equations with higher convergence order than the previously proposed frozen Jacobian methods for the same number of function evaluations. Another goal of the paper is to compare different pseudospectral collocation methods generated by J-GL-C methods by choosing different values for the parametersθ and ϕ.


## 2. New Multistep Iterative Method

To construct the higher order frozen Jacobian multistep iterative method we first construct method with unknown coefficients as(9)fori=1,3foriiF′y0ϕi=Fyi-1foriiyi=yi-1-αi,1ϕiendF′y0ϕ4=Fy3forj=5,qforiiF′y0ϕj=F′y3ϕj-1endy4=y3+∑k=1qα4,kϕkF′y0ϕq+1=Fy4forj=q+2,rforiiF′y0ϕj=F′y3ϕj-1endy5=y4+∑k=q+1rα5,kϕk.We haver+3 unknowns: namely, α1,1,α2,1,α3,1,α4,1,…,α4,q,α5,q+1,…,α5,r. To develop our method we find q, r, and these unknowns to obtain the maximum convergence order. Once we have solved that problem we write down the multistep part as (10)fors=5,mforiiF′y0ϕq+1=Fys-1foriiforj=q+2,rforhhhhiiF′y0ϕj=F′y3ϕj-1foriiendforiiys=ys-1+∑k=q+1rα5,kϕkend.Notice that α5,q+1,…,α5,r are already computed in the base method part and remain fixed in the multistep part. After finding the values of unknowns, we establish our proposed new multistep iterative method (NMIM) as (11)NMIM=Number  of  steps=m≥3Convergence  order=5m-11Function  evaluations=mJacobian  evaluations=2LU-factorization=1Matrix-vector  multiplications=4m-3Vector-vector  multiplications=5m-3Number  of  lower  and  upper  triangular  systems=5m-12Number  of  stepNumber  of  stepNumber  of  stepBasemethod⟶Number  of  stepNumber  of  stepNumber  of  stepNumber  of  stepNumber  of  stepNumber  of  stepNumber  of  stepNumber  of  stepMulti-step  part⟶Number  of  stepy0=initialguessforj=1,3forjF′y0ϕj=Fyj-1forjyj=yj-1-ϕjendF′y0ϕ4=Fy3forj=5,8forjF′y0ϕj=F′y3ϕj-1endy4=y3-214ϕ4+11ϕ5-232ϕ6+6ϕ7-54ϕ8forjsssssfors=5,mforsF′y0ϕ9=Fys-1forsforj=10,13forsforjF′y0ϕj=F′y3ϕj-1forsendforsys=ys-1-112ϕ9+12ϕ10-13ϕ11+7ϕ12-32ϕ13endy0=ym.In NMIM, the convergence order of the base method is nine, and we get an increment of five in the convergence order per additional function evaluation. The LU factors of F′(y0) are used in the multistep part to solve five lower and upper triangular systems of equations. NMIM requires only two Jacobian evaluations and the LU-factorization of one of the Jacobians.


## 3. Convergence Analysis

We use the method of mathematical induction to prove the convergence order of NMIM. We prove the convergence order form=4 and, then, using induction, we will prove the convergence order for m>4. In our analysis, we expand the system of nonlinear equations around its simple root by using Taylor’s series and use higher order Fréchet derivatives in the expansion. We use the symbolic algebra of Maple software to deal with symbolic computations. The Fréchet differentiability of the system of nonlinear equations is an essential constraint of our proof. Gateaux differentiability would not be enough, but that differentiability is neither enough to ensure that F(x) has a good local linear approximation, which is in soul of Newton-like methods. A function F(·) is said to be Fréchet differentiable at a point y if there exists a linear operator A∈L(Rn,Rq) such that (12)limh→0Fy+h-Fy-Ahh=0.The linear operator A is called the first-order Fréchet derivative and we denote it by F′(y). The higher order Fréchet derivatives can be computed recursively using (13)F′y=JacobianFy,Fjyvj-1=JacobianFj-1yvj-1,j≥2,where v is a vector independent from y.Theorem 1.
LetF:Γ⊆Rn→Rn be a sufficiently Fréchet differentiable function on an open convex neighborhood Γ of y∗∈Rn with F(y∗)=0 and det⁡(F′(y∗))≠0, where F′(y) denotes the Fréchet derivative of F(y). Let A1=F′(y∗) and Aj=1/j!F′(y∗)-1F(j)(y∗), for j≥2, where F(j)(y) denotes the j-order Fréchet derivative of F(y). Then, for m=4, with an initial guess in the neighborhood of y∗, the sequence {ym} generated by NMIM converges to y∗ with at least local order of convergence nine and error equation (14)e4=Le09+Oe010,where e0=y0-y∗, e0p=(e0,e0,…,e0)︷ptimes, and L=184A28+8A23A3A23-24A3A26 is a 9-linear function: that is, L∈LRn,Rn,Rn,…,Rn︷9times with Le09∈Rn.Proof.
We define the error at thenth step as en=yn-y∗. To complete the convergence proof, we perform detailed computations using Maple with the following results: (15)Fy0=A1e0+A2e02+A3e03+A4e04+A5e05+A6e06+A7e07+A8e08+Oe010,F′y0-1=I-2A2e0+-3A3+4A22e02+-4A4+6A3A2+6A2A3-8A23e03+-5A5-12A22A3-12A3A22-12A2A3A2+8A4A2+9A32+8A2A4+16A24e04+-6A6-16A22A4-32A25-16A4A22-18A32A2-16A2A4A2-18A3A2A3-18A2A32+24A3A23+24A2A3A22+24A22A3A2+10A5A2+12A4A3+12A3A4+10A2A5+24A23A3e05+⋯+Oe010A1-1,ϕ1=e0-A2e02+-2A3+2A22e03+-3A4+4A2A3+3A3A2-4A23e04+-4A5+6A2A4+4A4A2-6A2A3A2+6A32+8A24-6A3A22-8A22A3e05+-5A6-12A3A2A3-8A2A4A2+8A2A5+9A3A4+8A4A3+5A5A2+12A22A3A2+12A2A3A22+12A3A23-9A32A2-12A2A32-16A25-12A22A4-8A4A22+16A23A3e06+⋯+Oe010,e1=A2e02+2A3-2A22e03+3A4-4A2A3-3A3A2+4A23e04+4A5-6A2A4-4A4A2+6A2A3A2-6A32-8A24+6A3A22+8A22A3e05+5A6+12A3A2A3+8A2A4A2-8A2A5-9A3A4-8A4A3-5A5A2-12A22A3A2-12A2A3A22-12A3A23+9A32A2+12A2A32+16A25+12A22A4+8A4A22-16A23A3e06+⋯+Oe010,Fy1=A1A2e02+2A3-2A22e03+3A4-4A2A3-3A3A2+5A23e04+4A5-6A2A4-6A32-4A4A2-12A24+8A2A3A2+6A3A22+10A22A3e05+5A6+8A4A22+12A3A2A3-24A23A3-19A22A3A2-19A2A3A22-11A3A23+9A32A2+11A2A4A2+16A2A32+28A25-8A2A5-9A3A4-8A4A3-5A5A2+15A22A4e06+⋯+Oe010,ϕ2=A2e02+-4A22+2A3e03+3A4+13A23-8A2A3-6A3A2e04+4A5-12A2A4-12A32-8A4A2-38A24+20A2A3A2+18A3A22+26A22A3e05+5A6+24A4A22+36A3A2A3-76A23A3-59A22A3A2-55A2A3A22-50A3A23+27A32A2+27A2A4A2+40A2A32+104A25-16A2A5-18A3A4-16A4A3-10A5A2+39A22A4e06+⋯+Oe010,e2=2e03A22+-9A23+4A2A3+3A3A2e04+30A24+6A2A4+6A32+4A4A2-14A2A3A2-12A3A22-18A22A3e05+47A22A3A2+43A2A3A22+38A3A23-16A4A22-18A32A2-19A2A4A2-24A3A2A3-28A2A32+60A23A3-88A25+8A2A5+9A3A4+8A4A3+5A5A2-27A22A4e06+⋯+Oe010,Fy2=A12A22e03+-9A23+4A2A3+3A3A2e04+30A24+6A2A4+6A32+4A4A2-14A2A3A2-12A3A22-18A22A3e05+47A22A3A2+43A2A3A22+38A3A23-16A4A22-18A32A2-19A2A4A2-24A3A2A3-28A2A32+60A23A3-84A25+8A2A5+9A3A4+8A4A3+5A5A2-27A22A4e06+⋯+Oe010,ϕ3=2A22e03+-13A23+4A2A3+3A3A2e04+6A2A4+4A4A2-20A2A3A2+56A24+6A32-18A3A22-26A22A3e05+-27A2A4A2-36A3A2A3+8A2A5+9A3A4+8A4A3+5A5A2+87A22A3A2+79A2A3A22+77A3A23-24A4A22-27A32A2-40A2A32+112A23A3-196A25-39A22A4e06+⋯+Oe010,e3=4A23e04+-26A24+6A2A3A2+6A3A22+8A22A3e05+-40A22A3A2-36A2A3A22-39A3A23+8A4A22+9A32A2+8A2A4A2+12A3A2A3+12A2A32-52A23A3+108A25+12A22A4e06+-80A22A32-72A2A3A2A3-78A3A22A3+16A4A2A3+18A33+16A2A4A3+18A2A3A4+18A3A2A4+216A24A3-78A23A4-356A26+10A5A22+12A4A3A2+12A3A4A2+10A2A5A2+168A23A3A2+150A22A3A22+148A2A3A23+168A3A24+16A22A5-52A4A23-54A32A22-48A2A4A22-54A2A32A2-60A3A2A3A2-54A22A4A2e07+⋯+Oe010,Fy3=A14A23e04+-26A24+6A2A3A2+6A3A22+8A22A3e05+-40A22A3A2-36A2A3A22-39A3A23+8A4A22+9A32A2+8A2A4A2+12A3A2A3+12A2A32-52A23A3+108A25+12A22A4e06+-80A22A32-72A2A3A2A3-78A3A22A3+16A4A2A3+18A33+16A2A4A3+18A2A3A4+18A3A2A4+216A24A3-78A23A4-356A26+10A5A22+12A4A3A2+12A3A4A2+10A2A5A2+168A23A3A2+150A22A3A22+148A2A3A23+168A3A24+16A22A5-52A4A23-54A32A22-48A2A4A22-54A2A32A2-60A3A2A3A2-54A22A4A2e07+⋯+Oe010,F′y3=A1I+8A24e04+12A2A3A22+12A22A3A2-52A25+16A23A3e05+24A22A32+24A2A3A2A3-80A23A3A2-72A22A3A22-78A2A3A23+16A2A4A22+18A2A32A2+16A22A4A2+24A23A4+216A26-104A24A3e06+36A22A3A4+36A2A3A2A4+296A22A3A23+300A23A3A22+20A2A5A22+24A2A4A3A2+24A2A3A4A2+20A22A5A2+336A24A3A2-156A2A3A22A3-144A22A3A2A3-160A23A32-156A24A4+336A2A3A24+36A2A33+32A2A4A2A3-108A23A4A2+32A22A4A3-108A22A32A2-120A2A3A2A3A2-104A2A4A23-108A2A32A22-96A22A4A22+32A23A5+432A25A3-712A27e07+⋯+Oe010,ϕ4=4A23e04+-34A24+6A2A3A2+6A3A22+8A22A3e05+12A2A32-68A23A3+176A25+12A22A4+8A2A4A2-52A22A3A2-48A2A3A22-51A3A23+8A4A22+9A32A2+12A3A2A3e06+-104A22A32-96A2A3A2A3-102A3A22A3+16A4A2A3+18A33+16A2A4A3+18A2A3A4+18A3A2A4+352A24A3-102A23A4-708A26+10A5A22+12A4A3A2+12A3A4A2+10A2A5A2+272A23A3A2+246A22A3A22+250A2A3A23+270A3A24+16A22A5-68A4A23-72A32A22-64A2A4A22-72A2A32A2-78A3A2A3A2-70A22A4A2e07+⋯+Oe010,ϕ5=4A23e04+-42A24+6A2A3A2+6A3A22+8A22A3e05+8A2A4A2+12A3A2A3-64A22A3A2-60A2A3A22-63A3A23+8A4A22+9A32A2+12A2A32-84A23A3+260A25+12A22A4e06+-128A22A32-120A2A3A2A3-126A3A22A3+16A4A2A3+18A33+16A2A4A3+18A2A3A4+18A3A2A4+520A24A3-126A23A4-1228A26+10A5A22+12A4A3A2+12A3A4A2+10A2A5A2+400A23A3A2+366A22A3A22+376A2A3A23+396A3A24+16A22A5-84A4A23-90A32A22-80A2A4A22-90A2A32A2-96A3A2A3A2-86A22A4A2e07+⋯+Oe010ϕ6=4A23e04+-50A24+6A2A3A2+6A3A22+8A22A3e05+-76A22A3A2-72A2A3A22-75A3A23+8A4A22+9A32A2+8A2A4A2+12A3A2A3+12A2A32-100A23A3+360A25+12A22A4e06+-152A22A32-144A2A3A2A3-150A3A22A3+16A4A2A3+18A33+16A2A4A3+18A2A3A4+18A3A2A4+720A24A3-150A23A4-1948A26+10A5A22+12A4A3A2+12A3A4A2+10A2A5A2+552A23A3A2+510A22A3A22+526A2A3A23+546A3A24+16A22A5-100A4A23-108A32A22-96A2A4A22-108A2A32A2-114A3A2A3A2-102A22A4A2e07+⋯+Oe010,ϕ7=4A23e04+-58A24+6A2A3A2+6A3A22+8A22A3e05+12A3A2A3+8A2A4A2+476A25+12A22A4+12A2A32-116A23A3-88A22A3A2-84A2A3A22-87A3A23+8A4A22+9A32A2e06+-176A22A32-168A2A3A2A3-174A3A22A3+16A4A2A3+18A33+16A2A4A3+18A2A3A4+18A3A2A4+952A24A3-174A23A4-2900A26+10A5A22+12A4A3A2+12A3A4A2+10A2A5A2+728A23A3A2+678A22A3A22+700A2A3A23+720A3A24+16A22A5-116A4A23-126A32A22-112A2A4A22-126A2A32A2-132A3A2A3A2-118A22A4A2e07+⋯+Oe010,ϕ8=4A23e04+-66A24+6A2A3A2+6A3A22+8A22A3e05+12A3A2A3+12A2A32-132A23A3+608A25+12A22A4-100A22A3A2-96A2A3A22-99A3A23+8A4A22+9A32A2+8A2A4A2e06+-200A22A32-192A2A3A2A3-198A3A22A3+16A4A2A3+18A33+16A2A4A3+18A2A3A4+18A3A2A4+1216A24A3-198A23A4-4116A26+10A5A22+12A4A3A2+12A3A4A2+10A2A5A2+928A23A3A2+870A22A3A22+898A2A3A23+918A3A24+16A22A5-132A4A23-144A32A22-128A2A4A22-144A2A32A2-150A3A2A3A2-134A22A4A2e07+⋯+Oe010e4=184A28+8A23A3A23-24A3A26e09+Oe010.Now, we present the proof of convergence of NMIM via induction.Theorem 2.
The convergence order of NMIM form≥4 is 5m-11.Proof.
All the computations are made under the assumption of Theorem1. We know from Theorem 1 that the convergence order of NMIM method is nine for m=4. By performing the computations in a similar manner to that of Theorem 1, we get the following error equation for m=5: (16)e5=4A23A3-12A3A23+60A25184A28+8A23A3A23-24A3A26e014+Oe015.This proves that the convergence order of NMIM is 5m-11 for m=5. Now, assume that the convergence order of NMIM is 5m-11 for a number of steps’ m≥5 and and let us prove that the convergence order of NMIM for m steps is 5m-6 for (m)th step. If the convergence order of NMIM for m steps is 5m-11, then (17)em=ym-y∗~de05m-11,where d is the asymptotic constant and symbol ~ means “approximately equal to.” By using (17), we perform the following steps:(18)Fy0-1~I-2A2e0A1-1,F′y3~A1I+8A24e04,Fys~A1de05m-11,ϕ9~d-2A2de0e05m-11,ϕ10~d-4A2de0+4A22de02+8A24de04-32A25de05e05m-11,ϕ11~d-6A2de0+12A22de02-8A23de03+16A24de04-96A25de05e05m-11,ϕ12~d-8A2de0+24A22de02-32A23de03+40A24de04-192A25de05e05m-11,ϕ13~d-10A2de0+40A22de02-80A23de03+112A24de04-352A25de05e05m-11,em+1~48A25de05m-6.This completes the proof that the convergence order of NMIM is 5m-6 for m+1 steps.


## 4. Computational Cost

The computational costs in terms of multiplications of different iterative methods are shown in Table1, where m is the number of steps, n is the number of dimensions of the problem, l is the ratio between the computational cost of a division and the computational cost of a multiplication, α is the ratio between the computational cost of a function component evaluation and the computational cost of a multiplication, β is the ratio between the computational cost of a component of the Jacobian and the computational cost of a multiplication, and γ is the ratio between the computational cost of a second derivative and the computational cost of a multiplication. Table 2 displays the difference of computational costs of methods with respect to NMIM when the methods have the same convergence order as NMIM with m steps. Table 3 gives the conditions on α and bounds on m so that, for the same convergence order, the computational cost of NMIM is smaller than that of the other methods.Table 1
Computational costs of different iterative methods.
Methods
Computational cost
MNR
1
3
n
3
+
β
+
m
-
1
2
+
l
2
n
2
+
m
α
+
l
m
-
m
+
1
6
-
l
2
n
HJ
1
3
n
3
+
2
β
-
3
2
+
l
2
+
3
m
n
2
+
m
α
-
α
+
7
6
-
3
2
l
+
2
l
m
n
FTUC
1
3
n
3
+
2
β
-
7
2
+
l
2
+
3
m
n
2
+
m
α
-
α
+
19
6
-
5
2
l
-
m
+
2
l
m
n
MSF
1
3
n
3
+
2
β
-
7
2
+
l
2
+
5
m
n
2
+
m
α
+
γ
+
19
6
-
3
2
l
-
2
m
+
3
l
m
n
NMIM
1
3
n
3
+
2
β
-
49
2
+
l
2
+
9
m
n
2
+
α
m
-
17
6
-
25
2
l
+
5
l
m
nTable 2
Difference of computational costs of different iterative methods with respect to NMIM when the methods have the same convergence order as of NMIM withm steps.
Methods
Difference of computational costs
MNR-NMIM
4
n
α
m
-
β
n
2
-
4
n
2
m
-
10
n
α
+
2
l
n
-
5
n
m
+
14
n
2
+
13
n
HJ-NMIM
13
2
n
2
-
3
2
n
2
m
+
3
2
n
α
m
-
13
2
n
α
+
4
n
FTUC-NMIM
14
n
2
-
4
n
2
m
+
2
3
n
α
m
-
10
3
n
α
+
25
n
3
+
16
3
l
n
-
5
3
n
m
-
5
3
l
m
n
MSF-NMIM
-
2
3
n
2
m
+
8
3
n
2
+
2
3
n
α
m
-
11
3
n
α
+
γ
n
+
40
n
3
-
10
3
n
mTable 3
Conditions and bounds onm so that, for the same convergence order, the computational cost of NMIM is smaller.
Computational cost
Bounds onm
Conditions
C
M
N
R
>
C
N
M
I
M
m
<
n
β
+
10
α
-
2
l
-
14
n
-
13
-
4
n
+
4
α
-
5
α
<
n
+
5
4
C
H
J
>
C
N
M
I
M
m
<
1
3
-
13
n
+
13
α
-
8
-
n
+
α
α
<
n
C
F
T
U
C
>
C
N
M
I
M
m
<
-
42
n
-
16
l
+
10
α
-
25
-
12
n
-
5
l
+
2
α
-
5
α
<
6
n
+
5
2
+
5
2
l
C
M
S
F
>
C
N
M
I
M
m
<
-
1
2
8
n
-
11
α
+
3
γ
+
40
-
n
+
α
-
5
α
<
n
+
5


## 5. Numerical Simulations

In this section, we verify the convergence order of NMIM and solve fifteen nonlinear initial and boundary value problems to show the validity, accuracy, and efficiency of our proposed iterative method. For the discretization, we use four pseudospectral collocation methods. In all experiments but one, we apply NMIM once with several steps to reduce the absolute error in the computed numerical solution. To verify the convergence order computationally, we adopt the following definition of computational convergence order (CCO):(19)CCO=log⁡Fxk+1∞/Fxk∞log⁡Fxk∞/Fxk-1∞.


### 5.1. Computational Verification of Convergence Order

In general, it is hard to verify the convergence order by solving initial and boundary value problems. Therefore, we consider the following five systems of nonlinear equationsF(x)=0. (20)Problem1=F1x=x2x3+x4x2+x3,F2x=x1x3+x4x1+x3,F3x=x1x2+x4x1+x2,F4x=x1x2+x3x1+x2-1,Problem2=Fix=xi2xi+1-1,i=1,2,…,n-1,Fix=xnx1-1,i=n,Problem3=F1x=3-0.5x1x1-2x2+1,F2x=3-0.5xnxn-2xn-1+1,Fix=3-0.5xixi-xi-1+2xi+1+1,i=2,3,…,n-1,Problem4=Fix=2n+i1-cos⁡xi-sin⁡xi-∑j=1ncos⁡xj2sin⁡xi-cos⁡xi,i=1,2,…,n,Problem5=F1x=10x1+sin⁡x1+x2-1,F2x=8x2-cos2x3-x2-1,F3x=12x3+sin⁡x3-1.The system of nonlinear equations is solved by Mathematica using 11,000 digits. Table4 shows that CCOs agree with theoretical convergence orders.Table 4
Verification of convergence order.
F
(
x
k
)
∞
Iter.
m
=
4
m
=
5
m
=
6
m
=
7
Problem 1 (n=4)
1
6.56
e
-
8
1.09
e
-
11
1.82
e
-
15
3.04
e
-
19
2
6.31
e
-
68
5.74
e
-
158
8.66
e
-
286
2.17
e
-
451
3
4.44
e
-
608
6.90
e
-
2206
6.32
e
-
5422
6.67
e
-
10823
CCO
9
14
19
24
Problem 2 (n=10)
1
7.31
e
-
3
6.45
e
-
4
5.83
e
-
5
5.29
e
-
6
2
1.54
e
-
21
1.38
e
-
47
5.42
e
-
84
8.20
e
-
131
3
1.31
e
-
189
5.63
e
-
659
1.35
e
-
1585
3.06
e
-
3126
CCO
9
14
19
24
Problem 3 (n=100)
1
7.09
e
-
5
2.12
e
-
7
6.36
e
-
10
1.90
e
-
12
2
5.45
e
-
43
1.69
e
-
94
1.53
e
-
166
1.64
e
-
258
3
6.87
e
-
348
8.80
e
-
1143
7.95
e
-
2673
8.85
e
-
5184
CCO
8
12
16
20
Problem 4 (n=5)
2
1.24
e
-
3
8.47
e
-
6
8.17
e
-
9
9.64
e
-
13
3
2.73
e
-
19
2.08
e
-
54
1.45
e
-
120
2.27
e
-
229
4
6.73
e
-
148
1.79
e
-
645
7.35
e
-
1920
5.31
e
-
4575
CCO
8.21
12.2
16.1
20.1
Problem 5 (n=3)
1
3.77
e
-
6
8.57
e
-
10
1.94
e
-
13
4.39
e
-
17
2
2.90
e
-
67
1.71
e
-
151
6.21
e
-
269
1.46
e
-
419
3
3.09
e
-
581
2.79
e
-
1984
1.71
e
-
4737
9.99
e
-
9298
CCO
8.41
12.9
17.5
22.1
Theoretical convergence order (5m-11)
9
14
19
24


### 5.2. Solving Initial and Boundary Value Problems

In this section, we solve different initial and boundary value problems by using different J-GL-C discretization methods. In all computations, we make spatial and temporal discretizations to translate the entire IVPs and BVPs into systems of nonlinear equations. The systems of nonlinear equations are then solved by using the proposed iterative method NMIM. In most of the problems, we use vector0 with initial and boundary conditions introduced in it as the initial guess. The inclusion of initial and boundary conditions in 0 makes the initial guess nonsmooth but it works for most of the problems. When the nonsmooth initial guess does not work, we then make it smooth to get convergence.


#### 5.2.1. The Lane-Emden Equation

The Lane-Emden equation is(21)x′′t+2tx′t+xtp=0,x′0=0,x0=1.J-GL-C discretization methods translate the Lane-Emden equation into the following system of nonlinear equations:(22)Fx=Sxxx+diag2tSxx+xp=0,F′x=Sxx+diag2tSx+pdiagxp-1,where diag(·) denotes a diagonal matrix with values in the diagonal equal to the elements of its vector argument, Sx=2/bQx, Sxx=2/bQx2, and [0,b] is the domain of the problem. We have solved the Lane-Emden equation for different values of p ranging from two to five and b=3. Figure 1 shows the initial guess and the computed numerical solutions for θ=1/2, ϕ=-1/2, and 50 grid points. Table 5 gives the norm of the residue in the solution of the nonlinear systems by a single iteration of NMIM with varying m.Table 5
Norm of the residue of the system of nonlinear equation resulting from Lane-Emden equation with varyingm.
F
(
x
k
)
∞
m
p
2
3
4
5
10
2.69
e
-
17
3.70
e
-
14
1.04
e
-
12
1.03
e
-
11
20
1.36
e
-
46
7.06
e
-
41
9.33
e
-
39
1.22
e
-
37
30
9.33
e
-
76
1.68
e
-
67
7.56
e
-
65
1.97
-
63
40
6.16
e
-
105
4.53
e
-
94
6.19
e
-
91
2.26
e
-
89
50
4.08
e
-
134
1.21
e
-
120
5.11
e
-
117
2.18
e
-
115
60
2.71
e
-
163
3.22
e
-
147
4.24
e
-
143
2.00
e
-
141
70
1.79
e
-
192
8.555
e
-
174
3.52
e
-
169
1.83
e
-
167
80
1.19
e
-
221
2.27
e
-
200
2.92
e
-
195
1.69
e
-
193
90
7.89
e
-
251
6.03
e
-
227
2.43
e
-
221
1.58
e
-
219
100
5.23
e
-
280
1.60
e
-
253
2.02
e
-
247
1.47
e
-
245Figure 1
Initial guess and solutions of the Lane-Emden equation for different values ofp.


#### 5.2.2. Bratu Problem

Bratu problem is the boundary value problem:(23)x′′t+αext=0,x0=0,x1=0,where α is a parameter. The discretization is carried out using J-GL-C methods. The resulting system of nonlinear equations is (24)Fx=Sxxx+αex=0,F′x=Sxx+αdiagex,where Sxx=2Qx2. We considered several values of α, values θ=-1/2 and ϕ=-1/2 for the parameters of J-GL-C, 50 grid points, and solved the resulting system of nonlinear equations using a single application of NMIM with varying m. The initial guess and computed numerical solutions are given in Figure 2. Table 6 gives the norm of the residue of the solution of the system of nonlinear equations for varying m.Table 6
Norm of the residue of the solution of the system of nonlinear equations for Bratu problem and varyingm.
F
(
x
k
)
∞
m
α
1
2
3
10
7.68
e
-
41
4.10
e
-
31
1.97
e
-
28
20
2.20
e
-
93
1.01
e
-
70
9.65
e
-
64
30
7.02
e
-
146
2.50
e
-
110
4.72
e
-
99
40
2.23
e
-
198
6.17
e
-
150
2.30
e
-
134
50
7.13
e
-
251
1.52
e
-
189
1.12
e
-
169
60
2.27
e
-
303
3.75
e
-
229
5.52
e
-
205
70
7.30
e
-
356
9.27
e
-
269
2.69
e
-
240
80
2.30
e
-
408
2.28
e
-
308
1.32
e
-
275
90
7.40
e
-
461
5.70
e
-
348
6.50
e
-
311
100
2.30
e
-
513
1.40
e
-
387
3.20
e
-
346Figure 2
Initial guess and solutions of Bratu problem for different values ofα.


#### 5.2.3. Frank-Kamenetskii Problem

Frank-Kamenetskii problem is(25)x′′t+1tx′t+αext=0,x′0=0,x1=0,where α is a parameter. After discretization by J-GL-C methods we get the system of nonlinear equations:(26)Fx=Sxxx+diag1tSxx+αex=0,F′x=Sxx+diag1tSx+αdiagex.We considered several values of α, values θ=-1/2 and ϕ=-1/2 for the parameters of J-GL-C, 50 grid points, and solved the resulting system of nonlinear equations by using a single application of NMIM with varying m. Figure 3 plots the initial guess and computed solutions and Table 7 gives the norm of the residue of the solution of the system of nonlinear equations.Table 7
Norm of the residue of the solution of the system of nonlinear equations for Fran-Kamenetskii problem and varyingm.
F
(
x
k
)
∞
m
α
1
1.1
1.2
1.3
10
7.09
e
-
41
1.89
e
-
39
2.63
e
-
38
1.80
e
-
37
20
2.46
e
-
82
2.36
e
-
79
6.21
e
-
77
4.05
e
-
75
30
8.58
e
-
124
2.95
e
-
119
1.46
e
-
115
9.09
e
-
113
40
2.98
e
-
165
3.69
e
-
159
3.45
e
-
154
2.03
e
-
150
50
1.03
e
-
206
4.61
e
-
199
8.15
e
-
193
4.57
e
-
188
60
3.61
e
-
248
5.76
e
-
239
1.92
e
-
231
1.02
e
-
225
70
1.25
e
-
289
7.19
e
-
279
4.53
e
-
270
2.29
e
-
263
80
4.40
e
-
331
9.00
e
-
319
1.10
e
-
308
5.15
e
-
301
90
1.50
e
-
372
1.10
e
-
358
2.50
e
-
347
1.20
e
-
338
100
5.30
e
-
414
1.40
e
-
398
6.00
e
-
386
2.60
e
-
376Figure 3
Initial guess and solutions of Frank-Kamenetskii problem with different values ofα.


#### 5.2.4. 1 + 1 Nonlinear Schrödinger Equation

1 + 1 nonlinear Schrödinger equation is(27)i∂tψx,t+∂xxψx,t+2γψx,t2ψx,t-2δRx,tψx,t=0,x,t∈Dx×Dt,where(28)Dx=ax,bx,Dt=0,tf,with the initial and boundary conditions (29)ψax,t=η1t,ψbx,t=η2t,ψx,0=ζx.The function ψ(x,t) is a complex function and can be written as ψ(x,t)=ψ1(x,t)+iψ2(x,t). The Schrödinger equation can be rewritten in terms of the real functions ψ1(x,t) and ψ2(x,t) as(30)∂t∂xx∂xx-∂tψ1ψ2-2δRx,t0110ψ1ψ2+2γψ12+ψ220110ψ1ψ2=0,with initial and boundary conditions (31)ψ1ax,t+iψ2ax,t=η11t+iη12,ψ1bx,t+iψ2bx,t=η21t+iη22t,ψ1x,0+iψ2x,0=ζ1x+iζ2x.Using J-GL-C methods, we can discretize (30) as(32)Fψ=St⊗IxIt⊗SxxIt⊗Sxx-St⊗Ixψ1ψ2-2δOdiagRdiagROψ1ψ2+2γOdiagψ12+ψ22diagψ12+ψ22Oψ1ψ2=0,F′ψ=St⊗IxIt⊗SxxIt⊗Sxx-St⊗Ix-2δOdiagRdiagRO+2γ2diagψ1⊙ψ23diagψ223diagψ122diagψ1⊙ψ2,where ψj=[ψj1,1,ψj1,2,…,ψj1,nt,…,ψjnx,1,ψjnx,2, …,ψj(nx,nt)]T, j=1,2, R=[R1,1,R1,2,…,R1,nt, …,Rnx,1,Rnx,2,…,R(nx,nt)]T, ⊗ is the Kronecker product, ⊙ is the element-wise multiplication, O is a matrix of all zeros, I is the identity matrix, St=2/tfQt, Sxx=2/bx-axQx2, and Q is the J-GL-C operational matrix. The system of nonlinear equations (32) can be rewritten as (33)Fψ=Aψ+gψ-p,where A=St⊗IxIt⊗SxxIt⊗Sxx-St⊗Ix-2δOdiag(R)diag(R)O, g(ψ)=2γ2diag(ψ1⊙ψ2)3diag(ψ22)3diag(ψ12)2diag(ψ1⊙ψ2), and p is the vector incorporating the initial-boundary conditions. Four nonlinear Schrödinger equations are listed in Table 8 with their corresponding analytical solutions. We chose a domain (x,t)∈[-1,1]×[0,1] and solved the equations using different pair of values for the parameters θ and ϕ of J-GL-C and different grids for the domain. As initial guess, we used ψ1=0, ψ2=0 with boundary conditions set for problem 3. For problems 1, 2, and 4 NMIM does not converge with that initial guess and we made the initial guess smoother. We smoothed the initial guess by applying the iteration(34)ψ=-A-12gψ-ponce for problem 1 and twice for problems 2 and 4. Tables 9, 10, 11, and 12 give the norm of the error in the solution when the nonlinear system of equations is solved using a single application of NMIM with varying m. Figures 4–7 display the solutions and the error in the components of ψ1 and ψ2 for θ=-4/2, ϕ=-1/2, and a grid 15×30.Table 8
Nonlinear Schrödinger equations.
Schrödinger equation
Analytical Solution
i∂tψ+∂xxψ+2|ψ|2ψ=0
ψ
(
x
,
t
)
=
e
i
(
x
+
t
)
i
∂
t
ψ
+
∂
x
x
ψ
-
2
|
ψ
|
2
ψ
=
0
ψ
(
x
,
t
)
=
e
i
(
x
-
3
t
)
i∂tψ+12∂xxψ-|ψ|2ψ-cos⁡(x)2ψ=0
ψ
(
x
,
t
)
=
e
-
3
i
t
/
2
sin
⁡
(
x
)
i
∂
t
ψ
+
∂
x
x
ψ
+
2
|
ψ
|
2
ψ
-
2
δ
ψ
=
0
ψ
(
x
,
t
)
=
e
i
(
x
+
(
1
-
2
δ
)
t
)Table 9
Norm of the error in the solution for equation (17).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
5.49
e
-
01
6.21
e
-
01
6.41
e
-
01
4
5.36
e
-
02
2.35
e
-
01
2.44
e
-
01
8
2.85
e
-
04
4.55
e
-
04
4.97
e
-
04
12
2.85
e
-
04
5.81
e
-
08
6.38
e
-
08
16
2.85
e
-
04
5.69
e
-
10
4.38
e
-
12
20
2.85
e
-
04
6.13
e
-
11
5.95
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
6.86
e
-
01
Divergent
6.40
e
-
01
4
1.96
e
-
01
2.48
e
-
01
8
1.64
e
-
03
4.94
e
-
04
12
1.59
e
-
03
6.16
e
-
08
16
1.59
e
-
03
6.02
e
-
12
20
1.59
e
-
03
2.16
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
5.24
e
-
01
6.25
e
-
01
6.41
e
-
01
4
5.90
e
-
02
2.39
e
-
01
2.46
e
-
01
8
1.14
e
-
03
5.05
e
-
04
5.40
e
-
04
12
1.14
e
-
03
7.73
e
-
08
8.00
e
-
08
16
1.14
e
-
03
1.90
e
-
09
5.28
e
-
12
20
1.14
e
-
03
1.90
e
-
09
4.99
e
-
13
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
5.95
e
-
01
6.20
e
-
01
6.41
e
-
01
4
3.73
e
-
02
2.37
e
-
01
2.45
e
-
01
8
1.52
e
-
03
4.24
e
-
04
4.81
e
-
04
12
1.52
e
-
03
5.22
e
-
08
5.81
e
-
08
16
1.52
e
-
03
4.90
e
-
10
1.71
e
-
12
20
1.52
e
-
03
4.90
e
-
10
9.31
e
-
14Table 10
Norm of the error in the solution for equation (30).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
1.28
e
+
00
1.34
e
+
00
1.34
e
+
00
4
4.27
e
-
01
9.65
e
-
01
9.72
e
-
01
8
1.25
e
-
03
2.15
e
-
03
2.18
e
-
03
12
1.46
e
-
04
9.17
e
-
07
1.05
e
-
06
16
1.46
e
-
04
3.32
e
-
10
4.61
e
-
10
20
1.46
e
-
04
2.90
e
-
11
1.85
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
1.27
e
+
00
1.35
e
+
00
1.34
e
+
00
4
4.34
e
-
01
9.46
e
-
01
9.72
e
-
01
8
1.63
e
-
03
2.18
e
-
03
2.18
e
-
03
12
1.52
e
-
03
3.89
e
-
07
8.77
e
-
07
16
1.52
e
-
03
5.28
e
-
09
1.95
e
-
10
20
1.52
e
-
03
3.73
e
-
09
1.34
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
1.26
e
+
00
1.34
e
+
00
1.34
e
+
00
4
3.98
e
-
01
9.24
e
-
01
9.40
e
-
01
8
2.63
e
-
03
1.35
e
-
03
1.63
e
-
03
12
1.00
e
-
03
2.63
e
-
06
2.85
e
-
06
16
1.00
e
-
03
2.35
e
-
08
1.65
e
-
08
20
1.00
e
-
03
7.88
e
-
10
6.29
e
-
11
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
1.33
e
+
00
1.34
e
+
00
1.34
e
+
00
4
3.90
e
-
01
9.80
e
-
01
9.80
e
-
01
8
1.03
e
-
03
2.15
e
-
03
2.19
e
-
03
12
1.03
e
-
03
6.56
e
-
07
7.42
e
-
07
16
1.03
e
-
03
3.08
e
-
10
1.57
e
-
10
20
1.03
e
-
03
2.95
e
-
10
5.77
e
-
14Table 11
Norm of the error in the solution for equation (32).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomial(θ,ϕ)=(0,0)
9.98
e
-
02
8.73
e
-
02
9.66
e
-
02
4
7.62
e
-
04
3.16
e
-
04
3.25
e
-
04
8
4.27
e
-
05
3.64
e
-
11
3.12
e
-
11
12
4.27
e
-
05
3.16
e
-
11
1.89
e
-
13
16
4.27
e
-
05
3.16
e
-
11
1.88
e
-
13
20
4.27
e
-
05
3.16
e
-
11
1.89
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
9.71
e
-
02
9.33
e
-
02
9.67
e
-
02
4
1.18
e
-
03
3.37
e
-
04
3.36
e
-
04
8
7.68
e
-
04
3.23
e
-
08
2.93
e
-
11
12
7.68
e
-
04
1.89
e
-
10
1.07
e
-
13
16
7.68
e
-
04
1.87
e
-
10
1.06
e
-
13
20
7.68
e
-
04
1.87
e
-
10
1.07
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
9.47
e
-
02
9.17
e
-
02
9.65
e
-
02
4
9.17
e
-
04
3.14
e
-
04
3.38
e
-
04
8
6.70
e
-
04
9.13
e
-
11
8.39
e
-
11
12
6.70
e
-
04
8.92
e
-
11
1.19
e
-
13
16
6.70
e
-
04
8.92
e
-
11
1.19
e
-
13
20
6.70
e
-
04
8.92
e
-
11
1.19
e
-
13
1
Chebyshev polynomial of first kind(θ,ϕ)=(-1/2,-1/2)
1.01
e
-
01
9.01
e
-
02
9.61
e
-
02
4
1.31
e
-
03
2.99
e
-
04
3.26
e
-
04
8
8.53
e
-
04
4.90
e
-
11
2.27
e
-
11
12
8.53
e
-
04
5.59
e
-
11
3.28
e
-
14
16
8.53
e
-
04
5.59
e
-
11
3.24
e
-
14
20
8.53
e
-
04
5.59
e
-
11
3.20
e
-
14Table 12
Norm of the error in the solution for equation (35).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
6.44
e
-
01
7.84
e
-
01
8.01
e
-
01
4
6.52
e
-
02
4.18
e
-
02
4.18
e
-
02
8
2.77
e
-
04
1.36
e
-
06
1.33
e
-
06
12
2.75
e
-
04
1.70
e
-
10
1.81
e
-
10
16
2.75
e
-
04
4.14
e
-
11
4.87
e
-
13
20
2.75
e
-
04
4.14
e
-
11
4.67
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
6.51
e
-
01
Divergent
8.08
e
-
01
4
4.52
e
-
02
4.15
e
-
02
8
1.50
e
-
03
1.25
e
-
06
12
1.51
e
-
03
1.17
e
-
10
16
1.51
e
-
03
1.94
e
-
13
20
1.51
e
-
03
1.99
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
7.48
e
-
01
7.85
e
-
01
8.04
e
-
01
4
5.30
e
-
02
4.24
e
-
02
4.22
e
-
02
8
1.02
e
-
03
2.20
e
-
06
2.61
e
-
06
12
1.02
e
-
03
2.67
e
-
09
1.00
e
-
09
16
1.02
e
-
03
2.25
e
-
09
4.86
e
-
13
20
1.02
e
-
03
2.25
e
-
09
1.96
e
-
13
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
7.03
e
-
01
7.88
e
-
01
8.02
e
-
01
4
3.66
e
-
02
4.20
e
-
02
4.17
e
-
02
8
1.46
e
-
03
1.27
e
-
06
1.18
e
-
06
12
1.46
e
-
03
5.96
e
-
10
2.39
e
-
11
16
1.46
e
-
03
5.84
e
-
10
8.56
e
-
14
20
1.46
e
-
03
5.84
e
-
10
8.49
e
-
14Figure 4
Solution and error components for (17).Figure 5
Solution and error components for (30).Figure 6
Solution and error components for (32).Figure 7
Solution and error components for (35).


#### 5.2.5. Klein Gordon Equation

Klein Gordon equation is(35)uttx,t-c2uxxx,t+ku-γux,t3=0,x,t∈Dx×Dt,where(36)Dx=ax,bx,Dt=0,tf,with the initial and boundary conditions (37)uax,t=ϱ1t,ubx,t=ϱ2t,ux,0=ξ1x,utx,0=ξ2x.By discretization using J-GL-C methods we obtain (38)Fu=Sttu-c2Sxxu+ku-γu3=0,F′u=Stt-c2Sxx+kI-3γdiagu2,where u=[u(1,1),u(1,2),…,u(1,nt),…,u(nx,1),u(nx,2),…, u(nx,nt)]T, Stt=2/tfQt2, and Sxx=2/bx-axQx2. The analytical solution of (35) is u(x,t)=δsech⁡κ(x-vt), where κ=k/c2-v2 and δ=2k/γ. We set the parameters c=1, γ=1, v=0.5, and k=0.5, picked up the domain (x,t)∈[-10,10]×[0,1], and solved the nonlinear system of equations by a single application of NMIM with initial guess u=0, with the boundary conditions introduced. We used several pairs of the θ, ϕ parameters of J-GL-C. According to Table 13, the Chebyshev collocation method of first kind exhibits best accuracy when the grid is more refined. The solution and the error components for θ=-1/2, ϕ=-1/2, nx=120, and nt=20 are plotted in Figure 8.Table 13
Norm of the error in the Klein Gordon equation as a function ofm for several values for θ and ϕ and grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
30
×
30
60
×
30
120
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
3.80
e
-
01
4.10
e
-
01
4.09
e
-
01
2
7.37
e
-
02
8.39
e
-
02
8.34
e
-
02
3
1.18
e
-
02
8.73
e
-
03
8.65
e
-
03
4
8.53
e
-
03
1.42
e
-
05
1.40
e
-
06
5
8.53
e
-
03
1.42
e
-
05
4.37
e
-
09
6
8.53
e
-
03
1.42
e
-
05
4.56
e
-
09
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
3.98
e
-
01
4.10
e
-
01
4.07
e
-
01
2
7.64
e
-
02
8.40
e
-
02
8.28
e
-
02
3
1.01
e
-
02
8.76
e
-
03
8.55
e
-
03
4
7.67
e
-
03
3.55
e
-
05
1.38
e
-
06
5
7.67
e
-
03
3.55
e
-
05
1.60
e
-
09
6
7.67
e
-
03
3.55
e
-
05
1.45
e
-
09
1
Chebyshev polynomials of second Kind(θ,ϕ)=(1/2,1/2)
3.81
e
-
01
4.10
e
-
01
4.09
e
-
01
2
7.41
e
-
02
8.39
e
-
02
8.34
e
-
02
3
1.12
e
-
02
8.73
e
-
03
8.65
e
-
03
4
7.85
e
-
03
1.29
e
-
05
1.41
e
-
06
5
7.85
e
-
03
1.29
e
-
05
2.38
e
-
09
6
7.85
e
-
03
1.29
e
-
05
2.23
e
-
09
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
3.78
e
-
01
4.10
e
-
01
4.09
e
-
01
2
7.33
e
-
02
8.38
e
-
02
8.34
e
-
02
3
1.25
e
-
02
8.73
e
-
03
8.65
e
-
03
4
9.28
e
-
03
1.58
e
-
05
1.41
e
-
06
5
9.28
e
-
03
1.58
e
-
05
3.63
e
-
10
6
9.28
e
-
03
1.58
e
-
05
3.40
e
-
10Figure 8
Solution and error components for nonlinear Klein Gordon equation.


#### 5.2.6. 2D Nonlinear Wave Equation

The 2D nonlinear wave equation is(39)uttx,y,t-c2uxxx,y,t+uyyx,y,t+fu=pt,x,y,x,y,t∈Dx×Dy×Dt,where(40)Dx=ax,bx,Dy=ay,by,Dt=0,tf,with the initial and boundary conditions (41)uax,y,t=ϱ11y,t,ubx,y,t=ϱ12y,t,ux,ay,t=ϱ21x,t,ux,by,t=ϱ22x,t,ux,y,0=ζx,y.By discretization using J-GL-C methods we obtain (42)Fu=Stt⊗Ix⊗Iyu-c2It⊗Sxx⊗Iy+It⊗Ix⊗Syyu+fu-p=0,F′u=Stt⊗Ix⊗Iy-c2It⊗Sxx⊗Iy+It⊗Ix⊗Syy+diagf′u,where u=[u1,1,1,u1,1,2,…,u1,1,nt,…,unx,ny,1,unx,ny,2, …,u(nx,ny,nz)]T, Stt=2/tfQt2, Sxx=2/bx-axQx2, and Syy=2/by-ayQy2. The solution of (39) with c=1, f(u)=u3 and appropriate initial and boundary conditions is e-tsin⁡(x+y). We took a domain (x,y,t)∈[-1,1]×[-1,1]×[0,2] and solved the system of nonlinear equations using a single application of NMIM with initial guess u=0 with the initial and boundary conditions introduced. Table 14 gives the norm of the error as a function of m for several pairs of values for θ and ϕ and several grids. The collocation method with 20×20×20 which worked best was the Chebyshev collocation method of first kind. Figure 9 gives the error components for that collocation method and that grid.Table 14
Norm of the error for the 2D nonlinear wave equation as a function ofm for several pair of values for θ and ϕ and several grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
5
×
20
10
×
10
×
20
20
×
20
×
20
1
Legendre polynomials(θ,ϕ)=(0,0)
4.21
e
-
02
Divergent
Divergent
2
3.85
e
-
03
3
2.03
e
-
04
9
1.09
e
-
04
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
4.45
e
-
02
5.12
e
-
02
5.14
e
-
02
2
3.50
e
-
03
3.74
e
-
03
1.53
e
-
02
5
7.34e-04
5.74
e
-
08
2.16
e
-
04
7
7.34
e
-
04
1.40
e
-
09
3.76
e
-
06
9
7.34
e
-
04
1.40
e
-
09
6.56
e
-
08
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
4.40
e
-
02
5.13
e
-
02
Divergent
2
3.43
e
-
03
3.81
e
-
03
3
5.98
e
-
04
1.62
e
-
04
4
4.39
e
-
04
8.01
e
-
09
5
4.39
e
-
04
7.64
e
-
10
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
3.88
e
-
02
5.02
e
-
02
5.13
e
-
02
2
3.91
e
-
03
3.85
e
-
03
3.80
e
-
03
3
5.94
e
-
04
1.54
e
-
04
1.79
e
-
03
4
6.60
e
-
04
4.54
e
-
09
4.28
e
-
05
5
6.60
e
-
04
6.78
e
-
10
7.94
e
-
07
9
6.60
e
-
04
6.78
e
-
10
5.39
e
-
11Figure 9
Error components for the nonlinear 2D wave equation.


#### 5.2.7. 3D Nonlinear Wave Equation

The 3D nonlinear wave equation is(43)uttx,y,z,t-c2uxxx,y,z,t+uyyx,y,z,t+uzzx,y,z,t+fu=px,y,z,x,y,z,t∈Dx×Dy×Dz×Dt,where(44)Dx=ax,bx,Dy=ay,by,Dz=az,bz,Dt=0,tf,with the initial and boundary conditions (45)uax,y,z,t=ϱ11y,z,t,ubx,y,z,t=ϱ12y,z,t,ux,ay,z,t=ϱ21x,z,t,ux,by,z,t=ϱ22x,z,t,ux,y,az,t=ϱ31x,y,t,ux,y,bz,t=ϱ32x,y,t,ux,y,z,0=ζx,y,z.By discretization using J-GL-C methods we obtain (46)Fu=Stt⊗Ix⊗Iy⊗Izu-c2It⊗Sxx⊗Iy⊗Iz+It⊗Ix⊗Syy⊗Iz+It⊗Ix⊗Iy⊗Szzu+fu-p=0,F′u=Stt⊗Ix⊗Iy⊗Iz-c2It⊗Sxx⊗Iy⊗Iz+It⊗Ix⊗Syy⊗Iz+It⊗Ix⊗Iy⊗Szz+diagf′u,where u=[u(1,1,1,1),u(1,1,1,2),…,u(1,1,1,nt),…,u(nx,ny,nz,1), u(nx,ny,nz,2),…,u(nx,ny,nz,nt)]T, Stt=2/tfQt2, Sxx=2/bx-axQx2, Syy=2/by-ayQy2, and Szz=2/bz-azQz2. The solution of (43) with c=0.8, f(u)=u2 and appropriate initial and boundary conditions is e-tsin⁡(x+y+z). We took a domain (x,y,z,t)∈[0,1]×[0,1]×[0,1]×[0,1] and solved the system of nonlinear equations using a single application of NMIM with initial guess u=0 and with initial and boundary conditions introduced. Table 15 gives the norm of the error as a function of m for several pairs of values for θ,  ϕ and several grids. The collocation method which worked best for the grid 7×7×7×12 was the Legendre collocation method. Figure 10 gives the error components for that collocation method and that grid.Table 15
Norm of the error for the 3D nonlinear wave equation as a function ofm for several values for θ and ϕ and grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
5
×
5
×
12
6
×
6
×
6
×
12
7
×
7
×
7
×
12
1
Legendre polynomials(θ,ϕ)=(0,0)
9.95
e
-
01
1.00
e
+
00
1.14
e
+
00
3
1.14
e
-
02
8.53
e
-
03
1.09
e
-
02
5
1.25
e
-
01
1.12
e
-
01
1.73
e
-
01
7
5.53
e
-
06
4.12
e
-
07
1.73
e
-
06
8
5.54
e
-
06
7.54
e
-
08
4.24
e
-
09
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
1.00
e
+
00
1.00
e
+
00
divergent
2
1.00
e
-
01
8.68
e
-
02
4
1.31
e
+
00
1.29
e
+
00
7
2.24
e
-
05
1.99
e
-
06
8
2.24
e
-
05
2.11
e
-
06
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
9.98
e
-
01
9.99
e
-
01
1.32
e
+
00
2
1.02
e
-
01
8.06
e
-
02
1.08
e
-
01
5
1.22
e
-
01
1.18
e
-
01
2.46
e
-
01
7
1.81
e
-
05
1.84
e
-
06
3.72
e
-
05
8
1.81
e
-
05
1.76
e
-
06
4.87
e
-
08
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
9.95
e
-
01
divergent
1.09
e
+
00
2
9.98
e
-
02
1.04
e
-
01
5
1.31
e
-
01
1.51
e
-
01
7
1.92
e
-
05
4.37
e
-
07
8
1.92
e
-
05
3.27
e
-
08Figure 10
Error components for the nonlinear 3D wave equation.


#### 5.2.8. 3D Nonlinear Poisson Equation

The 3D nonlinear Poisson equation is(47)uxxx,y,z+uyyx,y,z+uzzx,y,z+fu=px,y,z,x,y,z∈Dx×Dy×Dz,where(48)Dx=ax,bx,Dy=ay,by,Dz=az,bz,with the initial and boundary conditions (49)uax,y,z=ϱ11y,z,ubx,y,z=ϱ12y,z,ux,ay,z=ϱ21x,z,ux,by,z=ϱ22x,z,ux,y,az=ϱ31x,y,ux,y,bz=ϱ32x,y.By discretization using J-GL-C methods we obtain (50)Fu=Sxx⊗Iy⊗Iz+Ix⊗Syy⊗Iz+Ix⊗Iy⊗Szzu+fu-p=0,F′u=Sxx⊗Iy⊗Iz+Ix⊗Syy⊗Iz+Ix⊗Iy⊗Szz+diagf′u,where u=[u1,1,1,u1,1,2,…,u1,1,nz,…,unx,ny,1,unx,ny,2, …,u(nx,ny,nz)]T, Sxx=2/bx-axQx2, Syy=2/by-ayQy2, and Szz=2/bz-azQz2. The solution of (47) with f(u)=u4 with appropriate initial and boundary conditions is sin⁡(x+y+z). We took a domain (x,y,z)∈[0,1]×[0,1]×[0,1] and solved the system of nonlinear equations using a single application of NMIM with initial guess u=0 and with initial and boundary conditions introduced. Table 16 gives the norm of the error as a function of m for several pairs of values for θ and ϕ and several grids. The Legendre collocation method offers best accuracy over different grids. Figure 11 gives the error components for that collocation method and that grid.Table 16
Norm of the error for the 3D Nonlinear Poisson equation as a function ofm for several values for θ and ϕ and grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
8
×
8
×
8
10
×
10
×
10
12
×
12
×
12
1
Legendre polynomials(θ,ϕ)=(0,0)
4.29
e
-
02
4.44
e
-
02
4.52
e
-
02
4
5.54
e
-
08
6.03
e
-
08
6.30
e
-
08
5
2.85
e
-
11
8.77
e
-
12
9.27
e
-
12
6
2.95
e
-
11
1.51
e
-
14
6.44
e
-
15
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
4.37
e
-
02
4.51
e
-
02
4.57
e
-
02
3
6.80
e
-
04
7.09
e
-
04
7.24
e
-
04
5
1.33
e
-
09
9.88
e
-
12
9.75
e
-
12
6
1.33
e
-
09
9.39
e
-
13
7.66
e
-
15
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
4.33
e
-
02
4.46
e
-
02
4.53
e
-
02
2
5.39
e
-
03
5.63
e
-
03
5.76
e
-
03
4
5.58
e
-
08
6.10
e
-
08
6.35
e
-
08
6
9.59
e
-
10
6.78
e
-
13
1.88
e
-
14
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
4.23
e
-
02
4.41
e
-
02
4.50
e
-
02
3
6.33
e
-
04
6.78
e
-
04
7.01
e
-
04
5
5.63
e
-
10
8.53
e
-
12
9.17
e
-
12
6
5.63
e
-
10
3.26
e
-
13
8.55
e
-
15Figure 11
Error components for the nonlinear 3D Poisson equation.


#### 5.2.9. The Nonlinear Complex Generalized Zakharov System

The nonlinear complex Zakharov system has importance in plasma physics [3]. The system includes two coupled nonlinear PDEs which can be written as(51)i∂tψx,t+δ1∂ttψx,t-δ2ψx,twx,t+δ3ψx,t2ψx,t=0,∂ttwx,t-cs2∂xxwx,t-δ4∂xxψx,t2=0,x,t∈ax,bx×at,bt,subject to the initial and boundary conditions(52)ψax,t=ψ1t,ψbx,t=ψ2t,ψx,0=ψ0x,wax,t=w1t,wbx,t=w2x,wx,0=w3x,x,t∈ax,bx×at,bt.Several numerical methods have been proposed recently for approximating the solution of (51) and (52) such as the homotopy method [24], the finite difference method [25, 26], and the variational iteration method [27]. Also, Bao et al. [28] suggested some high-accurate numerical methods for solving numerically (51) and (52). Bao and Sun [29] applied a new technique based on time-splitting discretization for approximating the solution of a variant of (51) and (52).One can split (51) using the real and imaginary parts of ψ(x,t), u(x,t), and v(x,t), as(53)∂tux,t+δ1∂xxvx,t-δ2vx,twx,t+δ3u2x,t+v2x,tvx,t=0,-∂tvx,t+δ1∂xxux,t-δ2ux,twx,t+δ3u2x,t+v2x,tux,t=0,∂ttwx,t-cs2∂xxwx,t-2δ4ux,t∂xxux,t+∂xux,t2+vx,t∂xxvx,t+∂xvx,t2=0,with the initial and boundary conditions(54)uax,t=α1t,ubx,t=α2t,vax,t=α3t,vbx,t=α4t,wax,t=α5t,wbx,t=α6t,ux,at=β1x,vx,at=β2x,wx,at=β3x,wtx,at=β4x,x,t∈ax,bx×at,bt.The matrix form of the nonlinear system (53) is(55)∂tδ1∂xx0δ1∂xx-∂t000∂tt-cs2∂xxuvw+q1u,v,wq2u,v,wq3u,v,w=000,where (56)q1=-δ2vx,twx,t+δ3u2x,t+v2x,tvx,t,q2=-δ2ux,twx,t+δ3u2x,t+v2x,tux,t,q3=-2δ4ux,t∂xxux,t+∂xux,t2+vx,t∂xxvproblemx,t+∂xvx,t2,where the constants δ1, δ2, δ3, δ4, and cs are given and the functions αi(t), 1≤i≤6 and βj(x), 1≤j≤4 are known. We use J-GL-C methods for discretizing (55) subject to the initial and boundary conditions (54) to reduce (55) to a system of nonlinear algebraic equations.


#### 5.2.10. Complex Generalized Zakharov Equation I

The first complex generalized Zakharov problem [3] is (57)i∂tψ-∂xxψ-ψw=0,∂ttw-∂xxw+∂xxψ2=0,with domain [-1,1]×[0,3.0]. The nonlinear problem (57) has the analytical solution (58)ψx,t=ux,t+ivx,t=3eix+ttanh⁡12x+2t,wx,t=1-tanh212x+2t.We define(59)E=maxU-Uanalytical∞,V-Vanalytical∞,W-Wanalytical∞to check the numerical accuracy in the computed solutions. Problem (57) is discretized by using Chebyshev collocation method of first kind. We took the initial guess 0 with initial and boundary conditions introduced and a grid 50×20. The achieved numerical accuracy is shown in Table 17. Our computed solutions are better than those in [3] because the maximum accuracy in [3] is of the order of 10-7. The solutions and error components are plotted in Figures 12, 13, and 14.Table 17
Maximum absolute error in the solution of the system of nonlinear equations for complex generalized Zakharov equation and varyingm.
m
E
1
8.56
e
-
1
2
2.72
e
0
3
4.16
e
0
4
8.92
e
0
5
1.12
e
1
6
3.40
e
1
7
2.04
e
1
8
2.83
e
1
9
4.05
e
0
10
2.34
e
0
11
7.77
e
-
1
12
3.25
e
-
1
13
1.14
e
-
1
14
6.20
e
-
2
15
2.79
e
-
2
16
1.85
e
-
2
17
8.11
e
-
3
18
4.44
e
-
3
19
1.79
e
-
3
20
8.86
e
-
4
21
3.91
e
-
4
22
1.66
e
-
4
23
5.98
e
-
5
24
1.97
e
-
5
25
5.32
e
-
6
26
1.20
e
-
6
27
4.12
e
-
7
28
1.36
e
-
7
29
6.95
e
-
8
30
2.95
e
-
8
31
9.29
e
-
9
32
3.42
e
-
9
33
1.43
e
-
9
34
5.57
e
-
10
35
2.62
e
-
10
36
1.26
e
-
10
37
8.65
e
-
11Figure 12
Solution and error components for complex generalized Zakharov equation (57).Figure 13
Solution and error components for complex generalized Zakharov equation (57).Figure 14
Solution and error components for complex generalized Zakharov equation (57).


#### 5.2.11. Complex Generalized Zakharov Equation II

The second complex generalized Zakharov equation [3] is(60)i∂tψ+∂xxψ+2ψw-2ψ2=0,∂ttw-∂xxw+∂xxψ2=0,with the domain [-1,1]×[0,1]. That problem has the analytical solution (61)ψx,t=ux,t+ivx,t=32e-ix+3ttanh⁡x+2t,wx,t=-14tanh2x+2t.Problem (60) is discretized by using Chebyshev collocation method of first kind. We took the initial guess 0 with initial and boundary conditions introduced and a grid 21×28. The achieved numerical accuracy is shown in Table 18. Our computed solutions are better than those in [3] because the ones in [3] achieved a maximum accuracy of the order of 10-6. The solutions and error components are plotted in Figures 15, 16, and 17.Table 18
Maximum absolute error in the solution of the system of nonlinear equations for complex generalized Zakharov equation and varyingm.
m
E
1
9.2329
e
-
1
2
7.0559
e
-
1
3
5.726
e
-
1
4
6.3466
e
-
2
5
4.3852
e
-
3
6
2.3221
e
-
4
7
1.5501
e
-
5
8
4.2396
e
-
7
9
2.8906
e
-
8
10
7.7742
e
-
10
11
8.6019
e
-
11Figure 15
Solution and error components for complex generalized Zakharov equation (60).Figure 16
Solution and error components for complex generalized Zakharov equation (60).Figure 17
Solution and error components for complex generalized Zakharov equation (60).


#### 5.2.12. Murray Equation

Murray equation [2] is (62)ut=uxx+λ1uux+λ2u-λ3u2,x,t∈-1,1×0,1.The initial and boundary condition for Murray equation can be computed from its analytical solution (63)ux,t=λ22λ31+tanh⁡λ28λ322λ1λ3x+λ12+4λ32t.The Chebyshev collocation method of first kind is used to discretize Murray equation over the grid 16×14. The achieved numerical accuracy is shown in Table 19. The computed accuracy in [2] is of the order of 10-8. The solution and error components are shown in Figure 18.Table 19
Maximum absolute error in the solution of the system of nonlinear equations for Murray equation and varyingm.
m
E
1
2.63
e
-
01
2
1.33
e
-
01
3
3.85
e
-
02
4
1.47
e
-
04
5
7.43
e
-
07
6
2.15
e
-
09
7
7.53
e
-
12
8
5.59
e
-
14Figure 18
Solution and error components for Murray equation.


#### 5.2.13. Nonlinear Reaction and Diffusion Equation

Nonlinear reaction and diffusion equation is(64)ut=uxx+u1+u-u2,x,t∈-1,1×0,1.The analytical solution is (65)ux,t=12+52tanh⁡104x+54t.The Chebyshev collocation method of first kind is used to discretize nonlinear reaction and diffusion equation over the grid 22×14. The initial guess is 0 with initial and boundary conditions introduced. The achieved numerical accuracy is shown in Table 20. The computed accuracy in [2] is of the order of 10-8. The solution and error components are shown in Figure 19.Table 20
Maximum absolute error in the solution of the system of nonlinear equations for nonlinear reaction and diffusion equation and varyingm.
m
E
1
7.28
e
-
02
2
2.12
e
-
02
3
4.99
e
-
03
4
1.09
e
-
06
5
9.20
e
-
10
6
6.53
e
-
13
7
2.24
e
-
14Figure 19
Solution and error components for nonlinear reaction and diffusion equation.


#### 5.2.14. Reaction-Diffusion Equation with Fischer-Kolmogorov Reaction Term and Density Dependent Diffusion

Reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion is(66)ut=u2xx+u1-u,x,t∈-1,1×0,1.The analytical solution is (67)ux,t=21+tanh⁡x+t/4.The Chebyshev collocation method of first kind is used to discretize reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion over the grid 30×30. The initial guess is 0 with initial and boundary conditions introduced. The achieved numerical accuracy after 1 and 2 application of NMIM is shown in Table 21. The computed accuracy in [2] is only of the order of 10-9. The solution and error components are shown in Figure 20.Table 21
Maximum absolute error in the solution of the system of nonlinear equations for reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion and varyingm.
m (1 Iter.)
E
1
1.49
e
+
00
2
1.12
e
+
00
3
1.07
e
+
00
4
4.92
e
-
01
5
3.73
e
-
01
m (2 Iter.)
E
1
1.51
e
-
02
2
7.84
e
-
03
3
6.72
e
-
04
4
5.17
e
-
07
5
1.04
e
-
09
6
2.35
e
-
12
7
5.04
e
-
14Figure 20
Solution and error components for reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion.


#### 5.2.15. Two-Dimensional Nonlinear Reaction-Diffusion Equation

The two-dimensional nonlinear reaction-diffusion equation is(68)ut=12∇2u+u2u-1,x,y,t∈-1,1×-1,1×0,20.Its analytical solution is (69)ux,y,t=11+e1/2x+y-1/2t.The Chebyshev collocation method of first kind is used to discretize the two-dimensional nonlinear reaction-diffusion equation over grid 30×30. The initial guess is 0 with initial and boundary conditions introduced. The achieved numerical accuracy is shown in Table 22. The computed accuracy in [2] is only of the order of 10-9. The error components are shown in Figure 21.Table 22
Maximum absolute error in the solution of the system of nonlinear equations for two-dimensional nonlinear reaction-diffusion equation and varyingm.
m
E
1
7.82
e
-
02
2
1.19
e
-
02
3
3.63
e
-
03
4
4.04
e
-
05
5
9.89
e
-
07
6
2.65
e
-
08
7
7.48
e
-
10Figure 21
Error components for the two-dimensional nonlinear reaction-diffusion equation.


#### 5.2.16. Ostrovsky-Hunter Equation

Ostrovsky-Hunter equation is(70)ut-βuxxx+αuuxx=γu,where α is a nonlinear coefficient and β∈R, γ>0 are the dispersion coefficients. We assume a periodic boundary condition for this problem; that is, (71)u0,t=u2π,t,ux0,t=ux2π,t,uxx0,t=uxx2π,t,uxxx0,t=uxxx2π,t.The Chebyshev collocation method of first kind is used to discretize Ostrovsky-Hunter equation over a grid30×100. The initial guess is 0 with initial and boundary conditions introduced. The residue is shown in Table 23. Our solution for the grid {0,6π/11,π,16π/11}×{0.1,0.2,0.3,0.4} is given in Table 24. We notice that our results do not agree with those presented in [4]. Our solution is plotted in Figure 22.Table 23
Residue of the system of nonlinear equations for Ostrovsky-Hunter equation and varyingm.
m
F
(
x
)
∞
1
2
2
2.5963
3
2.797
4
0.72034
5
0.15165
6
0.025873
7
0.002996
8
0.00020443
9
8.7968
e
-
06
10
1.4609
e
-
07
11
1.3542
e
-
08
12
8.8853
e
-
10
13
8.5753
e
-
11Table 24
Solution for Ostrovsky-Hunter equation.
x
t
0.1
0.2
0.3
0.4
0
0.094587
0.17222
0.20233
0.21608
6
π
/
11
0.96556
0.94696
0.8702
0.73399
π
−0.10438
−0.22921
−0.41929
−0.67358
16
π
/
11
−0.98195
−0.96996
−0.97341
−0.97748Figure 22
Solution for Ostrovsky-Hunter equation.


## 5.1. Computational Verification of Convergence Order

In general, it is hard to verify the convergence order by solving initial and boundary value problems. Therefore, we consider the following five systems of nonlinear equationsF(x)=0. (20)Problem1=F1x=x2x3+x4x2+x3,F2x=x1x3+x4x1+x3,F3x=x1x2+x4x1+x2,F4x=x1x2+x3x1+x2-1,Problem2=Fix=xi2xi+1-1,i=1,2,…,n-1,Fix=xnx1-1,i=n,Problem3=F1x=3-0.5x1x1-2x2+1,F2x=3-0.5xnxn-2xn-1+1,Fix=3-0.5xixi-xi-1+2xi+1+1,i=2,3,…,n-1,Problem4=Fix=2n+i1-cos⁡xi-sin⁡xi-∑j=1ncos⁡xj2sin⁡xi-cos⁡xi,i=1,2,…,n,Problem5=F1x=10x1+sin⁡x1+x2-1,F2x=8x2-cos2x3-x2-1,F3x=12x3+sin⁡x3-1.The system of nonlinear equations is solved by Mathematica using 11,000 digits. Table4 shows that CCOs agree with theoretical convergence orders.Table 4
Verification of convergence order.
F
(
x
k
)
∞
Iter.
m
=
4
m
=
5
m
=
6
m
=
7
Problem 1 (n=4)
1
6.56
e
-
8
1.09
e
-
11
1.82
e
-
15
3.04
e
-
19
2
6.31
e
-
68
5.74
e
-
158
8.66
e
-
286
2.17
e
-
451
3
4.44
e
-
608
6.90
e
-
2206
6.32
e
-
5422
6.67
e
-
10823
CCO
9
14
19
24
Problem 2 (n=10)
1
7.31
e
-
3
6.45
e
-
4
5.83
e
-
5
5.29
e
-
6
2
1.54
e
-
21
1.38
e
-
47
5.42
e
-
84
8.20
e
-
131
3
1.31
e
-
189
5.63
e
-
659
1.35
e
-
1585
3.06
e
-
3126
CCO
9
14
19
24
Problem 3 (n=100)
1
7.09
e
-
5
2.12
e
-
7
6.36
e
-
10
1.90
e
-
12
2
5.45
e
-
43
1.69
e
-
94
1.53
e
-
166
1.64
e
-
258
3
6.87
e
-
348
8.80
e
-
1143
7.95
e
-
2673
8.85
e
-
5184
CCO
8
12
16
20
Problem 4 (n=5)
2
1.24
e
-
3
8.47
e
-
6
8.17
e
-
9
9.64
e
-
13
3
2.73
e
-
19
2.08
e
-
54
1.45
e
-
120
2.27
e
-
229
4
6.73
e
-
148
1.79
e
-
645
7.35
e
-
1920
5.31
e
-
4575
CCO
8.21
12.2
16.1
20.1
Problem 5 (n=3)
1
3.77
e
-
6
8.57
e
-
10
1.94
e
-
13
4.39
e
-
17
2
2.90
e
-
67
1.71
e
-
151
6.21
e
-
269
1.46
e
-
419
3
3.09
e
-
581
2.79
e
-
1984
1.71
e
-
4737
9.99
e
-
9298
CCO
8.41
12.9
17.5
22.1
Theoretical convergence order (5m-11)
9
14
19
24


## 5.2. Solving Initial and Boundary Value Problems

In this section, we solve different initial and boundary value problems by using different J-GL-C discretization methods. In all computations, we make spatial and temporal discretizations to translate the entire IVPs and BVPs into systems of nonlinear equations. The systems of nonlinear equations are then solved by using the proposed iterative method NMIM. In most of the problems, we use vector0 with initial and boundary conditions introduced in it as the initial guess. The inclusion of initial and boundary conditions in 0 makes the initial guess nonsmooth but it works for most of the problems. When the nonsmooth initial guess does not work, we then make it smooth to get convergence.


### 5.2.1. The Lane-Emden Equation

The Lane-Emden equation is(21)x′′t+2tx′t+xtp=0,x′0=0,x0=1.J-GL-C discretization methods translate the Lane-Emden equation into the following system of nonlinear equations:(22)Fx=Sxxx+diag2tSxx+xp=0,F′x=Sxx+diag2tSx+pdiagxp-1,where diag(·) denotes a diagonal matrix with values in the diagonal equal to the elements of its vector argument, Sx=2/bQx, Sxx=2/bQx2, and [0,b] is the domain of the problem. We have solved the Lane-Emden equation for different values of p ranging from two to five and b=3. Figure 1 shows the initial guess and the computed numerical solutions for θ=1/2, ϕ=-1/2, and 50 grid points. Table 5 gives the norm of the residue in the solution of the nonlinear systems by a single iteration of NMIM with varying m.Table 5
Norm of the residue of the system of nonlinear equation resulting from Lane-Emden equation with varyingm.
F
(
x
k
)
∞
m
p
2
3
4
5
10
2.69
e
-
17
3.70
e
-
14
1.04
e
-
12
1.03
e
-
11
20
1.36
e
-
46
7.06
e
-
41
9.33
e
-
39
1.22
e
-
37
30
9.33
e
-
76
1.68
e
-
67
7.56
e
-
65
1.97
-
63
40
6.16
e
-
105
4.53
e
-
94
6.19
e
-
91
2.26
e
-
89
50
4.08
e
-
134
1.21
e
-
120
5.11
e
-
117
2.18
e
-
115
60
2.71
e
-
163
3.22
e
-
147
4.24
e
-
143
2.00
e
-
141
70
1.79
e
-
192
8.555
e
-
174
3.52
e
-
169
1.83
e
-
167
80
1.19
e
-
221
2.27
e
-
200
2.92
e
-
195
1.69
e
-
193
90
7.89
e
-
251
6.03
e
-
227
2.43
e
-
221
1.58
e
-
219
100
5.23
e
-
280
1.60
e
-
253
2.02
e
-
247
1.47
e
-
245Figure 1
Initial guess and solutions of the Lane-Emden equation for different values ofp.


### 5.2.2. Bratu Problem

Bratu problem is the boundary value problem:(23)x′′t+αext=0,x0=0,x1=0,where α is a parameter. The discretization is carried out using J-GL-C methods. The resulting system of nonlinear equations is (24)Fx=Sxxx+αex=0,F′x=Sxx+αdiagex,where Sxx=2Qx2. We considered several values of α, values θ=-1/2 and ϕ=-1/2 for the parameters of J-GL-C, 50 grid points, and solved the resulting system of nonlinear equations using a single application of NMIM with varying m. The initial guess and computed numerical solutions are given in Figure 2. Table 6 gives the norm of the residue of the solution of the system of nonlinear equations for varying m.Table 6
Norm of the residue of the solution of the system of nonlinear equations for Bratu problem and varyingm.
F
(
x
k
)
∞
m
α
1
2
3
10
7.68
e
-
41
4.10
e
-
31
1.97
e
-
28
20
2.20
e
-
93
1.01
e
-
70
9.65
e
-
64
30
7.02
e
-
146
2.50
e
-
110
4.72
e
-
99
40
2.23
e
-
198
6.17
e
-
150
2.30
e
-
134
50
7.13
e
-
251
1.52
e
-
189
1.12
e
-
169
60
2.27
e
-
303
3.75
e
-
229
5.52
e
-
205
70
7.30
e
-
356
9.27
e
-
269
2.69
e
-
240
80
2.30
e
-
408
2.28
e
-
308
1.32
e
-
275
90
7.40
e
-
461
5.70
e
-
348
6.50
e
-
311
100
2.30
e
-
513
1.40
e
-
387
3.20
e
-
346Figure 2
Initial guess and solutions of Bratu problem for different values ofα.


### 5.2.3. Frank-Kamenetskii Problem

Frank-Kamenetskii problem is(25)x′′t+1tx′t+αext=0,x′0=0,x1=0,where α is a parameter. After discretization by J-GL-C methods we get the system of nonlinear equations:(26)Fx=Sxxx+diag1tSxx+αex=0,F′x=Sxx+diag1tSx+αdiagex.We considered several values of α, values θ=-1/2 and ϕ=-1/2 for the parameters of J-GL-C, 50 grid points, and solved the resulting system of nonlinear equations by using a single application of NMIM with varying m. Figure 3 plots the initial guess and computed solutions and Table 7 gives the norm of the residue of the solution of the system of nonlinear equations.Table 7
Norm of the residue of the solution of the system of nonlinear equations for Fran-Kamenetskii problem and varyingm.
F
(
x
k
)
∞
m
α
1
1.1
1.2
1.3
10
7.09
e
-
41
1.89
e
-
39
2.63
e
-
38
1.80
e
-
37
20
2.46
e
-
82
2.36
e
-
79
6.21
e
-
77
4.05
e
-
75
30
8.58
e
-
124
2.95
e
-
119
1.46
e
-
115
9.09
e
-
113
40
2.98
e
-
165
3.69
e
-
159
3.45
e
-
154
2.03
e
-
150
50
1.03
e
-
206
4.61
e
-
199
8.15
e
-
193
4.57
e
-
188
60
3.61
e
-
248
5.76
e
-
239
1.92
e
-
231
1.02
e
-
225
70
1.25
e
-
289
7.19
e
-
279
4.53
e
-
270
2.29
e
-
263
80
4.40
e
-
331
9.00
e
-
319
1.10
e
-
308
5.15
e
-
301
90
1.50
e
-
372
1.10
e
-
358
2.50
e
-
347
1.20
e
-
338
100
5.30
e
-
414
1.40
e
-
398
6.00
e
-
386
2.60
e
-
376Figure 3
Initial guess and solutions of Frank-Kamenetskii problem with different values ofα.


### 5.2.4. 1 + 1 Nonlinear Schrödinger Equation

1 + 1 nonlinear Schrödinger equation is(27)i∂tψx,t+∂xxψx,t+2γψx,t2ψx,t-2δRx,tψx,t=0,x,t∈Dx×Dt,where(28)Dx=ax,bx,Dt=0,tf,with the initial and boundary conditions (29)ψax,t=η1t,ψbx,t=η2t,ψx,0=ζx.The function ψ(x,t) is a complex function and can be written as ψ(x,t)=ψ1(x,t)+iψ2(x,t). The Schrödinger equation can be rewritten in terms of the real functions ψ1(x,t) and ψ2(x,t) as(30)∂t∂xx∂xx-∂tψ1ψ2-2δRx,t0110ψ1ψ2+2γψ12+ψ220110ψ1ψ2=0,with initial and boundary conditions (31)ψ1ax,t+iψ2ax,t=η11t+iη12,ψ1bx,t+iψ2bx,t=η21t+iη22t,ψ1x,0+iψ2x,0=ζ1x+iζ2x.Using J-GL-C methods, we can discretize (30) as(32)Fψ=St⊗IxIt⊗SxxIt⊗Sxx-St⊗Ixψ1ψ2-2δOdiagRdiagROψ1ψ2+2γOdiagψ12+ψ22diagψ12+ψ22Oψ1ψ2=0,F′ψ=St⊗IxIt⊗SxxIt⊗Sxx-St⊗Ix-2δOdiagRdiagRO+2γ2diagψ1⊙ψ23diagψ223diagψ122diagψ1⊙ψ2,where ψj=[ψj1,1,ψj1,2,…,ψj1,nt,…,ψjnx,1,ψjnx,2, …,ψj(nx,nt)]T, j=1,2, R=[R1,1,R1,2,…,R1,nt, …,Rnx,1,Rnx,2,…,R(nx,nt)]T, ⊗ is the Kronecker product, ⊙ is the element-wise multiplication, O is a matrix of all zeros, I is the identity matrix, St=2/tfQt, Sxx=2/bx-axQx2, and Q is the J-GL-C operational matrix. The system of nonlinear equations (32) can be rewritten as (33)Fψ=Aψ+gψ-p,where A=St⊗IxIt⊗SxxIt⊗Sxx-St⊗Ix-2δOdiag(R)diag(R)O, g(ψ)=2γ2diag(ψ1⊙ψ2)3diag(ψ22)3diag(ψ12)2diag(ψ1⊙ψ2), and p is the vector incorporating the initial-boundary conditions. Four nonlinear Schrödinger equations are listed in Table 8 with their corresponding analytical solutions. We chose a domain (x,t)∈[-1,1]×[0,1] and solved the equations using different pair of values for the parameters θ and ϕ of J-GL-C and different grids for the domain. As initial guess, we used ψ1=0, ψ2=0 with boundary conditions set for problem 3. For problems 1, 2, and 4 NMIM does not converge with that initial guess and we made the initial guess smoother. We smoothed the initial guess by applying the iteration(34)ψ=-A-12gψ-ponce for problem 1 and twice for problems 2 and 4. Tables 9, 10, 11, and 12 give the norm of the error in the solution when the nonlinear system of equations is solved using a single application of NMIM with varying m. Figures 4–7 display the solutions and the error in the components of ψ1 and ψ2 for θ=-4/2, ϕ=-1/2, and a grid 15×30.Table 8
Nonlinear Schrödinger equations.
Schrödinger equation
Analytical Solution
i∂tψ+∂xxψ+2|ψ|2ψ=0
ψ
(
x
,
t
)
=
e
i
(
x
+
t
)
i
∂
t
ψ
+
∂
x
x
ψ
-
2
|
ψ
|
2
ψ
=
0
ψ
(
x
,
t
)
=
e
i
(
x
-
3
t
)
i∂tψ+12∂xxψ-|ψ|2ψ-cos⁡(x)2ψ=0
ψ
(
x
,
t
)
=
e
-
3
i
t
/
2
sin
⁡
(
x
)
i
∂
t
ψ
+
∂
x
x
ψ
+
2
|
ψ
|
2
ψ
-
2
δ
ψ
=
0
ψ
(
x
,
t
)
=
e
i
(
x
+
(
1
-
2
δ
)
t
)Table 9
Norm of the error in the solution for equation (17).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
5.49
e
-
01
6.21
e
-
01
6.41
e
-
01
4
5.36
e
-
02
2.35
e
-
01
2.44
e
-
01
8
2.85
e
-
04
4.55
e
-
04
4.97
e
-
04
12
2.85
e
-
04
5.81
e
-
08
6.38
e
-
08
16
2.85
e
-
04
5.69
e
-
10
4.38
e
-
12
20
2.85
e
-
04
6.13
e
-
11
5.95
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
6.86
e
-
01
Divergent
6.40
e
-
01
4
1.96
e
-
01
2.48
e
-
01
8
1.64
e
-
03
4.94
e
-
04
12
1.59
e
-
03
6.16
e
-
08
16
1.59
e
-
03
6.02
e
-
12
20
1.59
e
-
03
2.16
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
5.24
e
-
01
6.25
e
-
01
6.41
e
-
01
4
5.90
e
-
02
2.39
e
-
01
2.46
e
-
01
8
1.14
e
-
03
5.05
e
-
04
5.40
e
-
04
12
1.14
e
-
03
7.73
e
-
08
8.00
e
-
08
16
1.14
e
-
03
1.90
e
-
09
5.28
e
-
12
20
1.14
e
-
03
1.90
e
-
09
4.99
e
-
13
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
5.95
e
-
01
6.20
e
-
01
6.41
e
-
01
4
3.73
e
-
02
2.37
e
-
01
2.45
e
-
01
8
1.52
e
-
03
4.24
e
-
04
4.81
e
-
04
12
1.52
e
-
03
5.22
e
-
08
5.81
e
-
08
16
1.52
e
-
03
4.90
e
-
10
1.71
e
-
12
20
1.52
e
-
03
4.90
e
-
10
9.31
e
-
14Table 10
Norm of the error in the solution for equation (30).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
1.28
e
+
00
1.34
e
+
00
1.34
e
+
00
4
4.27
e
-
01
9.65
e
-
01
9.72
e
-
01
8
1.25
e
-
03
2.15
e
-
03
2.18
e
-
03
12
1.46
e
-
04
9.17
e
-
07
1.05
e
-
06
16
1.46
e
-
04
3.32
e
-
10
4.61
e
-
10
20
1.46
e
-
04
2.90
e
-
11
1.85
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
1.27
e
+
00
1.35
e
+
00
1.34
e
+
00
4
4.34
e
-
01
9.46
e
-
01
9.72
e
-
01
8
1.63
e
-
03
2.18
e
-
03
2.18
e
-
03
12
1.52
e
-
03
3.89
e
-
07
8.77
e
-
07
16
1.52
e
-
03
5.28
e
-
09
1.95
e
-
10
20
1.52
e
-
03
3.73
e
-
09
1.34
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
1.26
e
+
00
1.34
e
+
00
1.34
e
+
00
4
3.98
e
-
01
9.24
e
-
01
9.40
e
-
01
8
2.63
e
-
03
1.35
e
-
03
1.63
e
-
03
12
1.00
e
-
03
2.63
e
-
06
2.85
e
-
06
16
1.00
e
-
03
2.35
e
-
08
1.65
e
-
08
20
1.00
e
-
03
7.88
e
-
10
6.29
e
-
11
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
1.33
e
+
00
1.34
e
+
00
1.34
e
+
00
4
3.90
e
-
01
9.80
e
-
01
9.80
e
-
01
8
1.03
e
-
03
2.15
e
-
03
2.19
e
-
03
12
1.03
e
-
03
6.56
e
-
07
7.42
e
-
07
16
1.03
e
-
03
3.08
e
-
10
1.57
e
-
10
20
1.03
e
-
03
2.95
e
-
10
5.77
e
-
14Table 11
Norm of the error in the solution for equation (32).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomial(θ,ϕ)=(0,0)
9.98
e
-
02
8.73
e
-
02
9.66
e
-
02
4
7.62
e
-
04
3.16
e
-
04
3.25
e
-
04
8
4.27
e
-
05
3.64
e
-
11
3.12
e
-
11
12
4.27
e
-
05
3.16
e
-
11
1.89
e
-
13
16
4.27
e
-
05
3.16
e
-
11
1.88
e
-
13
20
4.27
e
-
05
3.16
e
-
11
1.89
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
9.71
e
-
02
9.33
e
-
02
9.67
e
-
02
4
1.18
e
-
03
3.37
e
-
04
3.36
e
-
04
8
7.68
e
-
04
3.23
e
-
08
2.93
e
-
11
12
7.68
e
-
04
1.89
e
-
10
1.07
e
-
13
16
7.68
e
-
04
1.87
e
-
10
1.06
e
-
13
20
7.68
e
-
04
1.87
e
-
10
1.07
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
9.47
e
-
02
9.17
e
-
02
9.65
e
-
02
4
9.17
e
-
04
3.14
e
-
04
3.38
e
-
04
8
6.70
e
-
04
9.13
e
-
11
8.39
e
-
11
12
6.70
e
-
04
8.92
e
-
11
1.19
e
-
13
16
6.70
e
-
04
8.92
e
-
11
1.19
e
-
13
20
6.70
e
-
04
8.92
e
-
11
1.19
e
-
13
1
Chebyshev polynomial of first kind(θ,ϕ)=(-1/2,-1/2)
1.01
e
-
01
9.01
e
-
02
9.61
e
-
02
4
1.31
e
-
03
2.99
e
-
04
3.26
e
-
04
8
8.53
e
-
04
4.90
e
-
11
2.27
e
-
11
12
8.53
e
-
04
5.59
e
-
11
3.28
e
-
14
16
8.53
e
-
04
5.59
e
-
11
3.24
e
-
14
20
8.53
e
-
04
5.59
e
-
11
3.20
e
-
14Table 12
Norm of the error in the solution for equation (35).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
6.44
e
-
01
7.84
e
-
01
8.01
e
-
01
4
6.52
e
-
02
4.18
e
-
02
4.18
e
-
02
8
2.77
e
-
04
1.36
e
-
06
1.33
e
-
06
12
2.75
e
-
04
1.70
e
-
10
1.81
e
-
10
16
2.75
e
-
04
4.14
e
-
11
4.87
e
-
13
20
2.75
e
-
04
4.14
e
-
11
4.67
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
6.51
e
-
01
Divergent
8.08
e
-
01
4
4.52
e
-
02
4.15
e
-
02
8
1.50
e
-
03
1.25
e
-
06
12
1.51
e
-
03
1.17
e
-
10
16
1.51
e
-
03
1.94
e
-
13
20
1.51
e
-
03
1.99
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
7.48
e
-
01
7.85
e
-
01
8.04
e
-
01
4
5.30
e
-
02
4.24
e
-
02
4.22
e
-
02
8
1.02
e
-
03
2.20
e
-
06
2.61
e
-
06
12
1.02
e
-
03
2.67
e
-
09
1.00
e
-
09
16
1.02
e
-
03
2.25
e
-
09
4.86
e
-
13
20
1.02
e
-
03
2.25
e
-
09
1.96
e
-
13
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
7.03
e
-
01
7.88
e
-
01
8.02
e
-
01
4
3.66
e
-
02
4.20
e
-
02
4.17
e
-
02
8
1.46
e
-
03
1.27
e
-
06
1.18
e
-
06
12
1.46
e
-
03
5.96
e
-
10
2.39
e
-
11
16
1.46
e
-
03
5.84
e
-
10
8.56
e
-
14
20
1.46
e
-
03
5.84
e
-
10
8.49
e
-
14Figure 4
Solution and error components for (17).Figure 5
Solution and error components for (30).Figure 6
Solution and error components for (32).Figure 7
Solution and error components for (35).


### 5.2.5. Klein Gordon Equation

Klein Gordon equation is(35)uttx,t-c2uxxx,t+ku-γux,t3=0,x,t∈Dx×Dt,where(36)Dx=ax,bx,Dt=0,tf,with the initial and boundary conditions (37)uax,t=ϱ1t,ubx,t=ϱ2t,ux,0=ξ1x,utx,0=ξ2x.By discretization using J-GL-C methods we obtain (38)Fu=Sttu-c2Sxxu+ku-γu3=0,F′u=Stt-c2Sxx+kI-3γdiagu2,where u=[u(1,1),u(1,2),…,u(1,nt),…,u(nx,1),u(nx,2),…, u(nx,nt)]T, Stt=2/tfQt2, and Sxx=2/bx-axQx2. The analytical solution of (35) is u(x,t)=δsech⁡κ(x-vt), where κ=k/c2-v2 and δ=2k/γ. We set the parameters c=1, γ=1, v=0.5, and k=0.5, picked up the domain (x,t)∈[-10,10]×[0,1], and solved the nonlinear system of equations by a single application of NMIM with initial guess u=0, with the boundary conditions introduced. We used several pairs of the θ, ϕ parameters of J-GL-C. According to Table 13, the Chebyshev collocation method of first kind exhibits best accuracy when the grid is more refined. The solution and the error components for θ=-1/2, ϕ=-1/2, nx=120, and nt=20 are plotted in Figure 8.Table 13
Norm of the error in the Klein Gordon equation as a function ofm for several values for θ and ϕ and grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
30
×
30
60
×
30
120
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
3.80
e
-
01
4.10
e
-
01
4.09
e
-
01
2
7.37
e
-
02
8.39
e
-
02
8.34
e
-
02
3
1.18
e
-
02
8.73
e
-
03
8.65
e
-
03
4
8.53
e
-
03
1.42
e
-
05
1.40
e
-
06
5
8.53
e
-
03
1.42
e
-
05
4.37
e
-
09
6
8.53
e
-
03
1.42
e
-
05
4.56
e
-
09
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
3.98
e
-
01
4.10
e
-
01
4.07
e
-
01
2
7.64
e
-
02
8.40
e
-
02
8.28
e
-
02
3
1.01
e
-
02
8.76
e
-
03
8.55
e
-
03
4
7.67
e
-
03
3.55
e
-
05
1.38
e
-
06
5
7.67
e
-
03
3.55
e
-
05
1.60
e
-
09
6
7.67
e
-
03
3.55
e
-
05
1.45
e
-
09
1
Chebyshev polynomials of second Kind(θ,ϕ)=(1/2,1/2)
3.81
e
-
01
4.10
e
-
01
4.09
e
-
01
2
7.41
e
-
02
8.39
e
-
02
8.34
e
-
02
3
1.12
e
-
02
8.73
e
-
03
8.65
e
-
03
4
7.85
e
-
03
1.29
e
-
05
1.41
e
-
06
5
7.85
e
-
03
1.29
e
-
05
2.38
e
-
09
6
7.85
e
-
03
1.29
e
-
05
2.23
e
-
09
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
3.78
e
-
01
4.10
e
-
01
4.09
e
-
01
2
7.33
e
-
02
8.38
e
-
02
8.34
e
-
02
3
1.25
e
-
02
8.73
e
-
03
8.65
e
-
03
4
9.28
e
-
03
1.58
e
-
05
1.41
e
-
06
5
9.28
e
-
03
1.58
e
-
05
3.63
e
-
10
6
9.28
e
-
03
1.58
e
-
05
3.40
e
-
10Figure 8
Solution and error components for nonlinear Klein Gordon equation.


### 5.2.6. 2D Nonlinear Wave Equation

The 2D nonlinear wave equation is(39)uttx,y,t-c2uxxx,y,t+uyyx,y,t+fu=pt,x,y,x,y,t∈Dx×Dy×Dt,where(40)Dx=ax,bx,Dy=ay,by,Dt=0,tf,with the initial and boundary conditions (41)uax,y,t=ϱ11y,t,ubx,y,t=ϱ12y,t,ux,ay,t=ϱ21x,t,ux,by,t=ϱ22x,t,ux,y,0=ζx,y.By discretization using J-GL-C methods we obtain (42)Fu=Stt⊗Ix⊗Iyu-c2It⊗Sxx⊗Iy+It⊗Ix⊗Syyu+fu-p=0,F′u=Stt⊗Ix⊗Iy-c2It⊗Sxx⊗Iy+It⊗Ix⊗Syy+diagf′u,where u=[u1,1,1,u1,1,2,…,u1,1,nt,…,unx,ny,1,unx,ny,2, …,u(nx,ny,nz)]T, Stt=2/tfQt2, Sxx=2/bx-axQx2, and Syy=2/by-ayQy2. The solution of (39) with c=1, f(u)=u3 and appropriate initial and boundary conditions is e-tsin⁡(x+y). We took a domain (x,y,t)∈[-1,1]×[-1,1]×[0,2] and solved the system of nonlinear equations using a single application of NMIM with initial guess u=0 with the initial and boundary conditions introduced. Table 14 gives the norm of the error as a function of m for several pairs of values for θ and ϕ and several grids. The collocation method with 20×20×20 which worked best was the Chebyshev collocation method of first kind. Figure 9 gives the error components for that collocation method and that grid.Table 14
Norm of the error for the 2D nonlinear wave equation as a function ofm for several pair of values for θ and ϕ and several grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
5
×
20
10
×
10
×
20
20
×
20
×
20
1
Legendre polynomials(θ,ϕ)=(0,0)
4.21
e
-
02
Divergent
Divergent
2
3.85
e
-
03
3
2.03
e
-
04
9
1.09
e
-
04
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
4.45
e
-
02
5.12
e
-
02
5.14
e
-
02
2
3.50
e
-
03
3.74
e
-
03
1.53
e
-
02
5
7.34e-04
5.74
e
-
08
2.16
e
-
04
7
7.34
e
-
04
1.40
e
-
09
3.76
e
-
06
9
7.34
e
-
04
1.40
e
-
09
6.56
e
-
08
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
4.40
e
-
02
5.13
e
-
02
Divergent
2
3.43
e
-
03
3.81
e
-
03
3
5.98
e
-
04
1.62
e
-
04
4
4.39
e
-
04
8.01
e
-
09
5
4.39
e
-
04
7.64
e
-
10
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
3.88
e
-
02
5.02
e
-
02
5.13
e
-
02
2
3.91
e
-
03
3.85
e
-
03
3.80
e
-
03
3
5.94
e
-
04
1.54
e
-
04
1.79
e
-
03
4
6.60
e
-
04
4.54
e
-
09
4.28
e
-
05
5
6.60
e
-
04
6.78
e
-
10
7.94
e
-
07
9
6.60
e
-
04
6.78
e
-
10
5.39
e
-
11Figure 9
Error components for the nonlinear 2D wave equation.


### 5.2.7. 3D Nonlinear Wave Equation

The 3D nonlinear wave equation is(43)uttx,y,z,t-c2uxxx,y,z,t+uyyx,y,z,t+uzzx,y,z,t+fu=px,y,z,x,y,z,t∈Dx×Dy×Dz×Dt,where(44)Dx=ax,bx,Dy=ay,by,Dz=az,bz,Dt=0,tf,with the initial and boundary conditions (45)uax,y,z,t=ϱ11y,z,t,ubx,y,z,t=ϱ12y,z,t,ux,ay,z,t=ϱ21x,z,t,ux,by,z,t=ϱ22x,z,t,ux,y,az,t=ϱ31x,y,t,ux,y,bz,t=ϱ32x,y,t,ux,y,z,0=ζx,y,z.By discretization using J-GL-C methods we obtain (46)Fu=Stt⊗Ix⊗Iy⊗Izu-c2It⊗Sxx⊗Iy⊗Iz+It⊗Ix⊗Syy⊗Iz+It⊗Ix⊗Iy⊗Szzu+fu-p=0,F′u=Stt⊗Ix⊗Iy⊗Iz-c2It⊗Sxx⊗Iy⊗Iz+It⊗Ix⊗Syy⊗Iz+It⊗Ix⊗Iy⊗Szz+diagf′u,where u=[u(1,1,1,1),u(1,1,1,2),…,u(1,1,1,nt),…,u(nx,ny,nz,1), u(nx,ny,nz,2),…,u(nx,ny,nz,nt)]T, Stt=2/tfQt2, Sxx=2/bx-axQx2, Syy=2/by-ayQy2, and Szz=2/bz-azQz2. The solution of (43) with c=0.8, f(u)=u2 and appropriate initial and boundary conditions is e-tsin⁡(x+y+z). We took a domain (x,y,z,t)∈[0,1]×[0,1]×[0,1]×[0,1] and solved the system of nonlinear equations using a single application of NMIM with initial guess u=0 and with initial and boundary conditions introduced. Table 15 gives the norm of the error as a function of m for several pairs of values for θ,  ϕ and several grids. The collocation method which worked best for the grid 7×7×7×12 was the Legendre collocation method. Figure 10 gives the error components for that collocation method and that grid.Table 15
Norm of the error for the 3D nonlinear wave equation as a function ofm for several values for θ and ϕ and grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
5
×
5
×
12
6
×
6
×
6
×
12
7
×
7
×
7
×
12
1
Legendre polynomials(θ,ϕ)=(0,0)
9.95
e
-
01
1.00
e
+
00
1.14
e
+
00
3
1.14
e
-
02
8.53
e
-
03
1.09
e
-
02
5
1.25
e
-
01
1.12
e
-
01
1.73
e
-
01
7
5.53
e
-
06
4.12
e
-
07
1.73
e
-
06
8
5.54
e
-
06
7.54
e
-
08
4.24
e
-
09
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
1.00
e
+
00
1.00
e
+
00
divergent
2
1.00
e
-
01
8.68
e
-
02
4
1.31
e
+
00
1.29
e
+
00
7
2.24
e
-
05
1.99
e
-
06
8
2.24
e
-
05
2.11
e
-
06
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
9.98
e
-
01
9.99
e
-
01
1.32
e
+
00
2
1.02
e
-
01
8.06
e
-
02
1.08
e
-
01
5
1.22
e
-
01
1.18
e
-
01
2.46
e
-
01
7
1.81
e
-
05
1.84
e
-
06
3.72
e
-
05
8
1.81
e
-
05
1.76
e
-
06
4.87
e
-
08
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
9.95
e
-
01
divergent
1.09
e
+
00
2
9.98
e
-
02
1.04
e
-
01
5
1.31
e
-
01
1.51
e
-
01
7
1.92
e
-
05
4.37
e
-
07
8
1.92
e
-
05
3.27
e
-
08Figure 10
Error components for the nonlinear 3D wave equation.


### 5.2.8. 3D Nonlinear Poisson Equation

The 3D nonlinear Poisson equation is(47)uxxx,y,z+uyyx,y,z+uzzx,y,z+fu=px,y,z,x,y,z∈Dx×Dy×Dz,where(48)Dx=ax,bx,Dy=ay,by,Dz=az,bz,with the initial and boundary conditions (49)uax,y,z=ϱ11y,z,ubx,y,z=ϱ12y,z,ux,ay,z=ϱ21x,z,ux,by,z=ϱ22x,z,ux,y,az=ϱ31x,y,ux,y,bz=ϱ32x,y.By discretization using J-GL-C methods we obtain (50)Fu=Sxx⊗Iy⊗Iz+Ix⊗Syy⊗Iz+Ix⊗Iy⊗Szzu+fu-p=0,F′u=Sxx⊗Iy⊗Iz+Ix⊗Syy⊗Iz+Ix⊗Iy⊗Szz+diagf′u,where u=[u1,1,1,u1,1,2,…,u1,1,nz,…,unx,ny,1,unx,ny,2, …,u(nx,ny,nz)]T, Sxx=2/bx-axQx2, Syy=2/by-ayQy2, and Szz=2/bz-azQz2. The solution of (47) with f(u)=u4 with appropriate initial and boundary conditions is sin⁡(x+y+z). We took a domain (x,y,z)∈[0,1]×[0,1]×[0,1] and solved the system of nonlinear equations using a single application of NMIM with initial guess u=0 and with initial and boundary conditions introduced. Table 16 gives the norm of the error as a function of m for several pairs of values for θ and ϕ and several grids. The Legendre collocation method offers best accuracy over different grids. Figure 11 gives the error components for that collocation method and that grid.Table 16
Norm of the error for the 3D Nonlinear Poisson equation as a function ofm for several values for θ and ϕ and grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
8
×
8
×
8
10
×
10
×
10
12
×
12
×
12
1
Legendre polynomials(θ,ϕ)=(0,0)
4.29
e
-
02
4.44
e
-
02
4.52
e
-
02
4
5.54
e
-
08
6.03
e
-
08
6.30
e
-
08
5
2.85
e
-
11
8.77
e
-
12
9.27
e
-
12
6
2.95
e
-
11
1.51
e
-
14
6.44
e
-
15
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
4.37
e
-
02
4.51
e
-
02
4.57
e
-
02
3
6.80
e
-
04
7.09
e
-
04
7.24
e
-
04
5
1.33
e
-
09
9.88
e
-
12
9.75
e
-
12
6
1.33
e
-
09
9.39
e
-
13
7.66
e
-
15
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
4.33
e
-
02
4.46
e
-
02
4.53
e
-
02
2
5.39
e
-
03
5.63
e
-
03
5.76
e
-
03
4
5.58
e
-
08
6.10
e
-
08
6.35
e
-
08
6
9.59
e
-
10
6.78
e
-
13
1.88
e
-
14
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
4.23
e
-
02
4.41
e
-
02
4.50
e
-
02
3
6.33
e
-
04
6.78
e
-
04
7.01
e
-
04
5
5.63
e
-
10
8.53
e
-
12
9.17
e
-
12
6
5.63
e
-
10
3.26
e
-
13
8.55
e
-
15Figure 11
Error components for the nonlinear 3D Poisson equation.


### 5.2.9. The Nonlinear Complex Generalized Zakharov System

The nonlinear complex Zakharov system has importance in plasma physics [3]. The system includes two coupled nonlinear PDEs which can be written as(51)i∂tψx,t+δ1∂ttψx,t-δ2ψx,twx,t+δ3ψx,t2ψx,t=0,∂ttwx,t-cs2∂xxwx,t-δ4∂xxψx,t2=0,x,t∈ax,bx×at,bt,subject to the initial and boundary conditions(52)ψax,t=ψ1t,ψbx,t=ψ2t,ψx,0=ψ0x,wax,t=w1t,wbx,t=w2x,wx,0=w3x,x,t∈ax,bx×at,bt.Several numerical methods have been proposed recently for approximating the solution of (51) and (52) such as the homotopy method [24], the finite difference method [25, 26], and the variational iteration method [27]. Also, Bao et al. [28] suggested some high-accurate numerical methods for solving numerically (51) and (52). Bao and Sun [29] applied a new technique based on time-splitting discretization for approximating the solution of a variant of (51) and (52).One can split (51) using the real and imaginary parts of ψ(x,t), u(x,t), and v(x,t), as(53)∂tux,t+δ1∂xxvx,t-δ2vx,twx,t+δ3u2x,t+v2x,tvx,t=0,-∂tvx,t+δ1∂xxux,t-δ2ux,twx,t+δ3u2x,t+v2x,tux,t=0,∂ttwx,t-cs2∂xxwx,t-2δ4ux,t∂xxux,t+∂xux,t2+vx,t∂xxvx,t+∂xvx,t2=0,with the initial and boundary conditions(54)uax,t=α1t,ubx,t=α2t,vax,t=α3t,vbx,t=α4t,wax,t=α5t,wbx,t=α6t,ux,at=β1x,vx,at=β2x,wx,at=β3x,wtx,at=β4x,x,t∈ax,bx×at,bt.The matrix form of the nonlinear system (53) is(55)∂tδ1∂xx0δ1∂xx-∂t000∂tt-cs2∂xxuvw+q1u,v,wq2u,v,wq3u,v,w=000,where (56)q1=-δ2vx,twx,t+δ3u2x,t+v2x,tvx,t,q2=-δ2ux,twx,t+δ3u2x,t+v2x,tux,t,q3=-2δ4ux,t∂xxux,t+∂xux,t2+vx,t∂xxvproblemx,t+∂xvx,t2,where the constants δ1, δ2, δ3, δ4, and cs are given and the functions αi(t), 1≤i≤6 and βj(x), 1≤j≤4 are known. We use J-GL-C methods for discretizing (55) subject to the initial and boundary conditions (54) to reduce (55) to a system of nonlinear algebraic equations.


### 5.2.10. Complex Generalized Zakharov Equation I

The first complex generalized Zakharov problem [3] is (57)i∂tψ-∂xxψ-ψw=0,∂ttw-∂xxw+∂xxψ2=0,with domain [-1,1]×[0,3.0]. The nonlinear problem (57) has the analytical solution (58)ψx,t=ux,t+ivx,t=3eix+ttanh⁡12x+2t,wx,t=1-tanh212x+2t.We define(59)E=maxU-Uanalytical∞,V-Vanalytical∞,W-Wanalytical∞to check the numerical accuracy in the computed solutions. Problem (57) is discretized by using Chebyshev collocation method of first kind. We took the initial guess 0 with initial and boundary conditions introduced and a grid 50×20. The achieved numerical accuracy is shown in Table 17. Our computed solutions are better than those in [3] because the maximum accuracy in [3] is of the order of 10-7. The solutions and error components are plotted in Figures 12, 13, and 14.Table 17
Maximum absolute error in the solution of the system of nonlinear equations for complex generalized Zakharov equation and varyingm.
m
E
1
8.56
e
-
1
2
2.72
e
0
3
4.16
e
0
4
8.92
e
0
5
1.12
e
1
6
3.40
e
1
7
2.04
e
1
8
2.83
e
1
9
4.05
e
0
10
2.34
e
0
11
7.77
e
-
1
12
3.25
e
-
1
13
1.14
e
-
1
14
6.20
e
-
2
15
2.79
e
-
2
16
1.85
e
-
2
17
8.11
e
-
3
18
4.44
e
-
3
19
1.79
e
-
3
20
8.86
e
-
4
21
3.91
e
-
4
22
1.66
e
-
4
23
5.98
e
-
5
24
1.97
e
-
5
25
5.32
e
-
6
26
1.20
e
-
6
27
4.12
e
-
7
28
1.36
e
-
7
29
6.95
e
-
8
30
2.95
e
-
8
31
9.29
e
-
9
32
3.42
e
-
9
33
1.43
e
-
9
34
5.57
e
-
10
35
2.62
e
-
10
36
1.26
e
-
10
37
8.65
e
-
11Figure 12
Solution and error components for complex generalized Zakharov equation (57).Figure 13
Solution and error components for complex generalized Zakharov equation (57).Figure 14
Solution and error components for complex generalized Zakharov equation (57).


### 5.2.11. Complex Generalized Zakharov Equation II

The second complex generalized Zakharov equation [3] is(60)i∂tψ+∂xxψ+2ψw-2ψ2=0,∂ttw-∂xxw+∂xxψ2=0,with the domain [-1,1]×[0,1]. That problem has the analytical solution (61)ψx,t=ux,t+ivx,t=32e-ix+3ttanh⁡x+2t,wx,t=-14tanh2x+2t.Problem (60) is discretized by using Chebyshev collocation method of first kind. We took the initial guess 0 with initial and boundary conditions introduced and a grid 21×28. The achieved numerical accuracy is shown in Table 18. Our computed solutions are better than those in [3] because the ones in [3] achieved a maximum accuracy of the order of 10-6. The solutions and error components are plotted in Figures 15, 16, and 17.Table 18
Maximum absolute error in the solution of the system of nonlinear equations for complex generalized Zakharov equation and varyingm.
m
E
1
9.2329
e
-
1
2
7.0559
e
-
1
3
5.726
e
-
1
4
6.3466
e
-
2
5
4.3852
e
-
3
6
2.3221
e
-
4
7
1.5501
e
-
5
8
4.2396
e
-
7
9
2.8906
e
-
8
10
7.7742
e
-
10
11
8.6019
e
-
11Figure 15
Solution and error components for complex generalized Zakharov equation (60).Figure 16
Solution and error components for complex generalized Zakharov equation (60).Figure 17
Solution and error components for complex generalized Zakharov equation (60).


### 5.2.12. Murray Equation

Murray equation [2] is (62)ut=uxx+λ1uux+λ2u-λ3u2,x,t∈-1,1×0,1.The initial and boundary condition for Murray equation can be computed from its analytical solution (63)ux,t=λ22λ31+tanh⁡λ28λ322λ1λ3x+λ12+4λ32t.The Chebyshev collocation method of first kind is used to discretize Murray equation over the grid 16×14. The achieved numerical accuracy is shown in Table 19. The computed accuracy in [2] is of the order of 10-8. The solution and error components are shown in Figure 18.Table 19
Maximum absolute error in the solution of the system of nonlinear equations for Murray equation and varyingm.
m
E
1
2.63
e
-
01
2
1.33
e
-
01
3
3.85
e
-
02
4
1.47
e
-
04
5
7.43
e
-
07
6
2.15
e
-
09
7
7.53
e
-
12
8
5.59
e
-
14Figure 18
Solution and error components for Murray equation.


### 5.2.13. Nonlinear Reaction and Diffusion Equation

Nonlinear reaction and diffusion equation is(64)ut=uxx+u1+u-u2,x,t∈-1,1×0,1.The analytical solution is (65)ux,t=12+52tanh⁡104x+54t.The Chebyshev collocation method of first kind is used to discretize nonlinear reaction and diffusion equation over the grid 22×14. The initial guess is 0 with initial and boundary conditions introduced. The achieved numerical accuracy is shown in Table 20. The computed accuracy in [2] is of the order of 10-8. The solution and error components are shown in Figure 19.Table 20
Maximum absolute error in the solution of the system of nonlinear equations for nonlinear reaction and diffusion equation and varyingm.
m
E
1
7.28
e
-
02
2
2.12
e
-
02
3
4.99
e
-
03
4
1.09
e
-
06
5
9.20
e
-
10
6
6.53
e
-
13
7
2.24
e
-
14Figure 19
Solution and error components for nonlinear reaction and diffusion equation.


### 5.2.14. Reaction-Diffusion Equation with Fischer-Kolmogorov Reaction Term and Density Dependent Diffusion

Reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion is(66)ut=u2xx+u1-u,x,t∈-1,1×0,1.The analytical solution is (67)ux,t=21+tanh⁡x+t/4.The Chebyshev collocation method of first kind is used to discretize reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion over the grid 30×30. The initial guess is 0 with initial and boundary conditions introduced. The achieved numerical accuracy after 1 and 2 application of NMIM is shown in Table 21. The computed accuracy in [2] is only of the order of 10-9. The solution and error components are shown in Figure 20.Table 21
Maximum absolute error in the solution of the system of nonlinear equations for reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion and varyingm.
m (1 Iter.)
E
1
1.49
e
+
00
2
1.12
e
+
00
3
1.07
e
+
00
4
4.92
e
-
01
5
3.73
e
-
01
m (2 Iter.)
E
1
1.51
e
-
02
2
7.84
e
-
03
3
6.72
e
-
04
4
5.17
e
-
07
5
1.04
e
-
09
6
2.35
e
-
12
7
5.04
e
-
14Figure 20
Solution and error components for reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion.


### 5.2.15. Two-Dimensional Nonlinear Reaction-Diffusion Equation

The two-dimensional nonlinear reaction-diffusion equation is(68)ut=12∇2u+u2u-1,x,y,t∈-1,1×-1,1×0,20.Its analytical solution is (69)ux,y,t=11+e1/2x+y-1/2t.The Chebyshev collocation method of first kind is used to discretize the two-dimensional nonlinear reaction-diffusion equation over grid 30×30. The initial guess is 0 with initial and boundary conditions introduced. The achieved numerical accuracy is shown in Table 22. The computed accuracy in [2] is only of the order of 10-9. The error components are shown in Figure 21.Table 22
Maximum absolute error in the solution of the system of nonlinear equations for two-dimensional nonlinear reaction-diffusion equation and varyingm.
m
E
1
7.82
e
-
02
2
1.19
e
-
02
3
3.63
e
-
03
4
4.04
e
-
05
5
9.89
e
-
07
6
2.65
e
-
08
7
7.48
e
-
10Figure 21
Error components for the two-dimensional nonlinear reaction-diffusion equation.


### 5.2.16. Ostrovsky-Hunter Equation

Ostrovsky-Hunter equation is(70)ut-βuxxx+αuuxx=γu,where α is a nonlinear coefficient and β∈R, γ>0 are the dispersion coefficients. We assume a periodic boundary condition for this problem; that is, (71)u0,t=u2π,t,ux0,t=ux2π,t,uxx0,t=uxx2π,t,uxxx0,t=uxxx2π,t.The Chebyshev collocation method of first kind is used to discretize Ostrovsky-Hunter equation over a grid30×100. The initial guess is 0 with initial and boundary conditions introduced. The residue is shown in Table 23. Our solution for the grid {0,6π/11,π,16π/11}×{0.1,0.2,0.3,0.4} is given in Table 24. We notice that our results do not agree with those presented in [4]. Our solution is plotted in Figure 22.Table 23
Residue of the system of nonlinear equations for Ostrovsky-Hunter equation and varyingm.
m
F
(
x
)
∞
1
2
2
2.5963
3
2.797
4
0.72034
5
0.15165
6
0.025873
7
0.002996
8
0.00020443
9
8.7968
e
-
06
10
1.4609
e
-
07
11
1.3542
e
-
08
12
8.8853
e
-
10
13
8.5753
e
-
11Table 24
Solution for Ostrovsky-Hunter equation.
x
t
0.1
0.2
0.3
0.4
0
0.094587
0.17222
0.20233
0.21608
6
π
/
11
0.96556
0.94696
0.8702
0.73399
π
−0.10438
−0.22921
−0.41929
−0.67358
16
π
/
11
−0.98195
−0.96996
−0.97341
−0.97748Figure 22
Solution for Ostrovsky-Hunter equation.


## 5.2.1. The Lane-Emden Equation

The Lane-Emden equation is(21)x′′t+2tx′t+xtp=0,x′0=0,x0=1.J-GL-C discretization methods translate the Lane-Emden equation into the following system of nonlinear equations:(22)Fx=Sxxx+diag2tSxx+xp=0,F′x=Sxx+diag2tSx+pdiagxp-1,where diag(·) denotes a diagonal matrix with values in the diagonal equal to the elements of its vector argument, Sx=2/bQx, Sxx=2/bQx2, and [0,b] is the domain of the problem. We have solved the Lane-Emden equation for different values of p ranging from two to five and b=3. Figure 1 shows the initial guess and the computed numerical solutions for θ=1/2, ϕ=-1/2, and 50 grid points. Table 5 gives the norm of the residue in the solution of the nonlinear systems by a single iteration of NMIM with varying m.Table 5
Norm of the residue of the system of nonlinear equation resulting from Lane-Emden equation with varyingm.
F
(
x
k
)
∞
m
p
2
3
4
5
10
2.69
e
-
17
3.70
e
-
14
1.04
e
-
12
1.03
e
-
11
20
1.36
e
-
46
7.06
e
-
41
9.33
e
-
39
1.22
e
-
37
30
9.33
e
-
76
1.68
e
-
67
7.56
e
-
65
1.97
-
63
40
6.16
e
-
105
4.53
e
-
94
6.19
e
-
91
2.26
e
-
89
50
4.08
e
-
134
1.21
e
-
120
5.11
e
-
117
2.18
e
-
115
60
2.71
e
-
163
3.22
e
-
147
4.24
e
-
143
2.00
e
-
141
70
1.79
e
-
192
8.555
e
-
174
3.52
e
-
169
1.83
e
-
167
80
1.19
e
-
221
2.27
e
-
200
2.92
e
-
195
1.69
e
-
193
90
7.89
e
-
251
6.03
e
-
227
2.43
e
-
221
1.58
e
-
219
100
5.23
e
-
280
1.60
e
-
253
2.02
e
-
247
1.47
e
-
245Figure 1
Initial guess and solutions of the Lane-Emden equation for different values ofp.


## 5.2.2. Bratu Problem

Bratu problem is the boundary value problem:(23)x′′t+αext=0,x0=0,x1=0,where α is a parameter. The discretization is carried out using J-GL-C methods. The resulting system of nonlinear equations is (24)Fx=Sxxx+αex=0,F′x=Sxx+αdiagex,where Sxx=2Qx2. We considered several values of α, values θ=-1/2 and ϕ=-1/2 for the parameters of J-GL-C, 50 grid points, and solved the resulting system of nonlinear equations using a single application of NMIM with varying m. The initial guess and computed numerical solutions are given in Figure 2. Table 6 gives the norm of the residue of the solution of the system of nonlinear equations for varying m.Table 6
Norm of the residue of the solution of the system of nonlinear equations for Bratu problem and varyingm.
F
(
x
k
)
∞
m
α
1
2
3
10
7.68
e
-
41
4.10
e
-
31
1.97
e
-
28
20
2.20
e
-
93
1.01
e
-
70
9.65
e
-
64
30
7.02
e
-
146
2.50
e
-
110
4.72
e
-
99
40
2.23
e
-
198
6.17
e
-
150
2.30
e
-
134
50
7.13
e
-
251
1.52
e
-
189
1.12
e
-
169
60
2.27
e
-
303
3.75
e
-
229
5.52
e
-
205
70
7.30
e
-
356
9.27
e
-
269
2.69
e
-
240
80
2.30
e
-
408
2.28
e
-
308
1.32
e
-
275
90
7.40
e
-
461
5.70
e
-
348
6.50
e
-
311
100
2.30
e
-
513
1.40
e
-
387
3.20
e
-
346Figure 2
Initial guess and solutions of Bratu problem for different values ofα.


## 5.2.3. Frank-Kamenetskii Problem

Frank-Kamenetskii problem is(25)x′′t+1tx′t+αext=0,x′0=0,x1=0,where α is a parameter. After discretization by J-GL-C methods we get the system of nonlinear equations:(26)Fx=Sxxx+diag1tSxx+αex=0,F′x=Sxx+diag1tSx+αdiagex.We considered several values of α, values θ=-1/2 and ϕ=-1/2 for the parameters of J-GL-C, 50 grid points, and solved the resulting system of nonlinear equations by using a single application of NMIM with varying m. Figure 3 plots the initial guess and computed solutions and Table 7 gives the norm of the residue of the solution of the system of nonlinear equations.Table 7
Norm of the residue of the solution of the system of nonlinear equations for Fran-Kamenetskii problem and varyingm.
F
(
x
k
)
∞
m
α
1
1.1
1.2
1.3
10
7.09
e
-
41
1.89
e
-
39
2.63
e
-
38
1.80
e
-
37
20
2.46
e
-
82
2.36
e
-
79
6.21
e
-
77
4.05
e
-
75
30
8.58
e
-
124
2.95
e
-
119
1.46
e
-
115
9.09
e
-
113
40
2.98
e
-
165
3.69
e
-
159
3.45
e
-
154
2.03
e
-
150
50
1.03
e
-
206
4.61
e
-
199
8.15
e
-
193
4.57
e
-
188
60
3.61
e
-
248
5.76
e
-
239
1.92
e
-
231
1.02
e
-
225
70
1.25
e
-
289
7.19
e
-
279
4.53
e
-
270
2.29
e
-
263
80
4.40
e
-
331
9.00
e
-
319
1.10
e
-
308
5.15
e
-
301
90
1.50
e
-
372
1.10
e
-
358
2.50
e
-
347
1.20
e
-
338
100
5.30
e
-
414
1.40
e
-
398
6.00
e
-
386
2.60
e
-
376Figure 3
Initial guess and solutions of Frank-Kamenetskii problem with different values ofα.


## 5.2.4. 1 + 1 Nonlinear Schrödinger Equation

1 + 1 nonlinear Schrödinger equation is(27)i∂tψx,t+∂xxψx,t+2γψx,t2ψx,t-2δRx,tψx,t=0,x,t∈Dx×Dt,where(28)Dx=ax,bx,Dt=0,tf,with the initial and boundary conditions (29)ψax,t=η1t,ψbx,t=η2t,ψx,0=ζx.The function ψ(x,t) is a complex function and can be written as ψ(x,t)=ψ1(x,t)+iψ2(x,t). The Schrödinger equation can be rewritten in terms of the real functions ψ1(x,t) and ψ2(x,t) as(30)∂t∂xx∂xx-∂tψ1ψ2-2δRx,t0110ψ1ψ2+2γψ12+ψ220110ψ1ψ2=0,with initial and boundary conditions (31)ψ1ax,t+iψ2ax,t=η11t+iη12,ψ1bx,t+iψ2bx,t=η21t+iη22t,ψ1x,0+iψ2x,0=ζ1x+iζ2x.Using J-GL-C methods, we can discretize (30) as(32)Fψ=St⊗IxIt⊗SxxIt⊗Sxx-St⊗Ixψ1ψ2-2δOdiagRdiagROψ1ψ2+2γOdiagψ12+ψ22diagψ12+ψ22Oψ1ψ2=0,F′ψ=St⊗IxIt⊗SxxIt⊗Sxx-St⊗Ix-2δOdiagRdiagRO+2γ2diagψ1⊙ψ23diagψ223diagψ122diagψ1⊙ψ2,where ψj=[ψj1,1,ψj1,2,…,ψj1,nt,…,ψjnx,1,ψjnx,2, …,ψj(nx,nt)]T, j=1,2, R=[R1,1,R1,2,…,R1,nt, …,Rnx,1,Rnx,2,…,R(nx,nt)]T, ⊗ is the Kronecker product, ⊙ is the element-wise multiplication, O is a matrix of all zeros, I is the identity matrix, St=2/tfQt, Sxx=2/bx-axQx2, and Q is the J-GL-C operational matrix. The system of nonlinear equations (32) can be rewritten as (33)Fψ=Aψ+gψ-p,where A=St⊗IxIt⊗SxxIt⊗Sxx-St⊗Ix-2δOdiag(R)diag(R)O, g(ψ)=2γ2diag(ψ1⊙ψ2)3diag(ψ22)3diag(ψ12)2diag(ψ1⊙ψ2), and p is the vector incorporating the initial-boundary conditions. Four nonlinear Schrödinger equations are listed in Table 8 with their corresponding analytical solutions. We chose a domain (x,t)∈[-1,1]×[0,1] and solved the equations using different pair of values for the parameters θ and ϕ of J-GL-C and different grids for the domain. As initial guess, we used ψ1=0, ψ2=0 with boundary conditions set for problem 3. For problems 1, 2, and 4 NMIM does not converge with that initial guess and we made the initial guess smoother. We smoothed the initial guess by applying the iteration(34)ψ=-A-12gψ-ponce for problem 1 and twice for problems 2 and 4. Tables 9, 10, 11, and 12 give the norm of the error in the solution when the nonlinear system of equations is solved using a single application of NMIM with varying m. Figures 4–7 display the solutions and the error in the components of ψ1 and ψ2 for θ=-4/2, ϕ=-1/2, and a grid 15×30.Table 8
Nonlinear Schrödinger equations.
Schrödinger equation
Analytical Solution
i∂tψ+∂xxψ+2|ψ|2ψ=0
ψ
(
x
,
t
)
=
e
i
(
x
+
t
)
i
∂
t
ψ
+
∂
x
x
ψ
-
2
|
ψ
|
2
ψ
=
0
ψ
(
x
,
t
)
=
e
i
(
x
-
3
t
)
i∂tψ+12∂xxψ-|ψ|2ψ-cos⁡(x)2ψ=0
ψ
(
x
,
t
)
=
e
-
3
i
t
/
2
sin
⁡
(
x
)
i
∂
t
ψ
+
∂
x
x
ψ
+
2
|
ψ
|
2
ψ
-
2
δ
ψ
=
0
ψ
(
x
,
t
)
=
e
i
(
x
+
(
1
-
2
δ
)
t
)Table 9
Norm of the error in the solution for equation (17).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
5.49
e
-
01
6.21
e
-
01
6.41
e
-
01
4
5.36
e
-
02
2.35
e
-
01
2.44
e
-
01
8
2.85
e
-
04
4.55
e
-
04
4.97
e
-
04
12
2.85
e
-
04
5.81
e
-
08
6.38
e
-
08
16
2.85
e
-
04
5.69
e
-
10
4.38
e
-
12
20
2.85
e
-
04
6.13
e
-
11
5.95
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
6.86
e
-
01
Divergent
6.40
e
-
01
4
1.96
e
-
01
2.48
e
-
01
8
1.64
e
-
03
4.94
e
-
04
12
1.59
e
-
03
6.16
e
-
08
16
1.59
e
-
03
6.02
e
-
12
20
1.59
e
-
03
2.16
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
5.24
e
-
01
6.25
e
-
01
6.41
e
-
01
4
5.90
e
-
02
2.39
e
-
01
2.46
e
-
01
8
1.14
e
-
03
5.05
e
-
04
5.40
e
-
04
12
1.14
e
-
03
7.73
e
-
08
8.00
e
-
08
16
1.14
e
-
03
1.90
e
-
09
5.28
e
-
12
20
1.14
e
-
03
1.90
e
-
09
4.99
e
-
13
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
5.95
e
-
01
6.20
e
-
01
6.41
e
-
01
4
3.73
e
-
02
2.37
e
-
01
2.45
e
-
01
8
1.52
e
-
03
4.24
e
-
04
4.81
e
-
04
12
1.52
e
-
03
5.22
e
-
08
5.81
e
-
08
16
1.52
e
-
03
4.90
e
-
10
1.71
e
-
12
20
1.52
e
-
03
4.90
e
-
10
9.31
e
-
14Table 10
Norm of the error in the solution for equation (30).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
1.28
e
+
00
1.34
e
+
00
1.34
e
+
00
4
4.27
e
-
01
9.65
e
-
01
9.72
e
-
01
8
1.25
e
-
03
2.15
e
-
03
2.18
e
-
03
12
1.46
e
-
04
9.17
e
-
07
1.05
e
-
06
16
1.46
e
-
04
3.32
e
-
10
4.61
e
-
10
20
1.46
e
-
04
2.90
e
-
11
1.85
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
1.27
e
+
00
1.35
e
+
00
1.34
e
+
00
4
4.34
e
-
01
9.46
e
-
01
9.72
e
-
01
8
1.63
e
-
03
2.18
e
-
03
2.18
e
-
03
12
1.52
e
-
03
3.89
e
-
07
8.77
e
-
07
16
1.52
e
-
03
5.28
e
-
09
1.95
e
-
10
20
1.52
e
-
03
3.73
e
-
09
1.34
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
1.26
e
+
00
1.34
e
+
00
1.34
e
+
00
4
3.98
e
-
01
9.24
e
-
01
9.40
e
-
01
8
2.63
e
-
03
1.35
e
-
03
1.63
e
-
03
12
1.00
e
-
03
2.63
e
-
06
2.85
e
-
06
16
1.00
e
-
03
2.35
e
-
08
1.65
e
-
08
20
1.00
e
-
03
7.88
e
-
10
6.29
e
-
11
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
1.33
e
+
00
1.34
e
+
00
1.34
e
+
00
4
3.90
e
-
01
9.80
e
-
01
9.80
e
-
01
8
1.03
e
-
03
2.15
e
-
03
2.19
e
-
03
12
1.03
e
-
03
6.56
e
-
07
7.42
e
-
07
16
1.03
e
-
03
3.08
e
-
10
1.57
e
-
10
20
1.03
e
-
03
2.95
e
-
10
5.77
e
-
14Table 11
Norm of the error in the solution for equation (32).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomial(θ,ϕ)=(0,0)
9.98
e
-
02
8.73
e
-
02
9.66
e
-
02
4
7.62
e
-
04
3.16
e
-
04
3.25
e
-
04
8
4.27
e
-
05
3.64
e
-
11
3.12
e
-
11
12
4.27
e
-
05
3.16
e
-
11
1.89
e
-
13
16
4.27
e
-
05
3.16
e
-
11
1.88
e
-
13
20
4.27
e
-
05
3.16
e
-
11
1.89
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
9.71
e
-
02
9.33
e
-
02
9.67
e
-
02
4
1.18
e
-
03
3.37
e
-
04
3.36
e
-
04
8
7.68
e
-
04
3.23
e
-
08
2.93
e
-
11
12
7.68
e
-
04
1.89
e
-
10
1.07
e
-
13
16
7.68
e
-
04
1.87
e
-
10
1.06
e
-
13
20
7.68
e
-
04
1.87
e
-
10
1.07
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
9.47
e
-
02
9.17
e
-
02
9.65
e
-
02
4
9.17
e
-
04
3.14
e
-
04
3.38
e
-
04
8
6.70
e
-
04
9.13
e
-
11
8.39
e
-
11
12
6.70
e
-
04
8.92
e
-
11
1.19
e
-
13
16
6.70
e
-
04
8.92
e
-
11
1.19
e
-
13
20
6.70
e
-
04
8.92
e
-
11
1.19
e
-
13
1
Chebyshev polynomial of first kind(θ,ϕ)=(-1/2,-1/2)
1.01
e
-
01
9.01
e
-
02
9.61
e
-
02
4
1.31
e
-
03
2.99
e
-
04
3.26
e
-
04
8
8.53
e
-
04
4.90
e
-
11
2.27
e
-
11
12
8.53
e
-
04
5.59
e
-
11
3.28
e
-
14
16
8.53
e
-
04
5.59
e
-
11
3.24
e
-
14
20
8.53
e
-
04
5.59
e
-
11
3.20
e
-
14Table 12
Norm of the error in the solution for equation (35).
ψ
-
ψ
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
30
10
×
30
15
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
6.44
e
-
01
7.84
e
-
01
8.01
e
-
01
4
6.52
e
-
02
4.18
e
-
02
4.18
e
-
02
8
2.77
e
-
04
1.36
e
-
06
1.33
e
-
06
12
2.75
e
-
04
1.70
e
-
10
1.81
e
-
10
16
2.75
e
-
04
4.14
e
-
11
4.87
e
-
13
20
2.75
e
-
04
4.14
e
-
11
4.67
e
-
13
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
6.51
e
-
01
Divergent
8.08
e
-
01
4
4.52
e
-
02
4.15
e
-
02
8
1.50
e
-
03
1.25
e
-
06
12
1.51
e
-
03
1.17
e
-
10
16
1.51
e
-
03
1.94
e
-
13
20
1.51
e
-
03
1.99
e
-
13
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
7.48
e
-
01
7.85
e
-
01
8.04
e
-
01
4
5.30
e
-
02
4.24
e
-
02
4.22
e
-
02
8
1.02
e
-
03
2.20
e
-
06
2.61
e
-
06
12
1.02
e
-
03
2.67
e
-
09
1.00
e
-
09
16
1.02
e
-
03
2.25
e
-
09
4.86
e
-
13
20
1.02
e
-
03
2.25
e
-
09
1.96
e
-
13
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
7.03
e
-
01
7.88
e
-
01
8.02
e
-
01
4
3.66
e
-
02
4.20
e
-
02
4.17
e
-
02
8
1.46
e
-
03
1.27
e
-
06
1.18
e
-
06
12
1.46
e
-
03
5.96
e
-
10
2.39
e
-
11
16
1.46
e
-
03
5.84
e
-
10
8.56
e
-
14
20
1.46
e
-
03
5.84
e
-
10
8.49
e
-
14Figure 4
Solution and error components for (17).Figure 5
Solution and error components for (30).Figure 6
Solution and error components for (32).Figure 7
Solution and error components for (35).


## 5.2.5. Klein Gordon Equation

Klein Gordon equation is(35)uttx,t-c2uxxx,t+ku-γux,t3=0,x,t∈Dx×Dt,where(36)Dx=ax,bx,Dt=0,tf,with the initial and boundary conditions (37)uax,t=ϱ1t,ubx,t=ϱ2t,ux,0=ξ1x,utx,0=ξ2x.By discretization using J-GL-C methods we obtain (38)Fu=Sttu-c2Sxxu+ku-γu3=0,F′u=Stt-c2Sxx+kI-3γdiagu2,where u=[u(1,1),u(1,2),…,u(1,nt),…,u(nx,1),u(nx,2),…, u(nx,nt)]T, Stt=2/tfQt2, and Sxx=2/bx-axQx2. The analytical solution of (35) is u(x,t)=δsech⁡κ(x-vt), where κ=k/c2-v2 and δ=2k/γ. We set the parameters c=1, γ=1, v=0.5, and k=0.5, picked up the domain (x,t)∈[-10,10]×[0,1], and solved the nonlinear system of equations by a single application of NMIM with initial guess u=0, with the boundary conditions introduced. We used several pairs of the θ, ϕ parameters of J-GL-C. According to Table 13, the Chebyshev collocation method of first kind exhibits best accuracy when the grid is more refined. The solution and the error components for θ=-1/2, ϕ=-1/2, nx=120, and nt=20 are plotted in Figure 8.Table 13
Norm of the error in the Klein Gordon equation as a function ofm for several values for θ and ϕ and grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
30
×
30
60
×
30
120
×
30
1
Legendre polynomials(θ,ϕ)=(0,0)
3.80
e
-
01
4.10
e
-
01
4.09
e
-
01
2
7.37
e
-
02
8.39
e
-
02
8.34
e
-
02
3
1.18
e
-
02
8.73
e
-
03
8.65
e
-
03
4
8.53
e
-
03
1.42
e
-
05
1.40
e
-
06
5
8.53
e
-
03
1.42
e
-
05
4.37
e
-
09
6
8.53
e
-
03
1.42
e
-
05
4.56
e
-
09
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
3.98
e
-
01
4.10
e
-
01
4.07
e
-
01
2
7.64
e
-
02
8.40
e
-
02
8.28
e
-
02
3
1.01
e
-
02
8.76
e
-
03
8.55
e
-
03
4
7.67
e
-
03
3.55
e
-
05
1.38
e
-
06
5
7.67
e
-
03
3.55
e
-
05
1.60
e
-
09
6
7.67
e
-
03
3.55
e
-
05
1.45
e
-
09
1
Chebyshev polynomials of second Kind(θ,ϕ)=(1/2,1/2)
3.81
e
-
01
4.10
e
-
01
4.09
e
-
01
2
7.41
e
-
02
8.39
e
-
02
8.34
e
-
02
3
1.12
e
-
02
8.73
e
-
03
8.65
e
-
03
4
7.85
e
-
03
1.29
e
-
05
1.41
e
-
06
5
7.85
e
-
03
1.29
e
-
05
2.38
e
-
09
6
7.85
e
-
03
1.29
e
-
05
2.23
e
-
09
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
3.78
e
-
01
4.10
e
-
01
4.09
e
-
01
2
7.33
e
-
02
8.38
e
-
02
8.34
e
-
02
3
1.25
e
-
02
8.73
e
-
03
8.65
e
-
03
4
9.28
e
-
03
1.58
e
-
05
1.41
e
-
06
5
9.28
e
-
03
1.58
e
-
05
3.63
e
-
10
6
9.28
e
-
03
1.58
e
-
05
3.40
e
-
10Figure 8
Solution and error components for nonlinear Klein Gordon equation.


## 5.2.6. 2D Nonlinear Wave Equation

The 2D nonlinear wave equation is(39)uttx,y,t-c2uxxx,y,t+uyyx,y,t+fu=pt,x,y,x,y,t∈Dx×Dy×Dt,where(40)Dx=ax,bx,Dy=ay,by,Dt=0,tf,with the initial and boundary conditions (41)uax,y,t=ϱ11y,t,ubx,y,t=ϱ12y,t,ux,ay,t=ϱ21x,t,ux,by,t=ϱ22x,t,ux,y,0=ζx,y.By discretization using J-GL-C methods we obtain (42)Fu=Stt⊗Ix⊗Iyu-c2It⊗Sxx⊗Iy+It⊗Ix⊗Syyu+fu-p=0,F′u=Stt⊗Ix⊗Iy-c2It⊗Sxx⊗Iy+It⊗Ix⊗Syy+diagf′u,where u=[u1,1,1,u1,1,2,…,u1,1,nt,…,unx,ny,1,unx,ny,2, …,u(nx,ny,nz)]T, Stt=2/tfQt2, Sxx=2/bx-axQx2, and Syy=2/by-ayQy2. The solution of (39) with c=1, f(u)=u3 and appropriate initial and boundary conditions is e-tsin⁡(x+y). We took a domain (x,y,t)∈[-1,1]×[-1,1]×[0,2] and solved the system of nonlinear equations using a single application of NMIM with initial guess u=0 with the initial and boundary conditions introduced. Table 14 gives the norm of the error as a function of m for several pairs of values for θ and ϕ and several grids. The collocation method with 20×20×20 which worked best was the Chebyshev collocation method of first kind. Figure 9 gives the error components for that collocation method and that grid.Table 14
Norm of the error for the 2D nonlinear wave equation as a function ofm for several pair of values for θ and ϕ and several grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
5
×
20
10
×
10
×
20
20
×
20
×
20
1
Legendre polynomials(θ,ϕ)=(0,0)
4.21
e
-
02
Divergent
Divergent
2
3.85
e
-
03
3
2.03
e
-
04
9
1.09
e
-
04
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
4.45
e
-
02
5.12
e
-
02
5.14
e
-
02
2
3.50
e
-
03
3.74
e
-
03
1.53
e
-
02
5
7.34e-04
5.74
e
-
08
2.16
e
-
04
7
7.34
e
-
04
1.40
e
-
09
3.76
e
-
06
9
7.34
e
-
04
1.40
e
-
09
6.56
e
-
08
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
4.40
e
-
02
5.13
e
-
02
Divergent
2
3.43
e
-
03
3.81
e
-
03
3
5.98
e
-
04
1.62
e
-
04
4
4.39
e
-
04
8.01
e
-
09
5
4.39
e
-
04
7.64
e
-
10
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
3.88
e
-
02
5.02
e
-
02
5.13
e
-
02
2
3.91
e
-
03
3.85
e
-
03
3.80
e
-
03
3
5.94
e
-
04
1.54
e
-
04
1.79
e
-
03
4
6.60
e
-
04
4.54
e
-
09
4.28
e
-
05
5
6.60
e
-
04
6.78
e
-
10
7.94
e
-
07
9
6.60
e
-
04
6.78
e
-
10
5.39
e
-
11Figure 9
Error components for the nonlinear 2D wave equation.


## 5.2.7. 3D Nonlinear Wave Equation

The 3D nonlinear wave equation is(43)uttx,y,z,t-c2uxxx,y,z,t+uyyx,y,z,t+uzzx,y,z,t+fu=px,y,z,x,y,z,t∈Dx×Dy×Dz×Dt,where(44)Dx=ax,bx,Dy=ay,by,Dz=az,bz,Dt=0,tf,with the initial and boundary conditions (45)uax,y,z,t=ϱ11y,z,t,ubx,y,z,t=ϱ12y,z,t,ux,ay,z,t=ϱ21x,z,t,ux,by,z,t=ϱ22x,z,t,ux,y,az,t=ϱ31x,y,t,ux,y,bz,t=ϱ32x,y,t,ux,y,z,0=ζx,y,z.By discretization using J-GL-C methods we obtain (46)Fu=Stt⊗Ix⊗Iy⊗Izu-c2It⊗Sxx⊗Iy⊗Iz+It⊗Ix⊗Syy⊗Iz+It⊗Ix⊗Iy⊗Szzu+fu-p=0,F′u=Stt⊗Ix⊗Iy⊗Iz-c2It⊗Sxx⊗Iy⊗Iz+It⊗Ix⊗Syy⊗Iz+It⊗Ix⊗Iy⊗Szz+diagf′u,where u=[u(1,1,1,1),u(1,1,1,2),…,u(1,1,1,nt),…,u(nx,ny,nz,1), u(nx,ny,nz,2),…,u(nx,ny,nz,nt)]T, Stt=2/tfQt2, Sxx=2/bx-axQx2, Syy=2/by-ayQy2, and Szz=2/bz-azQz2. The solution of (43) with c=0.8, f(u)=u2 and appropriate initial and boundary conditions is e-tsin⁡(x+y+z). We took a domain (x,y,z,t)∈[0,1]×[0,1]×[0,1]×[0,1] and solved the system of nonlinear equations using a single application of NMIM with initial guess u=0 and with initial and boundary conditions introduced. Table 15 gives the norm of the error as a function of m for several pairs of values for θ,  ϕ and several grids. The collocation method which worked best for the grid 7×7×7×12 was the Legendre collocation method. Figure 10 gives the error components for that collocation method and that grid.Table 15
Norm of the error for the 3D nonlinear wave equation as a function ofm for several values for θ and ϕ and grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
5
×
5
×
5
×
12
6
×
6
×
6
×
12
7
×
7
×
7
×
12
1
Legendre polynomials(θ,ϕ)=(0,0)
9.95
e
-
01
1.00
e
+
00
1.14
e
+
00
3
1.14
e
-
02
8.53
e
-
03
1.09
e
-
02
5
1.25
e
-
01
1.12
e
-
01
1.73
e
-
01
7
5.53
e
-
06
4.12
e
-
07
1.73
e
-
06
8
5.54
e
-
06
7.54
e
-
08
4.24
e
-
09
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
1.00
e
+
00
1.00
e
+
00
divergent
2
1.00
e
-
01
8.68
e
-
02
4
1.31
e
+
00
1.29
e
+
00
7
2.24
e
-
05
1.99
e
-
06
8
2.24
e
-
05
2.11
e
-
06
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
9.98
e
-
01
9.99
e
-
01
1.32
e
+
00
2
1.02
e
-
01
8.06
e
-
02
1.08
e
-
01
5
1.22
e
-
01
1.18
e
-
01
2.46
e
-
01
7
1.81
e
-
05
1.84
e
-
06
3.72
e
-
05
8
1.81
e
-
05
1.76
e
-
06
4.87
e
-
08
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
9.95
e
-
01
divergent
1.09
e
+
00
2
9.98
e
-
02
1.04
e
-
01
5
1.31
e
-
01
1.51
e
-
01
7
1.92
e
-
05
4.37
e
-
07
8
1.92
e
-
05
3.27
e
-
08Figure 10
Error components for the nonlinear 3D wave equation.


## 5.2.8. 3D Nonlinear Poisson Equation

The 3D nonlinear Poisson equation is(47)uxxx,y,z+uyyx,y,z+uzzx,y,z+fu=px,y,z,x,y,z∈Dx×Dy×Dz,where(48)Dx=ax,bx,Dy=ay,by,Dz=az,bz,with the initial and boundary conditions (49)uax,y,z=ϱ11y,z,ubx,y,z=ϱ12y,z,ux,ay,z=ϱ21x,z,ux,by,z=ϱ22x,z,ux,y,az=ϱ31x,y,ux,y,bz=ϱ32x,y.By discretization using J-GL-C methods we obtain (50)Fu=Sxx⊗Iy⊗Iz+Ix⊗Syy⊗Iz+Ix⊗Iy⊗Szzu+fu-p=0,F′u=Sxx⊗Iy⊗Iz+Ix⊗Syy⊗Iz+Ix⊗Iy⊗Szz+diagf′u,where u=[u1,1,1,u1,1,2,…,u1,1,nz,…,unx,ny,1,unx,ny,2, …,u(nx,ny,nz)]T, Sxx=2/bx-axQx2, Syy=2/by-ayQy2, and Szz=2/bz-azQz2. The solution of (47) with f(u)=u4 with appropriate initial and boundary conditions is sin⁡(x+y+z). We took a domain (x,y,z)∈[0,1]×[0,1]×[0,1] and solved the system of nonlinear equations using a single application of NMIM with initial guess u=0 and with initial and boundary conditions introduced. Table 16 gives the norm of the error as a function of m for several pairs of values for θ and ϕ and several grids. The Legendre collocation method offers best accuracy over different grids. Figure 11 gives the error components for that collocation method and that grid.Table 16
Norm of the error for the 3D Nonlinear Poisson equation as a function ofm for several values for θ and ϕ and grids.
u
-
u
a
n
a
l
y
t
i
c
a
l
∞
m
Grid size
8
×
8
×
8
10
×
10
×
10
12
×
12
×
12
1
Legendre polynomials(θ,ϕ)=(0,0)
4.29
e
-
02
4.44
e
-
02
4.52
e
-
02
4
5.54
e
-
08
6.03
e
-
08
6.30
e
-
08
5
2.85
e
-
11
8.77
e
-
12
9.27
e
-
12
6
2.95
e
-
11
1.51
e
-
14
6.44
e
-
15
1
(
θ
,
ϕ
)
=
(
1
/
2
,
0
)
4.37
e
-
02
4.51
e
-
02
4.57
e
-
02
3
6.80
e
-
04
7.09
e
-
04
7.24
e
-
04
5
1.33
e
-
09
9.88
e
-
12
9.75
e
-
12
6
1.33
e
-
09
9.39
e
-
13
7.66
e
-
15
1
Chebyshev polynomials of second kind(θ,ϕ)=(1/2,1/2)
4.33
e
-
02
4.46
e
-
02
4.53
e
-
02
2
5.39
e
-
03
5.63
e
-
03
5.76
e
-
03
4
5.58
e
-
08
6.10
e
-
08
6.35
e
-
08
6
9.59
e
-
10
6.78
e
-
13
1.88
e
-
14
1
Chebyshev polynomials of first kind(θ,ϕ)=(-1/2,-1/2)
4.23
e
-
02
4.41
e
-
02
4.50
e
-
02
3
6.33
e
-
04
6.78
e
-
04
7.01
e
-
04
5
5.63
e
-
10
8.53
e
-
12
9.17
e
-
12
6
5.63
e
-
10
3.26
e
-
13
8.55
e
-
15Figure 11
Error components for the nonlinear 3D Poisson equation.


## 5.2.9. The Nonlinear Complex Generalized Zakharov System

The nonlinear complex Zakharov system has importance in plasma physics [3]. The system includes two coupled nonlinear PDEs which can be written as(51)i∂tψx,t+δ1∂ttψx,t-δ2ψx,twx,t+δ3ψx,t2ψx,t=0,∂ttwx,t-cs2∂xxwx,t-δ4∂xxψx,t2=0,x,t∈ax,bx×at,bt,subject to the initial and boundary conditions(52)ψax,t=ψ1t,ψbx,t=ψ2t,ψx,0=ψ0x,wax,t=w1t,wbx,t=w2x,wx,0=w3x,x,t∈ax,bx×at,bt.Several numerical methods have been proposed recently for approximating the solution of (51) and (52) such as the homotopy method [24], the finite difference method [25, 26], and the variational iteration method [27]. Also, Bao et al. [28] suggested some high-accurate numerical methods for solving numerically (51) and (52). Bao and Sun [29] applied a new technique based on time-splitting discretization for approximating the solution of a variant of (51) and (52).One can split (51) using the real and imaginary parts of ψ(x,t), u(x,t), and v(x,t), as(53)∂tux,t+δ1∂xxvx,t-δ2vx,twx,t+δ3u2x,t+v2x,tvx,t=0,-∂tvx,t+δ1∂xxux,t-δ2ux,twx,t+δ3u2x,t+v2x,tux,t=0,∂ttwx,t-cs2∂xxwx,t-2δ4ux,t∂xxux,t+∂xux,t2+vx,t∂xxvx,t+∂xvx,t2=0,with the initial and boundary conditions(54)uax,t=α1t,ubx,t=α2t,vax,t=α3t,vbx,t=α4t,wax,t=α5t,wbx,t=α6t,ux,at=β1x,vx,at=β2x,wx,at=β3x,wtx,at=β4x,x,t∈ax,bx×at,bt.The matrix form of the nonlinear system (53) is(55)∂tδ1∂xx0δ1∂xx-∂t000∂tt-cs2∂xxuvw+q1u,v,wq2u,v,wq3u,v,w=000,where (56)q1=-δ2vx,twx,t+δ3u2x,t+v2x,tvx,t,q2=-δ2ux,twx,t+δ3u2x,t+v2x,tux,t,q3=-2δ4ux,t∂xxux,t+∂xux,t2+vx,t∂xxvproblemx,t+∂xvx,t2,where the constants δ1, δ2, δ3, δ4, and cs are given and the functions αi(t), 1≤i≤6 and βj(x), 1≤j≤4 are known. We use J-GL-C methods for discretizing (55) subject to the initial and boundary conditions (54) to reduce (55) to a system of nonlinear algebraic equations.


## 5.2.10. Complex Generalized Zakharov Equation I

The first complex generalized Zakharov problem [3] is (57)i∂tψ-∂xxψ-ψw=0,∂ttw-∂xxw+∂xxψ2=0,with domain [-1,1]×[0,3.0]. The nonlinear problem (57) has the analytical solution (58)ψx,t=ux,t+ivx,t=3eix+ttanh⁡12x+2t,wx,t=1-tanh212x+2t.We define(59)E=maxU-Uanalytical∞,V-Vanalytical∞,W-Wanalytical∞to check the numerical accuracy in the computed solutions. Problem (57) is discretized by using Chebyshev collocation method of first kind. We took the initial guess 0 with initial and boundary conditions introduced and a grid 50×20. The achieved numerical accuracy is shown in Table 17. Our computed solutions are better than those in [3] because the maximum accuracy in [3] is of the order of 10-7. The solutions and error components are plotted in Figures 12, 13, and 14.Table 17
Maximum absolute error in the solution of the system of nonlinear equations for complex generalized Zakharov equation and varyingm.
m
E
1
8.56
e
-
1
2
2.72
e
0
3
4.16
e
0
4
8.92
e
0
5
1.12
e
1
6
3.40
e
1
7
2.04
e
1
8
2.83
e
1
9
4.05
e
0
10
2.34
e
0
11
7.77
e
-
1
12
3.25
e
-
1
13
1.14
e
-
1
14
6.20
e
-
2
15
2.79
e
-
2
16
1.85
e
-
2
17
8.11
e
-
3
18
4.44
e
-
3
19
1.79
e
-
3
20
8.86
e
-
4
21
3.91
e
-
4
22
1.66
e
-
4
23
5.98
e
-
5
24
1.97
e
-
5
25
5.32
e
-
6
26
1.20
e
-
6
27
4.12
e
-
7
28
1.36
e
-
7
29
6.95
e
-
8
30
2.95
e
-
8
31
9.29
e
-
9
32
3.42
e
-
9
33
1.43
e
-
9
34
5.57
e
-
10
35
2.62
e
-
10
36
1.26
e
-
10
37
8.65
e
-
11Figure 12
Solution and error components for complex generalized Zakharov equation (57).Figure 13
Solution and error components for complex generalized Zakharov equation (57).Figure 14
Solution and error components for complex generalized Zakharov equation (57).


## 5.2.11. Complex Generalized Zakharov Equation II

The second complex generalized Zakharov equation [3] is(60)i∂tψ+∂xxψ+2ψw-2ψ2=0,∂ttw-∂xxw+∂xxψ2=0,with the domain [-1,1]×[0,1]. That problem has the analytical solution (61)ψx,t=ux,t+ivx,t=32e-ix+3ttanh⁡x+2t,wx,t=-14tanh2x+2t.Problem (60) is discretized by using Chebyshev collocation method of first kind. We took the initial guess 0 with initial and boundary conditions introduced and a grid 21×28. The achieved numerical accuracy is shown in Table 18. Our computed solutions are better than those in [3] because the ones in [3] achieved a maximum accuracy of the order of 10-6. The solutions and error components are plotted in Figures 15, 16, and 17.Table 18
Maximum absolute error in the solution of the system of nonlinear equations for complex generalized Zakharov equation and varyingm.
m
E
1
9.2329
e
-
1
2
7.0559
e
-
1
3
5.726
e
-
1
4
6.3466
e
-
2
5
4.3852
e
-
3
6
2.3221
e
-
4
7
1.5501
e
-
5
8
4.2396
e
-
7
9
2.8906
e
-
8
10
7.7742
e
-
10
11
8.6019
e
-
11Figure 15
Solution and error components for complex generalized Zakharov equation (60).Figure 16
Solution and error components for complex generalized Zakharov equation (60).Figure 17
Solution and error components for complex generalized Zakharov equation (60).


## 5.2.12. Murray Equation

Murray equation [2] is (62)ut=uxx+λ1uux+λ2u-λ3u2,x,t∈-1,1×0,1.The initial and boundary condition for Murray equation can be computed from its analytical solution (63)ux,t=λ22λ31+tanh⁡λ28λ322λ1λ3x+λ12+4λ32t.The Chebyshev collocation method of first kind is used to discretize Murray equation over the grid 16×14. The achieved numerical accuracy is shown in Table 19. The computed accuracy in [2] is of the order of 10-8. The solution and error components are shown in Figure 18.Table 19
Maximum absolute error in the solution of the system of nonlinear equations for Murray equation and varyingm.
m
E
1
2.63
e
-
01
2
1.33
e
-
01
3
3.85
e
-
02
4
1.47
e
-
04
5
7.43
e
-
07
6
2.15
e
-
09
7
7.53
e
-
12
8
5.59
e
-
14Figure 18
Solution and error components for Murray equation.


## 5.2.13. Nonlinear Reaction and Diffusion Equation

Nonlinear reaction and diffusion equation is(64)ut=uxx+u1+u-u2,x,t∈-1,1×0,1.The analytical solution is (65)ux,t=12+52tanh⁡104x+54t.The Chebyshev collocation method of first kind is used to discretize nonlinear reaction and diffusion equation over the grid 22×14. The initial guess is 0 with initial and boundary conditions introduced. The achieved numerical accuracy is shown in Table 20. The computed accuracy in [2] is of the order of 10-8. The solution and error components are shown in Figure 19.Table 20
Maximum absolute error in the solution of the system of nonlinear equations for nonlinear reaction and diffusion equation and varyingm.
m
E
1
7.28
e
-
02
2
2.12
e
-
02
3
4.99
e
-
03
4
1.09
e
-
06
5
9.20
e
-
10
6
6.53
e
-
13
7
2.24
e
-
14Figure 19
Solution and error components for nonlinear reaction and diffusion equation.


## 5.2.14. Reaction-Diffusion Equation with Fischer-Kolmogorov Reaction Term and Density Dependent Diffusion

Reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion is(66)ut=u2xx+u1-u,x,t∈-1,1×0,1.The analytical solution is (67)ux,t=21+tanh⁡x+t/4.The Chebyshev collocation method of first kind is used to discretize reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion over the grid 30×30. The initial guess is 0 with initial and boundary conditions introduced. The achieved numerical accuracy after 1 and 2 application of NMIM is shown in Table 21. The computed accuracy in [2] is only of the order of 10-9. The solution and error components are shown in Figure 20.Table 21
Maximum absolute error in the solution of the system of nonlinear equations for reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion and varyingm.
m (1 Iter.)
E
1
1.49
e
+
00
2
1.12
e
+
00
3
1.07
e
+
00
4
4.92
e
-
01
5
3.73
e
-
01
m (2 Iter.)
E
1
1.51
e
-
02
2
7.84
e
-
03
3
6.72
e
-
04
4
5.17
e
-
07
5
1.04
e
-
09
6
2.35
e
-
12
7
5.04
e
-
14Figure 20
Solution and error components for reaction-diffusion equation with Fischer-Kolmogorov reaction term and density dependent diffusion.


## 5.2.15. Two-Dimensional Nonlinear Reaction-Diffusion Equation

The two-dimensional nonlinear reaction-diffusion equation is(68)ut=12∇2u+u2u-1,x,y,t∈-1,1×-1,1×0,20.Its analytical solution is (69)ux,y,t=11+e1/2x+y-1/2t.The Chebyshev collocation method of first kind is used to discretize the two-dimensional nonlinear reaction-diffusion equation over grid 30×30. The initial guess is 0 with initial and boundary conditions introduced. The achieved numerical accuracy is shown in Table 22. The computed accuracy in [2] is only of the order of 10-9. The error components are shown in Figure 21.Table 22
Maximum absolute error in the solution of the system of nonlinear equations for two-dimensional nonlinear reaction-diffusion equation and varyingm.
m
E
1
7.82
e
-
02
2
1.19
e
-
02
3
3.63
e
-
03
4
4.04
e
-
05
5
9.89
e
-
07
6
2.65
e
-
08
7
7.48
e
-
10Figure 21
Error components for the two-dimensional nonlinear reaction-diffusion equation.


## 5.2.16. Ostrovsky-Hunter Equation

Ostrovsky-Hunter equation is(70)ut-βuxxx+αuuxx=γu,where α is a nonlinear coefficient and β∈R, γ>0 are the dispersion coefficients. We assume a periodic boundary condition for this problem; that is, (71)u0,t=u2π,t,ux0,t=ux2π,t,uxx0,t=uxx2π,t,uxxx0,t=uxxx2π,t.The Chebyshev collocation method of first kind is used to discretize Ostrovsky-Hunter equation over a grid30×100. The initial guess is 0 with initial and boundary conditions introduced. The residue is shown in Table 23. Our solution for the grid {0,6π/11,π,16π/11}×{0.1,0.2,0.3,0.4} is given in Table 24. We notice that our results do not agree with those presented in [4]. Our solution is plotted in Figure 22.Table 23
Residue of the system of nonlinear equations for Ostrovsky-Hunter equation and varyingm.
m
F
(
x
)
∞
1
2
2
2.5963
3
2.797
4
0.72034
5
0.15165
6
0.025873
7
0.002996
8
0.00020443
9
8.7968
e
-
06
10
1.4609
e
-
07
11
1.3542
e
-
08
12
8.8853
e
-
10
13
8.5753
e
-
11Table 24
Solution for Ostrovsky-Hunter equation.
x
t
0.1
0.2
0.3
0.4
0
0.094587
0.17222
0.20233
0.21608
6
π
/
11
0.96556
0.94696
0.8702
0.73399
π
−0.10438
−0.22921
−0.41929
−0.67358
16
π
/
11
−0.98195
−0.96996
−0.97341
−0.97748Figure 22
Solution for Ostrovsky-Hunter equation.


## 6. Conclusions

Frozen Jacobian multistep iterative methods are computationally efficient to solve systems of nonlinear equations. This is because the LU factors of the Jacobian can be reused to solve additional linear systems at a very small additional cost. We have developed a new frozen Jacobian multistep iterative method with increased order of convergence for a given number of function evaluations. The claimed order of convergence of NMIM has been confirmed computationally by solving several small systems of nonlinear equations. NMIM has been used to solve fifteen nonlinear IVPs and BVPs with different collocation methods. As a whole, the Chebyshev collocation method of the first kind seems to provide the best numerical accuracy in the set of fifteen IVPs and BVPs considered in the paper.



---

*Source: 9407656-2017-05-03.xml*