1 Introduction
The celebrated LaxRichtmeyer equivalence theorem ([21], [12], [29]) states that if the numerical scheme is stable then its global error is at least of the same order as its local truncation error. Indeed, it is the common experience that the frequently used methods produce solutions that have global error of the same order as the normalized local truncation error. Indeed, the Dahlquist’s Equivalence Theorem [35] states this fact for linear multistep methods. In fact, this behavior is so expected for all onestep and multistep methods that we typically define the order of a stable numerical method solely by the order conditions derived by Taylor series analysis of the local truncation error. This relationship between the order of the local truncation error and global error is also seen in finite difference schemes for partial differential equations (PDEs) [12, 29]. Work over the past decade reminds us that, while for a stable scheme the global error must at least of the same order as its local truncation error, it may in fact be of higher order.
In this paper we consider numerical solvers for ordinary differential equations (ODEs) of the form
(where we assume is at least continuous with respect to and Lipschitzcontinuous with respect to , and typically much smoother). Without loss of generality [14], we focus our attention on the autonomous ODE
(1)  
The canonical numerical method for this problem is the forward Euler method
where approximates the exact solution at time and we let . The forward Euler method has local truncation error and its approximation error at any given time defined by
and it produces a global error which is first order accurate
To develop methods that have higher order, we can add steps and define a linear multistep method [3], use multiple stages as in a Runge–Kutta methods [3], or include higher derivatives as we do in Taylor series methods [3]. We can also combine the approaches above: the combination of multiple steps and stages is commonly see as in the general linear methods described in [2, 14], and multiple derivatives and stages have also been used [4, 5, 10, 23, 24, 27, 30]. In all these cases the goal is to increase the order of the local truncation error and therefore of the global error. All these approaches typically produce methods that have global errors that are of the same order as their local truncation errors.
It is, however, possible in some cases to devise schemes that have global errors that are higher order than predicted by the local truncation errors. Kulikov (for Nordsiek methods) [18] and by Weiner and colleagues [40] (for explicit twostep peer methods) followed along the lines of the quasiconsistency theory first introduced by Skeel in 1978 [31] to find methods that have solution of order although their truncation errors are of order . In [6] Ditkowski and Gottlieb derived sufficient conditions for a family of general linear methods (GLMs) of any number of steps and stages under which one can control the accumulation of the local truncation error over time evolution and produce methods that are one order higher than expected from the truncation error alone. A similar theory was developed for implicitexplicit methods in [34]. In [19], [39], and [20]
, it was shown that under additional conditions on the method the precise form of the first surviving term in the global error (the vector multiplying
) can be computed explicitly and leveraged for error estimation. In
[7] Ditkowski, Gottlieb, and Grant showed that under less restrictive conditions the form of the first surviving term in the global error can be computed explicitly and exploited to enhance the accuracy to order , using postprocessing.In this paper, we extend the error inhibiting approach to methods that have multiple steps, multiple stages, and two derivatives. In Section 3.1 we show that the theory developed in [7] carries through easily to such methods, and by imposing the same EIS+ conditions as in [7] on on the twoderivative GLMs we consider we can not only produce solutions with accuracy one order higher than predicted by truncation error analysis but also precisely describe the leading coefficients of the error term. This allows us to postprocess the final time solution and obtain a global error that is two orders higher than the local truncation error. The construction of the postprocessor is exactly the same as in [7] and we review it in Section 3.2. We proceed to devise error inhibiting twostep GLMs this new EIS approach, present these methods in Section 4, and in Section 5 we test a selection of these methods on numerical examples to demonstrate their enhanced accuracy properties.
2 Twoderivative methods
Multiderivative Runge–Kutta methods, particularly twoderivative methods, were first considered in [25, 38, 36, 32, 33, 15, 16, 23, 27, 4], and later explored for use with partial differential equations (PDEs) [30, 37, 22, 28, 8]. Recent interest in exploring the strongstability properties of multiderivative Runge–Kutta methods, (also known as multistage multiderivative methods) is evident in [5, 10, 24]. A 2015 work by Okuonghae and Ikhile [26] discussed two and threederivatives methods with multiple stages and steps. These twoderivative general linear methods (GLMs) are similar to those in this work. (Also cite the paper Zack refereed and look in there for references)
2.1 Essentials
In this work we consider the class of twoderivative general linear methods (GLMs) of the form
(2) 
where is a vector of length that contains the numerical solution at times for :
(3) 
The function is defined as the componentwise function evaluation on the vector :
(4) 
and, similarly, the function is
(5) 
For convenience, we select so that the first element in the vector approximates the solution at time , and the abscissas are nondecreasing . To initialize these methods, we define the first element in the initial solution vector and the remaining elements are computed using some highly accurate method.
We define the corresponding projection of the exact solution of the ODE (1) onto the temporal grid:
(6) 
where and are the componentwise function evaluation on the vector . We define the global error as the difference between the vectors of the exact and the numerical solutions at some time
(7) 
Methods of the form (2) have a local truncation error at time that is related to the approximation error
(8) 
where
(9) 
where the truncation error vectors have the form equationparentequation
(10a)  
Here, is the vector of abscissas and is the vector of ones, and the terms are understood componentwise . For a method to have order it must satisfy the order conditions
Note that the form (2) includes both explicit and implicit schemes, as appears on both sides of the equation. However, if and are both strictly lower triangular the scheme is explicit. We are only interested in zerostable methods. A sufficient condition for this is that the coefficient matrix is a rank one matrix that has row sum one (i.e. satisfies the consistency condition). For simplicity we assume this to be the case in this the remainder of this work.
2.2 Preliminaries
In this subsection we make some observations that will be useful in the remainder of the paper. In all the following we assume sufficient smoothness of all the quantities such as etc. Furthermore, we assume that he order conditions hold for all , and that we are in the asymptotic regime so that the errors are small and we can say that .^{1}^{1}1We will verify that this is reasonable in Lemma 2.2
Observation 1 We observe that
(11) 
where .
Proof
This is simply due to the fact that is smooth enough that it can be expanded:
where the error vector is , and we use the notation . Each term can be expanded as
so that, exploiting the smoothness of we have
which, due to the fact that becomes
Observation 2 Given the smoothness of , we observe that
(12) 
where .
Proof
The proof of this observation is exactly the same as above, with the assumption that is smooth enough to be expanded.
We use these observations to develop an equation that describes the growth of the error:
Lemma
Given a zerostable method of the form (2) which satisfies the order conditions
and where the functions and are smooth, the evolution of the error can be described by:
(13) 
Proof
Recall that
and use this equation to subtract from
obtaining
using Observations 1 and 2 we have
so that
Note that and are all scalars. Now recall that since , and that terms of the form and are all bounded, so we have
Now move the terms to the left side:
Noting that we have the desired result
We now turn to verifying that the assumption we made in the first paragraph of this subsection, that we can assume that the error is small:
Lemma
The proof of this Lemma is given in [7], Lemma 3.
3 Designing multiderivative error inhibiting schemes that can be postprocessed to order
In this section we show how to construct twoderivative GLMs of the form (2) that are error inhibiting and so produce a solution of order even though one would only expect order from their truncation errors. We further show that under additional conditions we can express the exact form of the error and define an associated postprocessor that allow us to recover order from a scheme that would otherwise be only th order accurate.
3.1 Twoderivative error inhibiting schemes that have order
In this section we consider a multiderivative GLM of the form (2) where we assume that is a rank one matrix that satisfies the consistency condition so that this scheme is zerostable. Furthermore, we assume that the coefficient matrices , , , , are such that the order conditions are satisfied for all , so that the method will give us a numerical solution that has error that is guaranteed of order . In the following theorem we show that if the truncation error vector lives in the nullspace of the operator the order of the error can be shown to be of order . Furthermore, the theorem shows that under additional conditions on the coefficient matrices , , , , , we can determine precisely what the leading term of this error will look like and therefore remove it by postprocessing.
Theorem 1
Given a zerostable twoderivative general linear method of the form
where is a rank one matrix that satisfies the consistency condition , and the coefficient matrices , , , , satisfy the order conditions
where
if the error inhibiting condition equationparentequation
(14a)  
is satisfied, then the numerical solution produced by this method will have error  
Furthermore, if the conditions  
(14b)  
(14c) 
are also satisfied, then error vector will have the more precise form:
(15) 
Proof
Using Lemma 1 we obtain the equation for the evolution of the error
(16) 
the fact is smooth assures us that so that we can expand the first term to obtain
The discrete Duhamel principle (given as Lemma 5.1.1 in [12]) states that given an iterative process of the form
where is a linear operator, we have
(17) 
To analyze the error, we separate it into four parts:
and discuss each part separately:

The method is initialized so that the numerical solution vector is accurate enough to ensure that the initial error is negligible and we can ignore the first term.

The final term in the summation is . Recall that
so that this term contributes to the final time error the term

The term is a product of the operator
and the approximation error
so we obtain
due to the condition (14a) that states that .

Finally we look at the rest of the sum and use the boundedness of the operator to observe
Clearly, the first term here is a constant that depends only on the final time, so it is the product we need to bound. Using the definition of the operators
and the approximation error
we have
Using the fact that in our case and that we obtain
The first term and third terms disappear because of (14a) and (14b)
The second term is eliminated by (14c)
So that
All these parts together give us the result
(18)  
Remark
We note that generalizing this approach from our work in [7] to include a second derivative was quite simple. This is due to the fact that the second derivative appears with a term multiplying it, so that the operator is in fact exactly the same for a twoderivative method of the form (2) as for the methods considered in [7]. For this reason, expanding this approach to higher derivatives will be simple: the only change would be in the definition of the order conditions (10).
3.2 Designing a postprocesser to recover order
At every timestep the error vector has the particular form (15)
so that the leading order term can be removed by postprocessing at the end of the computation and a final time solution of order can be obtained, as we showed in [7]. We notice that the form of the error for multiderivative method (2) is exactly the same as that of the error of the EIS methods presented in [7], the only difference being the definition of the truncation error vector. The theoretical discussion of this postprocessor is presented in [7] and we refer the reader to the presentation there. We review one such construction below:

Select the number of computation steps that will be used, requiring that
where is the number of steps and the truncationerror order of the timestepping method.

Define the time vector of all the temporal grid points in the last computation steps:

Correspondingly, concatenate the final
Comments
There are no comments yet.