A note on error expressions for reflected and averaged implicit Runge-Kutta methods

In addition to their usefulness in the numerical solution of initial value ODE's, the implicit Runge-Kutta (IRK) methods are also important for the solution of two-point boundary value problems. Recently, several classes of modified IRK methods which improve significantly on the efficiency of the standard IRK methods in this application have been presented. One such class is the Averaged IRK methods; a member of the class is obtained by applying an averaging operation to a non-symmetric IRK method and its reflection. In this paper we investigate the forms of the error expressions for reflected and averaged IRK methods. Our first result relates the expression for the local error of the reflected method to that of the original method. The main result of this paper relates the error expression of an averaged method to that of the method upon which it is based. We apply these results to show that for each member of the class of the averaged methods, there exists an embedded lower order method which can be used for error estimation, in a formula-pair fashion.This work was supported by the Natural Science and Engineering Research Council of Canada.     

Runge-Kutta schemes that maintain invariant solutions for differential algebraic systems

Implicit Runge-Kutta (IRK) methods and projected IRK methods for the solution of semiexplicit index-2 systems of differential algebraic systems (DAEs) have been proposed by several authors. In this paper we prove that if a method satisfiesBA+A ^t B–bb ^t =0, it conserves quadratic invariants of DAEs.     

On implicit Runge-Kutta methods with high stage order

It is well known that high stage order is a desirable property for implicit Runge-Kutta methods. In this paper it is shown that it is always possible to construct ans-stage IRK method with a given stability function and stage orders−1 if the stability function is an approximation to the exponential function of at least orders. It is further indicated how to construct such methods as well as in which cases the constructed methods will be stiffly accurate.     

An improved approximate Newton method for implicit Runge-Kutta formulas

Implicit Runge-Kutta (IRK) methods (such as the s-stage Radau IIA method with s=3,5, or 7) for solving stiff ordinary differential equation systems have excellent stability properties and high solution accuracy orders, but their high computing costs in solving their nonlinear stage equations have seriously limited their applications to large scale problems. To reduce such a cost, several approximate Newton algorithms were developed, including a commonly used one called the simplified Newton method. In this paper, a new approximate Jacobian matrix and two new test rules for controlling the updating of approximate Jacobian matrices are proposed, yielding an improved approximate Newton method. Theoretical and numerical analysis show that the improved approximate Newton method can significantly improve the convergence and performance of the simplified Newton method.     

Butcher's simplifying assumption for symplectic integrators

Implicit Runge-Kutta methods with vanishingM matrix are discussed for preserving the symplectic structure of Hamiltonian systems. The number of the order conditions independent of the number of stages can be reduced considerably for the symplectic IRK method through the analysis utilizing the rooted tree and the corresponding elementary differentials. Butcher's simplifying condition further reduces the number of independent order conditions.     

A stability result for general linear methods with characteristic function having real poles only

It is shown that there exist A-stable multistep formulae, with a characteristic function havings poles, all of which are real, with orderp satisfyingp>s+1. This contradicts the widely held belief thatp=s+1 is the maximum possible order of such a method.     

Energy-preserving methods for Poisson systems

We present and analyze energy-conserving methods for the numerical integration of IVPs of Poisson type that are able to preserve some Casimirs. Their derivation and analysis is done following the ideas of Hamiltonian BVMs (HBVMs) (see Brugnano et al. [10] and references therein). It is seen that the proposed approach allows us to obtain the methods recently derived in Cohen and Hairer (2011) [17], giving an alternative derivation of such methods and a new proof of their order. Sufficient conditions that ensure the existence of a unique solution of the implicit equations defining the formulae are given. A study of the implementation of the methods is provided. In particular, order and preservation properties when the involved integrals are approximated by means of a quadrature formula, are derived.     

New embedded explicit pairs of exponentially fitted Runge-Kutta methods

Two new embedded pairs of exponentially fitted explicit Runge-Kutta methods with four and five stages for the numerical integration of initial value problems with oscillatory or periodic solutions are developed. In these methods, for a given fixed ω the coefficients of the formulae of the pair are selected so that they integrate exactly systems with solutions in the linear space generated by {sinh(ωt),cosh(ωt)}, the estimate of the local error behaves as O(h⁴) and the high-order formula has fourth-order accuracy when the stepsize h→0. These new pairs are compared with another one proposed by Franco [J.M. Franco, An embedded pair of exponentially fitted explicit Runge-Kutta methods, J. Comput. Appl. Math. 149 (2002) 407-414] on several problems to test the efficiency of the new methods.     

On the integration of stiff systems of O.D.E.s using extended backward differentiation formulae

Summary A class of extended backward differentiation formulae suitable for the approximate numerical integration of stiff systems of first order ordinary differential equations is derived. An algorithm is described whereby the required solution is predicted using a conventional backward differentiation scheme and then corrected using an extended backward differentiation scheme of higher order. This approach allows us to developL-stable schemes of order up to 4 andL()-stable schemes of order up to 9. An algorithm based on the integration formulae derived in this paper is illustrated by some numerical examples and it is shown that it is often superior to certain existing algorithms.     

On cubature formulae of degree 4k+1 attaining Möller's lower bound for integrals with circular symmetry

Summary The structure of cubature formulae of degree 4k+1 whose number of nodes is equal to Möller's lower bound is investigated for integrals with circular symmetry. A simple criterion is derived for the existence of such formulae. It shows that fork=1 Möller's lower bound can always be attained with Radon's formulae. It also allows to prove that for several integrals with circular symmetry and several values ofk>1, Möller's lower bound cannot be attained.     

A note on minimization problems and multistep methods

Summary. We use the qualitative properties of the solution flow of the gradient equation to compute a local minimum of a real-valued function . Under the regularity assumption of all equilibria we show a convergence result for bounded trajectories of a consistent, strictly stable linear multistep method applied to the gradient equation. Moreover, we compare the asymptotic features of the numerical and the exact solutions as done by Humphries, Stuart (1994) and Schropp (1995) for one-step methods. In the case of -stable formulae this leads to an efficient solver for stiff minimization problems. Received July 10, 1995 / Revised version received June 27, 1996     

High order methods for the numerical integration of ordinary differential equations

Summary High order implicit integration formulae with a large region of absolute stability are developed for the approximate numerical integration of both stiff and non-stiff systems of ordinary differential equations. The algorithms derived behave essentially like one step methods and are demonstrated by direct application to certain particular examples.     

Construction of Newton-like iteration methods for solving nonlinear equations

In this paper, we present a simple, and yet powerful and easily applicable scheme in constructing the Newton-like iteration formulae for the computation of the solutions of nonlinear equations. The new scheme is based on the homotopy analysis method applied to equations in general form equivalent to the nonlinear equations. It provides a tool to develop new Newton-like iteration methods or to improve the existing iteration methods which contains the well-known Newton iteration formula in logic; those all improve the Newton method. The orders of convergence and corresponding error equations of the obtained iteration formulae are derived analytically or with the help of Maple. Some numerical tests are given to support the theory developed in this paper.     

Global error estimation and control in linearly-implicit parallel two-step peer W-methods

The class of linearly-implicit parallel two-step peer W-methods has been designed recently for efficient numerical solutions of stiff ordinary differential equations. Those schemes allow for parallelism across the method, that is an important feature for implementation on modern computational devices. Most importantly, all stage values of those methods possess the same properties in terms of stability and accuracy of numerical integration. This property results in the fact that no order reduction occurs when they are applied to very stiff problems. In this paper, we develop parallel local and global error estimation schemes that allow the numerical solution to be computed for a user-supplied accuracy requirement in automatic mode. An algorithm of such global error control and other technical particulars are also discussed here. Numerical examples confirm efficiency of the presented error estimation and stepsize control algorithm on a number of test problems with known exact solutions, including nonstiff, stiff, very stiff and large-scale differential equations. A comparison with the well-known stiff solver RODAS is also shown.     

B-spline collocation for solution of two-point boundary value problems

A numerical method based on B-spline is developed to solve the general nonlinear two-point boundary value problems up to order 6. The standard formulation of sextic spline for the solution of boundary value problems leads to non-optimal approximations. In order to derive higher orders of accuracy, high order perturbations of the problem are generated and applied to construct the numerical algorithm. The error analysis and convergence properties of the method are studied via Green’s function approach. O(h⁶) global error estimates are obtained for numerical solution of these classes of problems. Numerical results are given to illustrate the efficiency of the proposed method. Results of numerical experiments verify the theoretical behavior of the orders of convergence.     

Positive definiteness in the numerical solution of Riccati differential equations

Summary. In this work we address the issue of integrating symmetric Riccati and Lyapunov matrix differential equations. In many cases -- typical in applications -- the solutions are positive definite matrices. Our goal is to study when and how this property is maintained for a numerically computed solution. There are two classes of solution methods: direct and indirect algorithms. The first class consists of the schemes resulting from direct discretization of the equations. The second class consists of algorithms which recover the solution by exploiting some special formulae that these solutions are known to satisfy. We show first that using a direct algorithm -- a one-step scheme or a strictly stable multistep scheme (explicit or implicit) -- limits the order of the numerical method to one if we want to guarantee that the computed solution stays positive definite. Then we show two ways to obtain positive definite higher order approximations by using indirect algorithms. The first is to apply a symplectic integrator to an associated Hamiltonian system. The other uses stepwise linearization. Received April 21, 1993     

Doubly quasi-consistent parallel explicit peer methods with built-in global error estimation

Recently, Kulikov presented the idea of double quasi-consistency, which facilitates global error estimation and control, considerably. More precisely, a local error control implemented in such methods plays a part of global error control at the same time. However, Kulikov studied only Nordsieck formulas and proved that there exists no doubly quasi-consistent scheme among those methods.Here, we prove that the class of doubly quasi-consistent formulas is not empty and present the first example of such sort. This scheme belongs to the family of superconvergent explicit two-step peer methods constructed by Weiner, Schmitt, Podhaisky and Jebens. We present a sample of s-stage doubly quasi-consistent parallel explicit peer methods of order s−1 when s=3. The notion of embedded formulas is utilized to evaluate efficiently the local error of the constructed doubly quasi-consistent peer method and, hence, its global error at the same time. Numerical examples of this paper confirm clearly that the usual local error control implemented in doubly quasi-consistent numerical integration techniques is capable of producing numerical solutions for user-supplied accuracy conditions in automatic mode.     

Explicit expressions for three-dimensional boundary integrals in linear elasticity

On employing isoparametric, piecewise linear shape functions over a flat triangle, exact formulae are derived for all surface potentials involved in the numerical treatment of three-dimensional singular and hyper-singular boundary integral equations in linear elasticity. These formulae are valid for an arbitrary source point in space and are represented as analytical expressions along the edges of the integration triangle. They can be employed to solve integral equations defined on triangulated surfaces via a collocation method or may be utilized as analytical expressions for the inner integrals in a Galerkin technique. A numerical example involving a unit triangle and a source point located at various distances above it, as well as sample problems solved by a collocation boundary element method for the Lamé equation are included to validate the proposed formulae.     

A fourth-order Runge–Kutta method based on BDF-type Chebyshev approximations

In this paper we consider a new fourth-order method of BDF-type for solving stiff initial-value problems, based on the interval approximation of the true solution by truncated Chebyshev series. It is shown that the method may be formulated in an equivalent way as a Runge–Kutta method having stage order four. The method thus obtained have good properties relatives to stability including an unbounded stability domain and large α$\alpha$-value concerning A(α)$A(\alpha )$-stability. A strategy for changing the step size, based on a pair of methods in a similar way to the embedding pair in the Runge–Kutta schemes, is presented. The numerical examples reveals that this method is very promising when it is used for solving stiff initial-value problems.     

On improving the absolute stability of local extrapolation

Summary A widely used technique for improving the accuracy of solutions of initial value problems in ordinary differential equations is local extrapolation. It is well known, however, that when using methods appropriate for solving stiff systems of ODES, the stability of the method can be seriously degraded if local extrapolation is employed. This is due to the fact that performing local extrapolation on a low order method is equivalent to using a higher order formula and this high order formula may not be suitable for solving stiff systems. In the present paper a general approach is proposed whereby the correction term added on in the process of local extrapolation is in a sense a rational, rather than a polynomial, function. This approach allows high order formulae with bounded growth functions to be developed. As an example we derive anA-stable rational correction algorithm based on the trapezoidal rule. This new algorithm is found to be efficient when low accuracy is requested (say a relative accuracy of about 1%) and its performance is compared with that of the more familiar Richardson extrapolation method on a large set of stiff test problems.