Mathematical Problems in Engineering
Volume 2009 (2009), Article ID 307298, 24 pages
doi:10.1155/2009/307298
Research Article

Modified Jacobian Newton Iterative Method: Theory and Applications

Jürgen Geiser

Department of Mathematics, Humboldt-Universität zu Berlin, Unter den Linden 6, D-10099 Berlin, Germany

Received 3 September 2008; Revised 10 November 2008; Accepted 4 December 2008

Academic Editor: José Roberto Castilho Piqueira

Copyright © 2009 Jürgen Geiser. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

This article proposes a new approach to the construction of a linearization method based on the iterative operator-splitting method for nonlinear differential equations. The convergence properties of such a method are studied. The main features of the proposed idea are the linearization of nonlinear equations and the application of iterative splitting methods. We present an iterative operator-splitting method with embedded Newton methods to solve nonlinearity. We confirm with numerical applications the effectiveness of the proposed iterative operator-splitting method in comparison with the classical Newton methods. We provide improved results and convergence rates.

1. Introduction

In this paper we propose a modified Jacobian-Newton iterative method to solve nonlinear differential equations. In the first paper we concentrate on ordinary differential equations, but numerical results are also obtained for partial differential equations. Basic studies of the operator-splitting methods are found in [1, 2]. Further important research was carried out to obtain a higher order for the splitting methods (see [3]). For this reason, the iterative splitting methods became more important for linear and nonlinear differential equations, while simple increasing of iteration steps affects the order of the scheme (see [4]). An interesting topic is Newton's methods for nonlinear problems with specifications for numerical implementations, (see [5]). Efficient modifications of Newton's methods, with regard to the computation of the Jacobian matrices are discussed in [6]. In our paper we discuss the benefit of the combination of splitting and linearization methods (see theoretical frameworks [7, 8]).

The outline of the paper is as follows. For our mathematical model we describe the convection-diffusion-reaction equation in Section 2. The fractional splitting is introduced in Section 3. We present the iterative splitting methods in Section 4. Section 5 discusses the Newton methods and their modifications. In Section 6 we present the numerical results from the solution of selected model problems. We end the article in Section 7 with a conclusion and comments.

2. Mathematical Model

The motivation for the study presented below originates from a computational simulation of heat-transfer [9] and convection-diffusion-reaction-equations [10–13].

In the present paper we concentrate on ordinary differential equations, given as 𝜕 𝑡 𝑢 ( 𝑡 ) = 𝐴 ( 𝑢 ( 𝑡 ) ) 𝑢 ( 𝑡 ) + 𝐵 ( 𝑢 ( 𝑡 ) ) 𝑢 ( 𝑡 ) , 𝑡 ∈ ( 0 , 𝑇 ) , ( 2 . 1 ) where the initial condition is 𝑢 ( 0 ) = 𝑢 0 . The operators 𝐴 ( 𝑢 ) and 𝐵 ( 𝑢 ) can be spatially discretized operators, that is, they can correspond to the discretized in space convection and diffusion operators (matrices). In the following, we deal with bounded nonlinear operators.

The aim of this paper is to present a new method based on Newton and iterative schemes.

In the next section we discuss the decoupling of the time-scale with a first-order fractional splitting method.

3. Fractional-Splitting Methods of First-Order for Linear Equations

First we describe the simplest operator-splitting, which is called sequential operator-splitting, for the following linear system of ordinary differential equations: 𝜕 𝑡 𝑢 ( 𝑡 ) = 𝐴 𝑢 ( 𝑡 ) + 𝐵 𝑢 ( 𝑡 ) , 𝑡 ∈ ( 0 , 𝑇 ) , ( 3 . 1 ) where the initial condition is 𝑢 ( 0 ) = 𝑢 0 . The operators 𝐴 and 𝐵 are linear and bounded operators in a Banach space (see also Section 2).

The sequential operator-splitting method is introduced as a method that solves two subproblems sequentially, where the different subproblems are connected via the initial conditions. This means that we replace the original problem (3.1) with the subproblems 𝜕 𝑢 ∗ ( 𝑡 ) 𝜕 𝑡 = 𝐴 𝑢 ∗ ( 𝑡 ) , w i t h 𝑢 ∗ 𝑡 𝑛 = 𝑢 𝑛 , 𝜕 𝑢 ∗ ∗ ( 𝑡 ) 𝜕 𝑡 = 𝐵 𝑢 ∗ ∗ ( 𝑡 ) , w i t h 𝑢 ∗ ∗ 𝑡 𝑛 = 𝑢 ∗ 𝑡 𝑛 + 1 , ( 3 . 2 ) where the splitting time-step is defined as 𝜏 𝑛 = 𝑡 𝑛 + 1 − 𝑡 𝑛 . The approximated solution is 𝑢 𝑛 + 1 = 𝑢 ∗ ∗ ( 𝑡 𝑛 + 1 ) .

Clearly, the replacement of the original problem with the subproblems usually results in an error, called splitting error. The splitting error of the sequential operator-splitting method can be derived as (cf., e.g., [1, 2]). 𝜌 𝑛 = 1 𝜏 𝑛 𝜏 e x p 𝑛 𝜏 ( 𝐴 + 𝐵 ) − e x p 𝑛 𝐵 𝜏 e x p 𝑛 𝐴 𝑢 𝑡 𝑛 = 𝑂 𝜏 0 , f o r [ 𝐴 , 𝐵 ] = 0 , 𝑛 , f o r [ 𝐴 , 𝐵 ] ≠ 0 , ( 3 . 3 ) where [ 𝐴 , 𝐵 ] ∶ = 𝐴 𝐵 − 𝐵 𝐴 is the commutator of 𝐴 and 𝐵 . Consequently, the splitting error is 𝑂 ( 𝜏 𝑛 ) when the operators 𝐴 and 𝐵 do not commute, otherwise the method is exact. Hence, by definition, the sequential operator-splitting is called the first-order splitting method.

4. The Iterative Splitting Method

The following algorithm is based on the iteration with fixed splitting discretization step-size 𝜏 . On the time interval [ 𝑡 𝑛 , 𝑡 𝑛 + 1 ] we solve the following subproblems consecutively for 𝑖 = 0 , 2 , … 2 𝑚 : 𝜕 𝑢 𝑖 ( 𝑥 , 𝑡 ) 𝜕 𝑡 = 𝐴 𝑢 𝑖 ( 𝑥 , 𝑡 ) + 𝐵 𝑢 𝑖 − 1 ( 𝑥 , 𝑡 ) , w i t h 𝑢 𝑖 ( 𝑡 𝑛 ) = 𝑢 𝑛 , 𝑢 0 ( 𝑥 , 𝑡 𝑛 ) = 𝑢 𝑛 , 𝑢 − 1 𝑢 = 0 , 𝑖 ( 𝑥 , 𝑡 ) = 𝑢 𝑖 − 1 ( 𝑥 , 𝑡 ) = 𝑢 1 , o n 𝜕 Ω × ( 0 , 𝑇 ) , 𝜕 𝑢 𝑖 + 1 ( 𝑥 , 𝑡 ) 𝜕 𝑡 = 𝐴 𝑢 𝑖 ( 𝑥 , 𝑡 ) + 𝐵 𝑢 𝑖 + 1 ( 𝑥 , 𝑡 ) , w i t h 𝑢 𝑖 + 1 ( 𝑥 , 𝑡 𝑛 ) = 𝑢 𝑛 , 𝑢 𝑖 ( 𝑥 , 𝑡 ) = 𝑢 𝑖 − 1 ( 𝑥 , 𝑡 ) = 𝑢 1 , o n 𝜕 Ω × ( 0 , 𝑇 ) , ( 4 . 1 ) where 𝑢 𝑛 is the known split approximation at the time level 𝑡 = 𝑡 𝑛 (see [14]).

Remark 4.1. We can generalize the iterative splitting method to a multi-iterative splitting method by introducing new splitting operators, for example, spatial operators. Then we obtain multi-indices to control the splitting process; each iterative splitting method can be solved independently, while connecting with further steps to the multi-splitting methods. In the following we introduce the multi-iterative splitting method for a combined time-space splitting method.

5. The Modified Jacobian-Newton Methods and Fixpoint-Iteration Methods

In this section we describe the modified Jacobian-Newton methods and Fixpoint-iteration methods.

We propose for weak nonlinearities, for example, quadratic nonlinearity, the fixpoint iteration method, where our iterative operator splitting method is one, see [4]. For stronger nonlinearities, for example, cubic or higher order polynomial nonlinearities, the modified Jacobian method with embedded iterative-splitting methods is suggested.

The point of embedding the splitting methods into the Newton methods is to decouple the equation systems into simpler equations. Such simple equation systems can be solved with scalar Newton methods.

5.1. The Altered Jacobian-Newton Iterative Methods with Embedded Sequential Splitting Methods

We confine our attention to time-dependent partial differential equations of the form 𝑑 𝑐 𝑡 𝑑 𝑡 = 𝐴 ( 𝑐 ( 𝑡 ) ) 𝑐 ( 𝑡 ) + 𝐵 ( 𝑐 ( 𝑡 ) ) 𝑐 ( 𝑡 ) , w i t h 𝑐 𝑛 = 𝑐 𝑛 , ( 5 . 1 ) where 𝐴 ( 𝑐 ) , 𝐵 ( 𝑐 ) ∶ 𝐗 → 𝐗 are linear and densely defined in the real Banach space 𝐗 , involving only spatial derivatives of 𝑐 , see [8]. We assume also that we have a weak nonlinear operator with 𝐴 ( 𝑐 ) 𝑐 = 𝜆 1 𝑐 and 𝐵 ( 𝑐 ) 𝑐 = 𝜆 2 𝑐 , where 𝜆 1 and 𝜆 2 are constant factors.

In the following we discuss the embedding of a sequential splitting method into the Newton method.

The altered Jacobian-Newton iterative method with an embedded iterative splitting method is given as follows.

Newton's Method
𝐹 ( 𝑐 ) = 𝑑 𝑐 / 𝑑 𝑡 − 𝐴 ( 𝑐 ( 𝑡 ) ) 𝑐 ( 𝑡 ) − 𝐵 ( 𝑐 ( 𝑡 ) ) 𝑐 ( 𝑡 ) and we can compute 𝑐 ( 𝑘 + 1 ) = 𝑐 ( 𝑘 ) − 𝐷 ( 𝐹 ( 𝑐 ( 𝑘 ) ) ) − 1 𝐹 ( 𝑐 ( 𝑘 ) ) , where 𝐷 ( 𝐹 ( 𝑐 ) ) is the Jacobian matrix and 𝑘 = 0 , 1 , … .

We stop the iterations when we obtain: | 𝑐 ( 𝑘 + 1 ) − 𝑐 ( 𝑘 ) | ≤ e r r , where e r r is an error bound, for example, e r r = 1 0 − 4 .

We assume the spatial discretization, with spatial grid points, 𝑖 = 1 , … , 𝑚 and obtain the differential equation system: ⎛ ⎜ ⎜ ⎜ ⎝ 𝐹 𝑐 𝐹 ( 𝑐 ) = 1 𝐹 𝑐 2 ⋮ 𝐹 𝑐 𝑚 ⎞ ⎟ ⎟ ⎟ ⎠ , ( 5 . 2 ) where 𝑐 = ( 𝑐 1 , … , 𝑐 𝑚 ) 𝑇 and 𝑚 is the number of spatial grid points.

The Jacobian matrix for the equation system is given as: ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 𝐷 𝐹 ( 𝑐 ) = 𝜕 𝐹 ( 𝑐 1 ) 𝑐 1 𝜕 𝐹 ( 𝑐 1 ) 𝑐 2 … 𝜕 𝐹 ( 𝑐 1 ) 𝑐 𝑚 𝜕 𝐹 ( 𝑐 2 ) 𝑐 1 𝜕 𝐹 ( 𝑐 2 ) 𝑐 2 … 𝜕 𝐹 ( 𝑐 2 ) 𝑐 𝑚 ⋮ 𝜕 𝐹 ( 𝑐 𝑚 ) 𝑐 1 𝜕 𝐹 ( 𝑐 𝑚 ) 𝑐 2 … 𝜕 𝐹 ( 𝑐 𝑚 ) 𝑐 𝑚 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ , ( 5 . 3 ) where 𝑐 = ( 𝑐 1 , … , 𝑐 𝑚 ) .

The modified Jacobian is given as: ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 𝐷 𝐹 ( 𝑐 ) = 𝜕 𝐹 ( 𝑐 1 ) 𝑐 1 + 𝐹 ( 𝑐 1 ) 𝜕 𝐹 ( 𝑐 1 ) 𝑐 2 … 𝜕 𝐹 ( 𝑐 1 ) 𝑐 𝑚 𝜕 𝐹 ( 𝑐 2 ) 𝑐 1 𝜕 𝐹 ( 𝑐 2 ) 𝑐 2 𝐹 ( 𝑐 2 ) … 𝜕 𝐹 ( 𝑐 2 ) 𝑐 𝑚 ⋮ 𝜕 𝐹 ( 𝑐 𝑚 ) 𝑐 1 𝜕 𝐹 ( 𝑐 𝑚 ) 𝑐 2 … 𝜕 𝐹 ( 𝑐 𝑚 ) 𝑐 𝑚 + 𝐹 ( 𝑐 𝑚 ) ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ , ( 5 . 4 ) where 𝑐 = ( 𝑐 1 , … , 𝑐 𝑛 ) .

By embedding the sequential splitting method we obtain the following algorithm. We decouple into two equation systems: 𝐹 1 𝑢 1 = 𝜕 𝑡 𝑢 1 𝑢 − 𝐴 1 𝑢 1 = 0 w i t h 𝑢 1 𝑡 𝑛 = 𝑐 𝑛 , 𝐹 2 𝑢 2 = 𝜕 𝑡 𝑢 2 𝑢 − 𝐵 2 𝑢 2 = 0 w i t h 𝑢 2 𝑡 𝑛 = 𝑢 1 𝑡 𝑛 + 1 , ( 5 . 5 ) where the results of the methods are 𝑐 ( 𝑡 𝑛 + 1 ) = 𝑢 2 ( 𝑡 𝑛 + 1 ) , and 𝑢 1 = ( 𝑢 1 1 , … , 𝑢 1 𝑛 ) , 𝑢 2 = ( 𝑢 2 1 , … , 𝑢 2 𝑛 ) .

Thus we have to solve two Newton methods, each in one equations system. The contribution is to reduce the Jacobian matrix into a diagonal entry, for example, with a weighted Newton method, see [15]. The splitting method with embedded Newton method is given as follows.

Algorithm 5.1. We assume the spatial operators 𝐴 and 𝐵 are discretized, for example, finite difference or finite element methods; further all initial conditions and boundary conditions are discrete given. Then we can apply the Newton's method in its discrete form as: 𝑢 1 ( 𝑘 + 1 ) = 𝑢 1 ( 𝑘 ) 𝐹 − 𝐷 1 𝑢 1 ( 𝑘 ) − 1 𝜕 𝑡 𝑢 1 ( 𝑘 ) 𝑢 − 𝐴 1 ( 𝑘 ) 𝑢 1 ( 𝑘 ) , 𝐹 w i t h 𝐷 1 𝑢 1 ( 𝑘 ) = 𝜕 𝜕 𝑢 1 ( 𝑘 ) 𝜕 𝑡 𝑢 1 ( 𝑘 ) 𝑢 − 𝐴 1 ( 𝑘 ) − 𝑢 𝜕 𝐴 1 ( 𝑘 ) 𝜕 𝑢 1 ( 𝑘 ) 𝑢 1 ( 𝑘 ) , 𝑢 1 ( 𝑘 ) 𝑡 𝑛 = 𝑐 𝑛 𝑢 , 𝑘 = 0 , 1 , 2 , … , 𝐾 , 2 ( 𝑙 + 1 ) = 𝑢 2 ( 𝑙 ) 𝐹 − 𝐷 2 𝑢 2 ( 𝑙 ) − 1 𝜕 𝑡 𝑢 2 ( 𝑙 ) 𝑢 − 𝐵 2 ( 𝑙 ) 𝑢 2 ( 𝑙 ) , 𝐹 w i t h 𝐷 2 𝑢 2 ( 𝑙 ) = 𝜕 𝜕 𝑢 1 ( 𝑘 ) 𝜕 𝑡 𝑢 2 ( 𝑘 ) 𝑢 − 𝐵 2 ( 𝑙 ) − 𝑢 𝜕 𝐵 2 ( 𝑙 ) 𝜕 𝑢 2 ( 𝑙 ) 𝑢 2 ( 𝑙 ) , 𝑢 2 ( 𝑙 ) 𝑡 𝑛 = 𝑢 𝐾 1 𝑡 𝑛 + 1 , 𝑙 = 0 , 1 , 2 , … , 𝐿 . ( 5 . 6 ) where 𝑘 and 𝑙 are the iteration indices, 𝐾 and 𝐿 the maximal iterative steps for each part of the Newton's method. The maximal iterative steps allow us to have at least an error of: | | 𝑢 ( 𝐾 ) ( 𝑡 𝑛 + 1 ) 1 − 𝑢 ( 𝐾 − 1 ) ( 𝑡 𝑛 + 1 ) 1 | | | | 𝑢 ≤ e r r , ( 𝐿 ) ( 𝑡 𝑛 + 1 ) 2 − 𝑢 ( 𝐿 − 1 ) ( 𝑡 𝑛 + 1 ) 2 | | ≤ e r r , ( 5 . 7 ) where e r r is the error bound, for example, e r r = 1 0 − 6 .
The approximated solution is given as: 𝑢 𝑡 𝑛 + 1 = 𝑢 ( 𝐿 ) ( 𝑡 𝑛 + 1 ) 2 . ( 5 . 8 )

For the improvement method, we can apply the weighted Newton method. We try to skip the delicate outer diagonals in the Jacobian matrix and apply: 𝑢 1 ( 𝑘 + 1 ) = 𝑢 1 ( 𝑘 ) − 𝐷 𝐹 1 𝑢 1 ( 𝑘 ) + 𝛿 1 𝑢 1 ( 𝑘 ) − 1 𝐹 1 𝑢 1 ( 𝑘 ) + 𝜖 𝑢 1 ( 𝑘 ) , ( 5 . 8 ) where the function 𝛿 can be applied as a scalar, for example, 𝛿 = 1 0 − 6 , and the same with 𝜖 . It is important to ensure that 𝛿 is small enough to preserve the convergence.

Remark 5.2. If we assume that we discretize (5.5) with the backward-Euler method, for example, 𝐹 1 𝑢 1 𝑡 𝑛 + 1 = 𝑢 1 𝑡 𝑛 + 1 − 𝑢 1 𝑡 𝑛 𝑢 − Δ 𝑡 𝐴 1 𝑡 𝑛 + 1 𝑢 1 𝑡 𝑛 + 1 = 0 w i t h 𝑢 1 𝑡 𝑛 = 𝑐 𝑛 , 𝐹 2 𝑢 2 = 𝑢 2 𝑡 𝑛 + 1 − 𝑢 2 𝑡 𝑛 𝑢 − Δ 𝑡 𝐵 2 𝑡 𝑛 + 1 𝑢 2 𝑡 𝑛 + 1 = 0 w i t h 𝑢 2 𝑡 𝑛 = 𝑢 1 𝑡 𝑛 + 1 , ( 5 . 1 0 ) then we obtain the derivations 𝐷 ( 𝐹 1 ( 𝑢 1 ( 𝑡 𝑛 + 1 ) ) ) and 𝐷 ( 𝐹 2 ( 𝑢 2 ( 𝑡 𝑛 + 1 ) ) ) 𝐷 𝐹 1 𝑢 1 𝑡 𝑛 + 1 𝐴 𝑢 = 1 − Δ 𝑡 1 𝑡 𝑛 + 1 + 𝑢 𝜕 𝐴 1 𝑡 𝑛 + 1 𝜕 𝑢 1 𝑡 𝑛 + 1 𝑢 1 𝑡 𝑛 + 1 , 𝐷 𝐹 2 𝑢 2 𝐵 𝑢 = 1 − Δ 𝑡 2 𝑡 𝑛 + 1 + 𝑢 𝜕 𝐵 2 𝑡 𝑛 + 1 𝜕 𝑢 2 𝑡 𝑛 + 1 𝑢 2 𝑡 𝑛 + 1 . ( 5 . 1 1 )

We can apply the equation (5.9) analogously 𝑢 2 ( 𝑙 + 1 ) .

5.2. Iterative Operator-Splitting Method as a Fixpoint Scheme

The iterative operator-splitting method is used as a fixpoint scheme to linearize the nonlinear operators, see [4, 16].

We confine our attention to time-dependent partial differential equations of the form: 𝑑 𝑢 𝑑 𝑡 = 𝐴 ( 𝑢 ( 𝑡 ) ) 𝑢 ( 𝑡 ) + 𝐵 ( 𝑢 ( 𝑡 ) ) 𝑢 ( 𝑡 ) , w i t h 𝑢 ( 𝑡 𝑛 ) = 𝑐 𝑛 , ( 5 . 1 2 ) where 𝐴 ( 𝑢 ) , 𝐵 ( 𝑢 ) ∶ 𝐗 → 𝐗 are linear and densely defined in the real Banach space 𝐗 , involving only spatial derivatives of 𝑐 , see [8]. In the following we discuss the standard iterative operator-splitting methods as a fixpoint iteration method to linearize the operators.

We split our nonlinear differential equation (5.12) by applying: 𝑑 𝑢 𝑖 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐴 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 𝑢 ( 𝑡 ) + 𝐵 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 − 1 ( 𝑡 ) , w i t h 𝑢 𝑖 𝑡 𝑛 = 𝑐 𝑛 , 𝑑 𝑢 𝑖 + 1 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐴 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 𝑢 ( 𝑡 ) + 𝐵 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 + 1 , w i t h 𝑢 𝑖 + 1 𝑡 𝑛 = 𝑐 𝑛 , ( 5 . 1 3 ) where the time-step is 𝜏 = 𝑡 𝑛 + 1 − 𝑡 𝑛 . The iterations are 𝑖 = 1 , 3 , … , 2 𝑚 + 1 . 𝑢 0 ( 𝑡 ) = 𝑐 𝑛 is the starting solution, where we assume the solution 𝑐 𝑛 + 1 is near 𝑐 𝑛 , or 𝑢 0 ( 𝑡 ) = 0 . So we have to solve the local fixpoint problem. 𝑐 𝑛 is the known split approximation at the time level 𝑡 = 𝑡 𝑛 .

The split approximation at time level 𝑡 = 𝑡 𝑛 + 1 is defined as 𝑐 𝑛 + 1 = 𝑢 2 𝑚 + 2 ( 𝑡 𝑛 + 1 ) . We assume the operators 𝐴 ( 𝑢 𝑖 − 1 ) , 𝐵 ( 𝑢 𝑖 − 1 ) ∶ 𝐗 → 𝐗 to be linear and densely defined on the real Banach space 𝐗 , for 𝑖 = 1 , 3 , … , 2 𝑚 + 1 .

Here the linearization is done with respect to the iterations, such that 𝐴 ( 𝑢 𝑖 − 1 ) , 𝐵 ( 𝑢 𝑖 − 1 ) are at least non-dependent operators in the iterative equations, and we can apply the linear theory.

The linearization is at least in the first equation 𝐴 ( 𝑢 𝑖 − 1 ) ≈ 𝐴 ( 𝑢 𝑖 ) , and in the second equation 𝐵 ( 𝑢 𝑖 − 1 ) ≈ 𝐵 ( 𝑢 𝑖 + 1 ) .

We have ‖ ‖ 𝐴 𝑢 𝑖 − 1 𝑡 𝑛 + 1 𝑢 𝑖 𝑡 𝑛 + 1 𝑢 − 𝐴 𝑛 + 1 𝑢 𝑡 𝑛 + 1 ‖ ‖ ≤ 𝜖 , ( 5 . 1 4 ) with sufficient iterations (5.14) 𝑖 = { 1 , 3 , … , 2 𝑚 + 1 } .

RemarkThe linearization with the fixpoint scheme can be used for smooth or weak nonlinear operators, otherwise we lose the convergence behavior, while we did not converge to the local fixpoint, see [].

The second idea is based on the Newton method.

5.3. Jacobian-Newton Iterative Method with Embedded Operator-Splitting Method

The Newton method is used to solve the nonlinear parts of the iterative operator-splitting method (see the linearization techniques in [, ]).

Newton Method
The function is given as: 𝐹 ( 𝑐 ) = 𝜕 𝑐 𝜕 𝑡 − 𝐴 ( 𝑐 ( 𝑡 ) ) 𝑐 ( 𝑡 ) − 𝐵 ( 𝑐 ( 𝑡 ) ) 𝑐 ( 𝑡 ) = 0 , ( 5 . 1 4 )
The iteration can be computed as: 𝑐 ( 𝑘 + 1 ) = 𝑐 ( 𝑘 ) 𝐹 𝑐 − 𝐷 ( 𝑘 ) − 1 𝐹 𝑐 ( 𝑘 ) , ( 5 . 1 4 ) where (5.14) 𝐷 ( 𝐹 ( 𝑐 ) ) is the Jacobian matrix and (5.14) 𝑘 = 0 , 1 , … . and (5.14) 𝑐 = ( 𝑐 1 , … , 𝑐 𝑚 ) is the solution vector of the spatial discretized nonlinear equation.

We then have to apply the iterative operator-splitting method and obtain: 𝐹 1 𝑢 𝑖 = 𝜕 𝑡 𝑢 𝑖 𝑢 − 𝐴 𝑖 𝑢 𝑖 𝑢 − 𝐵 𝑖 − 1 𝑢 𝑖 − 1 = 0 , w i t h 𝑢 𝑖 𝑡 𝑛 = 𝑐 𝑛 , 𝐹 2 𝑢 𝑖 + 2 = 𝜕 𝑡 𝑢 𝑖 + 1 𝑢 − 𝐴 𝑖 𝑢 𝑖 𝑢 − 𝐵 𝑖 + 1 𝑢 𝑖 + 1 = 0 , w i t h 𝑢 𝑖 + 1 𝑡 𝑛 = 𝑐 𝑛 , ( 5 . 1 4 ) where the time-step is 𝜏 = 𝑡 𝑛 + 1 − 𝑡 𝑛 . The iterations are 𝑖 = 1 , 3 , … , 2 𝑚 + 1 . 𝑐 0 ( 𝑡 ) = 0 is the starting solution and 𝑐 𝑛 is the known split approximation at the time-level 𝑡 = 𝑡 𝑛 . The results of the methods are 𝑐 ( 𝑡 𝑛 + 1 ) = 𝑢 2 𝑚 + 2 ( 𝑡 𝑛 + 1 ) .

Thus we have to solve two Newton methods and the contribution will be to reduce the Jacobian matrix, for example, skip the diagonal entries. The splitting method with the embedded Newton method is given as: 𝑢 𝑖 ( 𝑘 + 1 ) = 𝑢 𝑖 ( 𝑘 ) 𝐹 − 𝐷 1 𝑢 𝑖 ( 𝑘 ) − 1 𝜕 𝑡 𝑢 𝑖 ( 𝑘 ) 𝑢 − 𝐴 𝑖 ( 𝑘 ) 𝑢 𝑖 ( 𝑘 ) 𝑢 − 𝐵 ( 𝑘 ) 𝑖 − 1 𝑢 ( 𝑘 ) 𝑖 − 1 , 𝐹 w i t h 𝐷 1 𝑢 𝑖 ( 𝑘 ) 𝐴 𝑢 = − 𝑖 ( 𝑘 ) + 𝑢 𝜕 𝐴 𝑖 ( 𝑘 ) 𝜕 𝑢 𝑖 ( 𝑘 ) 𝑢 𝑖 ( 𝑘 ) , 𝑘 = 0 , 1 , 2 , … , 𝐾 , w i t h 𝑢 𝑖 𝑡 𝑛 = 𝑐 𝑛 , 𝑢 ( 𝑙 + 1 ) 𝑖 + 1 = 𝑢 ( 𝑙 ) 𝑖 + 1 𝐹 − 𝐷 2 𝑢 ( 𝑙 ) 𝑖 + 1 − 1 𝜕 𝑡 𝑢 ( 𝑙 ) 𝑖 + 1 𝑢 − 𝐴 𝑖 ( 𝑘 ) 𝑢 𝑖 ( 𝑘 ) 𝑢 − 𝐵 ( 𝑘 ) 𝑖 + 1 𝑢 ( 𝑘 ) 𝑖 + 1 𝑐 2 ( 𝑙 ) , 𝐹 w i t h 𝐷 2 𝑢 ( 𝑙 ) 𝑖 + 1 𝐵 𝑢 = − ( 𝑙 ) 𝑖 + 1 + 𝑢 𝜕 𝐵 ( 𝑙 ) 𝑖 + 1 𝜕 𝑢 ( 𝑙 ) 𝑖 + 1 𝑢 ( 𝑙 ) 𝑖 + 1 , 𝑙 = 0 , 1 , 2 , … , 𝐿 , w i t h 𝑢 𝑖 + 1 𝑡 𝑛 = 𝑐 𝑛 , ( 5 . 1 8 ) where the time-step is 𝜏 = 𝑡 𝑛 + 1 − 𝑡 𝑛 . The iterations are: 𝑖 = 1 , 3 , … , 2 𝑚 + 1 . 𝑐 0 ( 𝑡 ) = 0 is the starting solution and 𝑐 𝑛 is the known split approximation at the time-level 𝑡 = 𝑡 𝑛 . The results of the methods are 𝑐 ( 𝑡 𝑛 + 1 ) = 𝑢 2 𝑚 + 2 ( 𝑡 𝑛 + 1 ) .

For the improvement by skipping the delicate outer diagonals in the Jacobian matrix, we apply 𝑢 𝑖 ( 𝑘 + 1 ) = 𝑢 𝑖 ( 𝑘 ) − ( 𝐷 ( 𝐹 1 ( 𝑢 𝑖 ( 𝑘 ) ) ) + 𝛿 1 ( 𝑢 𝑖 ( 𝑘 ) ) ) − 1 ( 𝐹 1 ( 𝑢 𝑖 ( 𝑘 ) ) + 𝜖 𝑢 𝑖 ( 𝑘 ) ) , and analogously 𝑢 ( 𝑙 + 1 ) 𝑖 + 1 .

Remark 5.4. For the iterative operator-splitting method with the Newton iteration we have two iteration procedures. The first iteration is the Newton method to compute the solution of the nonlinear equations, and the second iteration is the iterative splitting method, which computes the resulting solution of the coupled equation systems. The embedded method is used for strong nonlinearities.

6. Numerical Results

In this section, we present the numerical results for nonlinear differential equation using several variations of the proposed Newton and iterative schemes as solvers.

6.1. First Numerical Example: Bernoulli Equation

As a nonlinear differential example, we choose the Bernoulli equation: 𝜕 𝑢 ( 𝑡 ) = 𝜆 𝜕 𝑡 1 + 𝜆 3 𝜆 𝑢 ( 𝑡 ) + 2 + 𝜆 4 ( 𝑢 ( 𝑡 ) ) 𝑝 , 𝑡 ∈ [ 0 , 𝑇 ] , w i t h 𝑢 ( 0 ) = 1 , ( 6 . 1 ) where the analytical solution can be derived as (see also [16]): 𝜆 𝑢 ( 𝑡 ) = e x p 1 + 𝜆 3 𝑡 − 𝜆 2 + 𝜆 4 𝜆 1 + 𝜆 3 𝜆 e x p 1 + 𝜆 3 ( 𝑝 − 1 ) 𝑡 + 𝑐 1 / ( 1 − 𝑝 ) . ( 6 . 2 )

Using 𝑢 ( 0 ) = 1 we find that 𝑐 = 1 + ( 𝜆 2 + 𝜆 4 ) / ( 𝜆 1 + 𝜆 3 ) , so 𝜆 𝑢 ( 𝑡 ) = e x p 1 + 𝜆 3 𝑡 𝜆 1 + 2 + 𝜆 4 𝜆 1 + 𝜆 3 𝜆 1 − e x p 1 + 𝜆 3 ( 𝑝 − 1 ) 𝑡 1 / ( 1 − 𝑝 ) . ( 6 . 3 )

We choose 𝑝 = 2 , 𝜆 1 = − 1 , 𝜆 2 = − 0 . 5 , 𝜆 3 = − 1 0 0 , 𝜆 4 = − 2 0 and, for example, Δ 𝑡 = 1 0 − 2 .

The analytical solutions can be given as: 𝑢 ( 𝑡 ) 1 − 𝑝 = 𝑢 0 𝜆 e x p ( 1 − 𝑝 ) 1 + 𝜆 3 𝑡 + 𝜆 2 + 𝜆 4 𝜆 1 + 𝜆 3 𝜆 e x p ( 1 − 𝑝 ) 1 + 𝜆 3 𝑡 − 1 ) . ( 6 . 4 )

We divide the time interval [ 0 , 𝑇 ] , with 𝑇 = 1 , in 𝑛 intervals with length 𝜏 𝑛 = 𝑇 / 𝑛 .

(1) The sequential operator-splitting method with analytical solutions is given as follows.

We apply the quasilinear iterative operator-splitting method: 𝑑 𝑢 1 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐴 1 𝑢 ( 𝑡 ) 1 ( 𝑡 ) , w i t h 𝑢 1 𝑡 𝑛 = 𝑢 𝑛 , 𝑑 𝑢 2 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐵 2 𝑢 ( 𝑡 ) 2 , w i t h 𝑢 2 𝑡 𝑛 = 𝑢 1 𝑛 + 1 , ( 6 . 5 ) with the nonlinear operators 𝐴 ( 𝑢 ) 𝑢 = 𝜆 1 𝑢 ( 𝑡 ) + 𝜆 2 ( 𝑢 ( 𝑡 ) ) 𝑝 − 1 𝑢 , 𝐵 ( 𝑢 ) 𝑢 = 𝜆 3 𝑢 ( 𝑡 ) + 𝜆 4 ( 𝑢 ( 𝑡 ) ) 𝑝 − 1 𝑢 . The result is given as 𝑢 2 ( 𝑡 𝑛 + 1 ) = 𝑢 𝑛 + 1 .

We apply the Newton method and discretize the operators with time discretization methods such as backward-Euler or higher Runge-Kutta methods.

The analytical result for each equation part is given as: 𝑢 1 ( 𝑡 ) 1 − 𝑝 𝑡 = 𝑢 𝑛 𝜆 e x p ( 1 − 𝑝 ) 1 𝑡 + 𝜆 2 𝜆 1 𝜆 e x p ( 1 − 𝑝 ) 1 𝑡 , 𝑢 − 1 2 𝑡 𝑛 + 1 1 − 𝑝 = 𝑢 1 𝑡 𝑛 + 1 𝜆 e x p ( 1 − 𝑝 ) 3 𝑡 + 𝜆 4 𝜆 3 𝜆 e x p ( 1 − 𝑝 ) 3 𝑡 , − 1 ( 6 . 6 ) where the result is given as 𝑢 ( 𝑡 𝑛 + 1 ) = 𝑢 2 ( 𝑡 𝑛 + 1 ) .

We can apply the simpler equations and solve the sequential operator-splitting method.

(2) The sequential operator-splitting method with embedded Newton method is given as follows.

We apply the quasilinear iterative operator-splitting method: 𝑑 𝑢 1 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐴 1 𝑢 ( 𝑡 ) 1 ( 𝑡 ) , w i t h 𝑢 1 𝑡 𝑛 = 𝑢 𝑛 , 𝑑 𝑢 2 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐵 2 𝑢 ( 𝑡 ) 2 , w i t h 𝑢 2 𝑡 𝑛 = 𝑢 1 𝑛 + 1 , ( 6 . 7 ) with the nonlinear operators 𝐴 ( 𝑢 ) 𝑢 = 𝜆 1 𝑢 ( 𝑡 ) + 𝜆 2 ( 𝑢 ( 𝑡 ) ) 𝑝 − 1 𝑢 , 𝐵 ( 𝑢 ) 𝑢 = 𝜆 3 𝑢 ( 𝑡 ) + 𝜆 4 ( 𝑢 ( 𝑡 ) ) 𝑝 − 1 𝑢 . The result is given as 𝑢 2 ( 𝑡 𝑛 + 1 ) = 𝑢 𝑛 + 1 .

We apply the Newton method and discretize the operators with time discretization methods such as backward-Euler or higher Runge-Kutta methods.

The splitting method with embedded Newton's method is given as 𝑢 1 ( 𝑘 + 1 ) = 𝑢 1 ( 𝑘 ) 𝐹 − 𝐷 1 𝑢 1 ( 𝑘 ) − 1 𝜕 𝑡 𝑢 1 ( 𝑘 ) 𝑢 − 𝐴 1 ( 𝑘 ) 𝑢 1 ( 𝑘 ) , 𝐹 w i t h 𝐷 1 𝑢 1 ( 𝑘 ) 𝐴 𝑢 = − 1 ( 𝑘 ) + 𝑢 𝜕 𝐴 1 ( 𝑘 ) 𝜕 𝑢 1 ( 𝑘 ) 𝑢 1 ( 𝑘 ) , 𝑢 1 ( 𝑘 ) 𝑡 𝑛 = 𝑐 𝑛 𝑢 , 𝑘 = 0 , 1 , 2 , … , 𝐾 , 2 ( 𝑙 + 1 ) = 𝑢 2 ( 𝑙 ) 𝐹 − 𝐷 2 𝑢 2 ( 𝑙 ) − 1 𝜕 𝑡 𝑢 2 ( 𝑙 ) 𝑢 − 𝐵 2 ( 𝑙 ) 𝑢 2 ( 𝑙 ) , 𝐹 w i t h 𝐷 2 𝑢 2 ( 𝑙 ) 𝐵 𝑢 = − 2 ( 𝑙 ) + 𝑢 𝜕 𝐵 2 ( 𝑙 ) 𝜕 𝑢 2 ( 𝑙 ) 𝑢 2 ( 𝑙 ) , 𝑢 2 ( 𝑙 ) 𝑡 𝑛 = 𝑢 𝐾 1 𝑡 𝑛 + 1 , 𝑙 = 0 , 1 , 2 , … , 𝐿 , ( 6 . 8 ) where we discretize the equations and obtain the discretized operators: 𝜕 𝑡 𝑢 1 ( 𝑘 ) 𝑢 − 𝐴 1 ( 𝑘 ) 𝑢 1 ( 𝑘 ) = 0 , ( 6 . 9 ) as 𝐹 1 𝑢 1 𝑡 𝑛 + 1 𝑢 1 ( 𝑘 ) 𝑡 𝑛 + 1 − 𝑢 1 𝑡 𝑛 𝑢 − Δ 𝑡 𝐴 1 ( 𝑘 ) 𝑡 𝑛 + 1 𝑢 1 ( 𝑘 ) 𝑡 𝑛 + 1 = 0 , ( 6 . 1 0 ) where we have the initialization of the Newton's method as 𝑢 1 ( 0 ) ( 𝑡 𝑛 + 1 ) = 0 or 𝑢 1 ( 0 ) ( 𝑡 𝑛 + 1 ) = 𝑢 1 ( 𝑡 𝑛 ) .

For the second iteration equation we have 𝜕 𝑡 𝑢 2 ( 𝑙 ) 𝑢 − 𝐵 2 ( 𝑙 ) 𝑢 2 ( 𝑙 ) , ( 6 . 1 1 ) as 𝐹 2 𝑢 2 𝑡 𝑛 + 1 = 𝑢 2 ( 𝑙 ) 𝑡 𝑛 + 1 − 𝑢 2 𝑡 𝑛 𝑢 − Δ 𝑡 𝐵 2 ( 𝑙 ) 𝑡 𝑛 + 1 𝑢 2 ( 𝑙 ) 𝑡 𝑛 + 1 = 0 , ( 6 . 1 2 ) where we have the initialization of the Newton's method as 𝑢 2 ( 0 ) ( 𝑡 𝑛 + 1 ) = 0 or 𝑢 2 ( 0 ) ( 𝑡 𝑛 + 1 ) = 𝑢 1 ( 𝑡 𝑛 ) .

The derivations are given as: 𝐷 𝐹 1 𝑢 1 𝑡 𝑛 + 1 𝐴 𝑢 = 1 − Δ 𝑡 1 𝑡 𝑛 + 1 + 𝑢 𝜕 𝐴 1 𝑡 𝑛 + 1 𝜕 𝑢 1 𝑡 𝑛 + 1 𝑢 1 𝑡 𝑛 + 1 , 𝐷 𝐹 2 𝑢 2 𝐵 𝑢 = 1 − Δ 𝑡 2 𝑡 𝑛 + 1 + 𝑢 𝜕 𝐵 2 𝑡 𝑛 + 1 𝜕 𝑢 2 𝑡 𝑛 + 1 𝑢 2 𝑡 𝑛 + 1 . ( 6 . 1 3 )

(3) The standard iterative operator-splitting method is given as follows.

We apply the quasilinear iterative operator-splitting method: 𝑑 𝑢 𝑖 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐴 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 𝑢 ( 𝑡 ) + 𝐵 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 − 1 ( 𝑡 ) , w i t h 𝑢 𝑖 𝑡 𝑛 = 𝑢 𝑛 , 𝑑 𝑢 𝑖 + 1 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐴 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 𝑢 ( 𝑡 ) + 𝐵 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 + 1 , w i t h 𝑢 𝑖 + 1 𝑡 𝑛 = 𝑢 𝑛 , ( 6 . 1 4 ) with the nonlinear operators 𝐴 ( 𝑢 ) 𝑢 = 𝜆 1 𝑢 ( 𝑡 ) + 𝜆 2 ( 𝑢 ( 𝑡 ) ) 𝑝 − 1 𝑢 , 𝐵 ( 𝑢 ) 𝑢 = 𝜆 3 𝑢 ( 𝑡 ) + 𝜆 4 ( 𝑢 ( 𝑡 ) ) 𝑝 − 1 𝑢 . The initialization of the fixpoint iteration is 𝑢 0 = 𝑢 𝑛 or 𝑢 0 = 0 with 𝐴 ( 𝑢 0 ) = 𝜆 1 and 𝐵 ( 𝑢 0 ) = 𝜆 3 .

For the iterations we can apply the analytical solution of each equation: 𝑢 𝑖 ( 𝑡 ) = 𝑢 𝑛 𝐴 𝑢 e x p 𝑖 − 1 𝑡 + 𝐴 𝑢 ( 𝑡 ) 𝑖 − 1 ( 𝑡 ) − 1 𝐵 𝑢 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 − 1 𝐴 𝑢 ( 𝑡 ) 1 − e x p 𝑖 − 1 𝑡 , 𝑢 ( 𝑡 ) 𝑖 + 1 ( 𝑡 ) = 𝑢 𝑛 𝐵 𝑢 e x p 𝑖 − 1 𝑡 + 𝐵 𝑢 ( 𝑡 ) 𝑖 − 1 ( 𝑡 ) − 1 𝐴 ( 𝑢 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 − 1 𝐵 𝑢 ( 𝑡 ) 1 − e x p 𝑖 − 1 𝑡 . ( 𝑡 ) ( 6 . 1 5 )

Further the iterative steps can be done.

(4) The Newton iterative method with embedded iterative operator-splitting method is given as follows.

We apply the quasilinear iterative operator-splitting method: 𝑑 𝑢 𝑖 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐴 𝑖 𝑢 ( 𝑡 ) 𝑖 𝑢 ( 𝑡 ) + 𝐵 𝑖 − 1 𝑢 ( 𝑡 ) 𝑖 − 1 ( 𝑡 ) , w i t h 𝑢 𝑖 𝑡 𝑛 = 𝑢 𝑛 , 𝑑 𝑢 𝑖 + 1 ( 𝑡 ) 𝑢 𝑑 𝑡 = 𝐴 𝑖 𝑢 ( 𝑡 ) 𝑖 𝑢 ( 𝑡 ) + 𝐵 𝑖 + 1 𝑢 ( 𝑡 ) 𝑖 + 1 , w i t h 𝑢 𝑖 + 1 𝑡 𝑛 = 𝑢 𝑛 , ( 6 . 1 6 ) with the nonlinear operators 𝐴 ( 𝑢 ) 𝑢 = 𝜆 1 𝑢 ( 𝑡 ) + 𝜆 2 ( 𝑢 ( 𝑡 ) ) 𝑝 − 1 𝑢 , 𝐵 ( 𝑢 ) 𝑢 = 𝜆 3 𝑢 ( 𝑡 ) + 𝜆 4 ( 𝑢 ( 𝑡 ) ) 𝑝 − 1 𝑢 . The initialization of the fixpoint iteration is 𝑢 0 ( 𝑡 𝑛 + 1 ) = 𝑢 𝑛 or 𝑢 0 ( 𝑡 𝑛 + 1 ) = 0 .

The discretization of the nonlinear ordinary differential equation is performed with higher-order Runge-Kutta methods.

The Newton method is applied as: 𝑢 𝑖 ( 𝑘 + 1 ) = 𝑢 𝑖 ( 𝑘 ) 𝐹 − 𝐷 1 𝑢 𝑖 ( 𝑘 ) − 1 𝜕 𝑡 𝑢 𝑖 ( 𝑘 ) 𝑢 − 𝐴 𝑖 ( 𝑘 ) 𝑢 𝑖 ( 𝑘 ) 𝑢 − 𝐵 𝑖 − 1 𝑢 𝑖 − 1 , 𝐷 𝐹 1 𝑢 𝑖 ( 𝑘 ) 𝐴 𝑢 = − 𝑖 ( 𝑘 ) + 𝑢 𝜕 𝐴 𝑖 ( 𝑘 ) 𝜕 𝑢 𝑖 ( 𝑘 ) 𝑢 𝑖 ( 𝑘 ) 𝑢 , 𝑘 = 0 , 1 , 2 , … , 𝐾 , 𝑖 𝑡 𝑛 = 𝑐 𝑛 , 𝑢 𝑖 𝑡 𝑛 + 1 = 𝑢 𝑖 𝑡 𝑛 + 1 𝐾 + 1 , w h e r e | 𝑢 𝑖 𝑡 𝑛 + 1 𝐾 + 1 − 𝑢 𝑖 𝑡 𝑛 + 1 𝐾 𝑢 | ≤ e r r , ( 𝑙 + 1 ) 𝑖 + 1 = 𝑢 ( 𝑙 ) 𝑖 + 1 𝐹 − 𝐷 2 𝑢 ( 𝑙 ) 𝑖 + 1 − 1 𝜕 𝑡 𝑢 ( 𝑙 ) 𝑖 + 1 𝑢 − 𝐴 𝑖 𝑢 𝑖 𝑢 − 𝐵 ( 𝑙 ) 𝑖 + 1 𝑢 ( 𝑙 ) 𝑖 + 1 𝑐 2 ( 𝑙 ) , 𝐷 𝐹 2 𝑢 ( 𝑙 ) 𝑖 + 1 𝐵 𝑢 = − ( 𝑙 ) 𝑖 + 1 + 𝑢 𝜕 𝐵 ( 𝑙 ) 𝑖 + 1 𝜕 𝑢 ( 𝑙 ) 𝑖 + 1 𝑢 ( 𝑙 ) 𝑖 + 1 𝑢 , 𝑙 = 0 , 1 , 2 , … , 𝐿 , 𝑖 + 1 𝑡 𝑛 = 𝑐 𝑛 , 𝑢 𝑖 + 1 𝑡 𝑛 + 1 = 𝑢 𝑖 + 1 𝑡 𝑛 + 1 𝐿 + 1 , w h e r e | 𝑢 𝑖 + 1 𝑡 𝑛 + 1 𝐿 + 1 − 𝑢 𝑖 + 1 𝑡 𝑛 + 1 𝐿 | ≤ e r r . ( 6 . 1 7 ) Here the time-step is 𝜏 = 𝑡 𝑛 + 1 − 𝑡 𝑛 . The iterations are 𝑖 = 1 , 3 , … , 2 𝑚 + 1 . 𝑢 0 ( 𝑡 ) = 0 is the starting solution and 𝑐 𝑛 is the known split approximation at the time-level 𝑡 = 𝑡 𝑛 . The results of the methods are 𝑢 ( 𝑡 𝑛 + 1 ) = 𝑢 2 𝑚 + 2 ( 𝑡 𝑛 + 1 ) .

We apply the discretization methods for the iteration steps.

We discretize the equations 𝜕 𝑡 𝑢 𝑖 ( 𝑘 ) 𝑢 − 𝐴 𝑖 ( 𝑘 ) 𝑢 𝑖 ( 𝑘 ) 𝑢 − 𝐵 𝑖 − 1 𝑢 𝑖 − 1 = 0 , ( 6 . 1 8 ) as 𝐹 1 𝑢 𝑖 ( 𝑘 ) 𝑡 𝑛 + 1 = 𝑢 𝑖 ( 𝑘 ) 𝑡 𝑛 + 1 − 𝑢 𝑖 𝑡 𝑛 𝐴 𝑢 − Δ 𝑡 𝑖 ( 𝑘 ) 𝑡 𝑛 + 1 𝑢 𝑖 ( 𝑘 ) 𝑡 𝑛 + 1 𝑢 + 𝐵 𝑖 − 1 𝑡 𝑛 + 1 𝑢 𝑖 − 1 𝑡 𝑛 + 1 , ( 6 . 1 9 ) where we have the initialization of the Newton's method as 𝑢 𝑖 ( 0 ) ( 𝑡 𝑛 + 1 ) = 0 or 𝑢 𝑖 ( 0 ) ( 𝑡 𝑛 + 1 ) = 𝑢 ( 𝑡 𝑛 ) .

For the second iteration equation we have: 𝜕 𝑡 𝑢 2 ( 𝑙 ) 𝑢 − 𝐴 𝑖 𝑢 𝑖 𝑢 − 𝐵 2 ( 𝑙 ) 𝑢 2 ( 𝑙 ) = 0 , ( 6 . 2 0 ) as 𝐹 2 𝑢 ( 𝑙 ) 𝑖 + 1 𝑡 𝑛 + 1 = 𝑢 ( 𝑙 ) 𝑖 + 1 𝑡 𝑛 + 1 − 𝑢 𝑖 + 1 𝑡 𝑛 𝐴 𝑢 − Δ 𝑡 𝑖 𝑡 𝑛 + 1 𝑢 𝑖 𝑡 𝑛 + 1 𝑢 + 𝐵 ( 𝑙 ) 𝑖 + 1 𝑡 𝑛 + 1 𝑢 ( 𝑙 ) 𝑖 + 1 𝑡 𝑛 + 1 , ( 6 . 2 1 ) where we have the initialization of the Newton's method as 𝑢 ( 0 ) 𝑖 + 1 ( 𝑡 𝑛 + 1 ) = 0 or 𝑢 ( 0 ) 𝑖 + 1 ( 𝑡 𝑛 + 1 ) = 𝑢 ( 𝑡 𝑛 ) .

The derivations are given as: 𝐷 𝐹 1 𝑢 𝑖 ( 𝑘 ) 𝑡 𝑛 + 1 𝐴 𝑢 = 1 − Δ 𝑡 ( 𝑘 ) 𝑖 + 1 𝑡 𝑛 + 1 + 𝑢 𝜕 𝐴 ( 𝑘 ) 𝑖 + 1 𝑡 𝑛 + 1 𝜕 𝑢 ( 𝑘 ) 𝑖 + 1 𝑡 𝑛 + 1 𝑢 ( 𝑘 ) 𝑖 + 1 𝑡 𝑛 + 1 , 𝐷 𝐹 2 𝑢 ( 𝑙 ) 𝑖 + 1 𝑡 𝑛 + 1 𝐵 𝑢 = 1 − Δ 𝑡 ( 𝑙 ) 𝑖 + 1 𝑡 𝑛 + 1 + 𝑢 𝜕 𝐵 ( 𝑙 ) 𝑖 + 1 𝑡 𝑛 + 1 𝜕 𝑢 ( 𝑙 ) 𝑖 + 1 𝑡 𝑛 + 1 𝑢 ( 𝑙 ) 𝑖 + 1 𝑡 𝑛 + 1 . ( 6 . 2 2 )

Our numerical results for the different methods are presented in Tables 1, 2, 3, and 4. The errors of the methods are shown in Figures 1, 2, and 3. We chose different iteration steps and time partitions. The error between the analytical and numerical solution is shown with the supremum norm at time 𝑇 = 1 . 0 .

tab1
Table 1: Numerical results for the Bernoulli equation with sequential operator-splitting method.
tab2
Table 2: Numerical results for the Bernoulli equation with sequential operator-splitting method with embedded Newton's method.
tab3
Table 3: Numerical results for the Bernoulli equation with iterative operator-splitting method.
tab4
Table 4: Numerical results for the Bernoulli equation with iterative operator-splitting method with embedded Newton's method.
307298.fig.001
Figure 1: Analytical and approximated solution with sequential operator-splitting method.
307298.fig.002
Figure 2: Analytical and approximated solution with sequential operator-splitting method with embedded Newton's method.
307298.fig.003
Figure 3: Analytical and approximated solution with iterative operator-splitting method.

The experiments show the reduced errors for more iteration steps and more time partitions. Because of the time-discretization method for ODEs, we restrict the number of iteration steps to a maximum of five. If we restrict the error bound to 1 0 − 3 , two iteration steps and five time partitions give the most effective combination.

6.2. Second Numerical Example: Mixed Convection-Diffusion and Burgers Equation

We deal with a 2D example which is a mixture of a convection-diffusion and Burgers equation. We can derive an analytical solution: 𝜕 𝑡 1 𝑢 = − 2 𝑢 𝜕 𝑥 1 𝑢 − 2 𝑢 𝜕 𝑦 1 𝑢 − 2 𝜕 𝑥 1 𝑢 − 2 𝜕 𝑦 𝑢 𝜕 + 𝜇 𝑥 𝑥 𝑢 + 𝜕 𝑦 𝑦 𝑢 + 𝑓 ( 𝑥 , 𝑦 , 𝑡 ) , ( 𝑥 , 𝑦 , 𝑡 ) ∈ Ω × [ 0 , 𝑇 ] , 𝑢 ( 𝑥 , 𝑦 , 0 ) = 𝑢 a n a ( 𝑥 , 𝑦 , 0 ) , ( 𝑥 , 𝑦 ) ∈ Ω , 𝑢 ( 𝑥 , 𝑦 , 𝑡 ) = 𝑢 a n a ( 𝑥 , 𝑦 , 𝑡 ) o n 𝜕 Ω × [ 0 , 𝑇 ] , ( 6 . 2 3 ) where Ω = [ 0 , 1 ] × [ 0 , 1 ] , 𝑇 = 1 . 2 5 , and 𝜇 is the viscosity.

The analytical solution is given as 𝑢 a n a ( 𝑥 , 𝑦 , 𝑡 ) = 1 + e x p 𝑥 + 𝑦 − 𝑡 2 𝜇 − 1 + e x p 𝑥 + 𝑦 − 𝑡 2 𝜇 , ( 6 . 2 4 ) where we compute 𝑓 ( 𝑥 , 𝑦 , 𝑡 ) accordingly.

We split the convection-diffusion and the Burgers equation. The operators are given as: 1 𝐴 ( 𝑢 ) 𝑢 = − 2 𝑢 𝜕 𝑥 1 𝑢 − 2 𝑢 𝜕 𝑦 1 𝑢 + 2 𝜇 𝜕 𝑥 𝑥 𝑢 + 𝜕 𝑦 𝑦 𝑢 , ( 6 . 2 5 ) hence 1 𝐴 ( 𝑢 ) = 2 − 𝑢 𝜕 𝑥 − 𝑢 𝜕 𝑦 𝜕 + 𝜇 𝑥 𝑥 + 𝜕 𝑦 𝑦 1 ( t h e B u r g e r s t e r m ) , 𝐵 𝑢 = − 2 𝜕 𝑥 1 𝑢 − 2 𝜕 𝑦 1 𝑢 + 2 𝜇 𝜕 𝑥 𝑥 𝑢 + 𝜕 𝑦 𝑦 𝑢 + 𝑓 ( 𝑥 , 𝑦 , 𝑡 ) ( t h e c o n v e c t i o n - d i ff u s i o n t e r m ) . ( 6 . 2 6 )

For the first equation we apply the nonlinear Algorithm 5.1 and obtain 𝐴 𝑢 𝑖 − 1 𝑢 𝑖 1 = − 2 𝑢 𝑖 − 1 𝜕 𝑥 𝑢 𝑖 − 1 2 𝑢 𝑖 − 1 𝜕 𝑦 𝑢 𝑖 + 1 2 𝜇 𝜕 𝑥 𝑥 𝑢 𝑖 + 𝜕 𝑦 𝑦 𝑢 𝑖 , 𝐵 𝑢 𝑖 − 1 = 1 2 − 𝜕 𝑥 − 𝜕 𝑦 𝜕 + 𝜇 𝑥 𝑥 + 𝜕 𝑦 𝑦 𝑢 𝑖 − 1 , ( 6 . 2 7 ) and we obtain linear operators, because 𝑢 𝑖 − 1 is known from the previous time-step.

In the second equation we obtain by using Algorithm 5.1: 𝐴 𝑢 𝑖 − 1 𝑢 𝑖 1 = − 2 𝑢 𝑖 − 1 𝜕 𝑥 𝑢 𝑖 − 1 2 𝑢 𝑖 − 1 𝜕 𝑦 𝑢 𝑖 + 1 2 𝜇 𝜕 𝑥 𝑥 𝑢 𝑖 + 𝜕 𝑦 𝑦 𝑢 𝑖 , 𝐵 𝑢 𝑖 + 1 = 1 2 − 𝜕 𝑥 − 𝜕 𝑦 𝜕 + 𝜇 𝑥 𝑥 + 𝜕 𝑦 𝑦 𝑢 𝑖 + 1 , ( 6 . 2 8 ) and we have linear operators.

We deal with different viscosities 𝜇 as well as different step-sizes in time and space. We have the following results (see Tables 5 and 6).

tab5
Table 5: Numerical results for the mixed convection-diffusion and Burgers equation with viscosit 𝜇 = 0 . 5 , initial conditio 𝑢 0 ( 𝑡 ) = 𝑐 𝑛 , and four iterations per time-step.
tab6
Table 6: Numerical results for the mixed convection-diffusion and Burgers equation with viscosit 𝜇 = 5 , initial conditio 𝑢 0 ( 𝑡 ) = 𝑐 𝑛 , and two iterations per time-step.

Figure 4 presents the profile of the 2D linear and nonlinear convection-diffusion equation.

fig4
Figure 4: Mixed convection-diffusion and Burgers equation at initial time 𝑡 = 0 . 0 (a) and end time 𝑡 = 1 . 2 5 (b) for viscosity 𝜇 = 0 . 5 .

Remark 6.1. In the examples, we deal with more iteration steps to obtain higher-order convergence results. In the first test we have four iterative steps but a smaller viscosity ( 𝜇 = 0 . 5 ) such that we can reach at least a second-order method. In the second test we use a high viscosity about 𝜇 = 5 and get the second-order result with two iteration steps. Here we see the loss of differentiability, that becomes stiff equation parts. To obtain the same results, we have to increase the number of iteration steps. Therefore we can show an improvement of the convergence order with respect to the iteration steps.

6.3. Third Numerical Example: Momentum Equation (Molecular Flow)

We deal with an example of a momentum equation, that is used to model the viscous flow of a fluid. 𝜕 𝑡 1 𝐮 = − 𝐮 ⋅ ∇ 𝐮 + 2 𝜇 ∇ 𝐷 ( 𝐮 ) + 3 ∇ 𝐮 + 𝐟 ( 𝑥 , 𝑦 , 𝑡 ) , ( 𝑥 , 𝑦 , 𝑡 ) ∈ Ω × [ 0 , 𝑇 ] , 𝐮 ( 𝑥 , 𝑦 , 0 ) = 𝐠 1 ( 𝑥 , 𝑦 ) , ( 𝑥 , 𝑦 ) ∈ Ω , 𝐮 ( 𝑥 , 𝑦 , 𝑡 ) = 𝐠 2 ( 𝑥 , 𝑦 , 𝑡 ) o n 𝜕 Ω × [ 0 , 𝑇 ] , ( e n c l o s e d fl o w ) , ( 6 . 2 9 ) where 𝐮 = ( 𝑢 1 , 𝑢 2 ) 𝑡 is the solution and Ω = [ 0 , 1 ] × [ 0 , 1 ] , 𝑇 = 1 . 2 5 , 𝜇 = 5 , and 𝐯 = ( 0 . 0 0 1 , 0 . 0 0 1 ) 𝑡 are the parameters and 𝐼 is the unit matrix.

The nonlinear function 𝐷 ( 𝐮 ) = 𝐮 ⋅ 𝐮 + 𝐯 ⋅ 𝐮 is the viscosity flow, and 𝐯 is a constant velocity.

We can derive the analytical solution with respect to the first two test examples with the functions: 𝑢 1 ( 𝑥 , 𝑦 , 𝑡 ) = 1 + e x p 𝑥 + 𝑦 − 𝑡 2 𝜇 − 1 + e x p 𝑥 + 𝑦 − 𝑡 , 𝑢 2 𝜇 2 ( 𝑥 , 𝑦 , 𝑡 ) = 1 + e x p 𝑥 + 𝑦 − 𝑡 2 𝜇 − 1 + e x p 𝑥 + 𝑦 − 𝑡 . 2 𝜇 ( 6 . 3 0 ) For the splitting method our operators are given as: 2 𝐴 ( 𝐮 ) 𝐮 = − 𝐮 ∇ 𝐮 + 2 𝜇 ∇ 𝐷 ( 𝐮 ) ( t h e n o n l i n e a r o p e r a t o r ) , 𝐵 𝐮 = 3 𝜇 Δ 𝐮 ( t h e l i n e a r o p e r a t o r ) . ( 6 . 3 1 )

We deal first with the one-dimensional case, 𝜕 𝑡 𝑢 = − 𝑢 ⋅ 𝜕 𝑥 𝑢 + 2 𝜇 𝜕 𝑥 1 𝐷 ( 𝑢 ) + 3 𝜕 𝑥 𝑢 + 𝑓 ( 𝑥 , 𝑡 ) , ( 𝑥 , 𝑡 ) ∈ Ω × [ 0 , 𝑇 ] , 𝑢 ( 𝑥 , 0 ) = 𝑔 1 ( 𝑥 ) , ( 𝑥 ) ∈ Ω , 𝑢 ( 𝑥 , 𝑡 ) = 𝑔 2 ( 𝑥 , 𝑡 ) o n 𝜕 Ω × [ 0 , 𝑇 ] , ( e n c l o s e d fl o w ) , ( 6 . 3 2 ) where 𝑢 is the solution and Ω = [ 0 , 1 ] , 𝑇 = 1 . 2 5 , 𝜇 = 5 , and 𝑣 = 0 . 0 0 1 are the parameters.

Then the operators are given as: 𝐴 ( 𝑢 ) 𝑢 = − 𝑢 𝜕 𝑥 𝑢 + 2 𝜇 𝜕 𝑥 2 𝐷 ( 𝑢 ) ( t h e n o n l i n e a r o p e r a t o r ) , 𝐵 𝑢 = 3 𝜇 𝜕 𝑥 𝑥 𝑢 ( t h e l i n e a r o p e r a t o r ) . ( 6 . 3 3 ) For the iterative operator-splitting as fixed point scheme, we have the following results (see Tables 7 and 9).

tab7
Table 7: Numerical results for the 1D momentum equation wit 𝜇 = 5 𝑣 = 0 . 0 0 1 , initial conditio 𝑢 0 ( 𝑡 ) = 𝑐 𝑛 , and two iterations per time-step.

Figure 5 presents the profile of the 1D momentum equation.

fig5
Figure 5: 1D momentum equation at initial time 𝑡 = 0 . 0 (a) and end time 𝑡 = 1 . 2 5 (b) for 𝜇 = 5 and 𝑣 = 0 . 0 0 1 .

We have the following results for the 2D case (see Tables 10, 11, and 12).

Figure 6 presents the profile of the 2D momentum equation.

fig6
Figure 6: 2D momentum equation at initial time 𝑡 = 0 . 0 (a), (c)) and end time 𝑡 = 1 . 2 5 (b), (d)) for 𝜇 = 0 . 5 and 𝑣 = ( 1 , 1 ) 𝑡 for the first and second component of the numerical solution.

For the Newton operator-splitting method we obtain the following functional matrices for the one-dimensional case: 𝐷 𝐹 ( 𝑢 ) = ( 4 𝜇 − 1 ) 𝜕 𝑥 𝜕 𝑢 , 𝐷 ( 𝐹 ( 𝑢 ) ) = − 𝑢 1 𝐹 1 ( 𝑢 ) 𝜕 𝑢 2 𝐹 1 𝜕 ( 𝑢 ) 𝑢 1 𝐹 2 ( 𝑢 ) 𝜕 𝑢 2 𝐹 2 ( 𝑢 ) = − − 𝜕 𝑥 𝑢 1 + 4 𝜇 𝜕 𝑥 𝑢 1 − 𝜕 𝑥 𝑢 2 + 4 𝜇 𝜕 𝑥 𝑢 2 − 𝜕 𝑦 𝑢 1 + 4 𝜇 𝜕 𝑦 𝑢 1 − 𝜕 𝑦 𝑢 2 + 4 𝜇 𝜕 𝑦 𝑢 2 = ( 4 𝜇 − 1 ) ∇ 𝑢 . ( 6 . 3 4 ) For the two-dimensional case, we use: 𝑢 𝐴 ( 𝐮 ) 𝐮 = − 𝐮 ∇ 𝐮 + 2 𝜇 ∇ 𝐷 ( 𝐮 ) = − 1 𝜕 𝑥 𝑢 1 + 𝑢 2 𝜕 𝑥 𝑢 2 𝑢 1 𝜕 𝑦 𝑢 1 + 𝑢 2 𝜕 𝑦 𝑢 2 + 2 𝜇 2 𝑢 1 𝜕 𝑥 𝑢 1 + 2 𝑢 2 𝜕 𝑥 𝑢 2 + 𝑣 1 𝜕 𝑥 𝑢 1 + 𝑣 2 𝜕 𝑥 𝑢 2 2 𝑢 1 𝜕 𝑦 𝑢 1 + 2 𝑢 2 𝜕 𝑦 𝑢 2 + 𝑣 1 𝜕 𝑦 𝑢 1 + 𝑣 2 𝜕 𝑦 𝑢 2 . ( 6 . 3 5 )

Here, we do not need the linearization and apply the standard iterative splitting method.

We only linearize the first split step and therefore we can relax this step with the second linear split step. Therefore we obtain stable methods, see [15]. For the Newton iterative method, we have the following results, see Table 8.

tab8
Table 8: Numerical results for the 1D momentum equation wit 𝜇 = 5 𝑣 = 0 . 0 0 1 , initial conditio 𝑢 0 ( 𝑡 ) = 𝑐 𝑛 , two iterations per time-step an 𝐾 = 1 using Newton iterative method.
tab9
Table 9: Numerical results for the 1D momentum equation wit 𝜇 = 5 0 𝑣 = 0 . 1 , initial conditio 𝑢 0 ( 𝑡 ) = 𝑐 𝑛 , and two iterations per time-step.
tab10
Table 10: Numerical results for the 2D momentum equation with 𝜇 = 2 , 𝑣 = ( 1 , 1 ) 𝑡 , initial condition 𝑢 0 ( 𝑡 ) = 𝑐 𝑛 , and two iterations per time-step.
tab11
Table 11: Numerical results for the 2D momentum equation for the first component with 𝜇 = 5 0 , 𝑣 = ( 1 0 0 , 0 . 0 1 ) 𝑡 , initial condition 𝑢 0 ( 𝑡 ) = 𝑐 𝑛 , and two iterations per time-step.
tab12
Table 12: Numerical results for the 2D momentum equation for the second component with 𝜇 = 5 0 , 𝑣 = ( 1 0 0 , 0 . 0 1 ) 𝑡 , initial condition 𝑢 0 ( 𝑡 ) = 𝑐 𝑛 , and two iterations per time-step.

Remark 6.2. In the more realistic examples of 1D and 2D momentum equations, we can also observe the stiffness problem, which we obtain with a more hyperbolic behavior. In the 1D experiments we deal with a more hyperbolic behavior and can obtain at least first-order convergence with two iterative steps. In the 2D experiments we obtain nearly second-order convergence results with two iterative steps, if we increase the parabolic behavior, for example, larger 𝜇 and 𝐯 values. For such methods, we have to balance the usage of the iterative steps with refinement in time and space with respect to the hyperbolicity of the equations. At least we can obtain a second-order method with more than two iterative steps. Therefore the stiffness influences the number of iterative steps.

7. Conclusion and Discussion

We present decomposition methods for differential equations based on iterative and non-iterative methods. The nonlinear equations are solved with embedded Newton's methods. We present new ideas on linearization to obtain more accurate results. The superiority of the new embedded Newton's methods over the traditional sequential methods is demonstrated in examples, especially through their simple implementation. Further, we have the smoothing properties of the iterative scheme that allow a balance between the nonlinear and the linear terms. The results show more accurate solutions with respect to time decomposition. In the future the iterative operator-splitting method can be generalized for multi-dimensional problems and also for non-smooth and nonlinear problems in time and space. In next paper we discuss error analysis of nonlinear methods.

References

  1. G. I. Marchuk, “Some applicatons of splitting-up methods to the solution of problems in mathematical physics,” Aplikace Matematiky, vol. 1, pp. 103–132, 1968.
  2. G. Strang, “On the construction and comparison of difference schemes,” SIAM Journal on Numerical Analysis, vol. 5, no. 3, pp. 506–517, 1968. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  3. H. Yoshida, “Construction of higher order symplectic integrators,” Physics Letters A, vol. 150, no. 5–7, pp. 262–268, 1990. View at Publisher · View at Google Scholar · View at MathSciNet
  4. J. F. Kanney, C. T. Miller, and C. T. Kelley, “Convergence of iterative split-operator approaches for approximating nonlinear reactive problems,” Advances in Water Resources, vol. 26, no. 3, pp. 247–261, 2003. View at Publisher · View at Google Scholar
  5. P. Deuflhard, Newton Methods for Nonlinear Problems: Affine Invariance and Adaptive Algorithm, vol. 35 of Springer Series in Computational Mathematics, Springer, Berlin, Germany, 2004. View at Zentralblatt MATH · View at MathSciNet
  6. P. Deuflhard, “A modified Newton method for the solution of ill-conditioned systems of nonlinear equations with application to multiple shooting,” Numerische Mathematik, vol. 22, no. 4, pp. 289–315, 1974. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  7. M. Ruzicka, Nichtlineare Funktionalanalysis, Springer, Berlin, Germany, 2004. View at Zentralblatt MATH
  8. E. Zeidler, Nonlinear Functional Analysis and Its Applications. II/B: Nonlinear Monotone Operators, Springer, New York, NY, USA, 1990. View at Zentralblatt MATH · View at MathSciNet
  9. J. Geiser, O. Klein, and P. Philip, “Anisotropic thermal conductivity in apparatus insulation: numerical study of effects on the temperature field during sublimation growth of silicon carbide single crystals,” preprint no. 1034 of Weierstrass Institute for Applied Analysis and Stochastics, Berlin, Germany, 2005.
  10. J. Geiser, “Discretization methods with embedded analytical solutions for convection dominated transport in porous media,” in Proceedings of the 3rd International Conference on Numerical Analysis and Its Applications (NAA '04), vol. 3401 of Lecture Notes in Computer Science, pp. 288–295, Springer, Rousse, Bulgaria, June-July 2004.
  11. W. Hundsdorfer and J. Verwer, Numerical Solution of Time-Dependent Advection-Diffusion-Reaction Equations, vol. 33 of Springer Series in Computational Mathematics, Springer, Berlin, Germany, 2003. View at Zentralblatt MATH · View at MathSciNet
  12. K. H. Karlsen, K.-A. Lie, J. R. Natvig, H. F. Nordhaug, and H. K. Dahle, “Operator splitting methods for systems of convection-diffusion equations: nonlinear error mechanisms and correction strategies,” Journal of Computational Physics, vol. 173, no. 2, pp. 636–663, 2001. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  13. K. H. Karlsen and N. H. Risebro, “Corrected operator splitting for nonlinear parabolic equations,” SIAM Journal on Numerical Analysis, vol. 37, no. 3, pp. 980–1003, 2000. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  14. I. Faragó and J. Geiser, “Iterative operator-splitting methods for linear problems,” International Journal of Computational Science and Engineering, vol. 3, no. 4, pp. 255–263, 2007. View at Publisher · View at Google Scholar
  15. C. T. Kelly, Iterative Methods for Linear and Nonlinear Equations, Frontiers in Applied Mathematics, SIAM, Philadelphia, Pa, USA, 1995. View at Zentralblatt MATH · View at MathSciNet
  16. J. Geiser, “Iterative operator-splitting methods with higher-order time integration methods and applications for parabolic partial differential equations,” Journal of Computational and Applied Mathematics, vol. 217, no. 1, pp. 227–242, 2008. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  17. K. H. Karlsen and N. H. Risebro, “An operator splitting method for nonlinear convection-diffusion equations,” Numerische Mathematik, vol. 77, no. 3, pp. 365–382, 1997. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet