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.     

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.     

On quasi-consistent integration by Nordsieck methods

In this paper we study advantages of numerical integration by quasi-consistent Nordsieck formulas. All quasi-consistent numerical methods possess at least one important property for practical use, which has not attracted attention yet, i.e. the global error of a quasi-consistent method has the same order as its local error. This means that the usual local error control will produce a numerical solution for the prescribed accuracy requirement if the principal term of the local error dominates strongly over remaining terms. In other words, the global error control can be as cheap as the local error control in the methods under discussion.     

Asymptotically correct error estimation for collocation methods applied to singular boundary value problems

We discuss an a posteriori error estimate for collocation methods applied to boundary value problems in ordinary differential equations with a singularity of the first kind. As an extension of previous results we show the asymptotical correctness of our error estimate for the most general class of singular problems where the coefficient matrix is allowed to have eigenvalues with positive real parts. This requires a new representation of the global error for the numerical solution obtained by piecewise polynomial collocation when applied to our problem class.     

Local error estimation for multistep collocation methods

Multistep collocation methods for initial value problems in ordinary differential equations are known to be a subclass of multistep Runge-Kutta methods and a generalisation of the well-known class of one-step collocation methods as well as of the one-leg methods of Dahlquist. In this paper we derive an error estimation method of embedded type for multistep collocation methods based on perturbed multistep collocation methods. This parallels and generalizes the results for one-step collocation methods by Nørsett and Wanner. Simple numerical experiments show that this error estimator agrees well with a theoretical error estimate which is a generalisation of an error estimate first derived by Dahlquist for one-leg methods.     

Boundary value methods for the solution of differential-algebraic equations

Summary Boundary value techniques for the solution of initial value problems of ODEs, despite their apparent higher cost, present some important advantages over initial value methods. Among them, there is the possibility to have greater accuracy, to control the global error, and to have an efficient parallel implementation.In this paper, the same techniques are applied to the solution of linear initial value problems of DAEs. We have considered three term numerical methods (Midpoint, Simpson, and an Adams type method) in order to obtain a block tridiagonal linear system as a discrete problem.Convergence results are stated in the case of constant coefficients, and numerical examples are given on linear time-varying problems.Work supported by the Ministero della Ricerca Scientifica, 40% project, and by the C.N.R. (contract of research # 92.00535.CT01)     

Family of symplectic implicit Runge-Kutta formulae

A family of formulae for the sympletic IRK method is investigated. Specifically, focus is given to general solutions for formula parameters of IRK under the symplectic and the order conditions. Examples of such formulae are constructed for up to three stages.     

Time-stepping and preserving orthonormality

Certain applications produce initial value ODEs whose solutions, regarded as time-dependent matrices, preserve orthonormality. Such systems arise in the computation of Lyapunov exponents and the construction of smooth singular value decompositions of parametrized matrices. For some special problem classes, there exist time-stepping methods that automatically inherit the orthonormality preservation. However, a more widely applicable approach is to apply a standard integrator and regularly replace the approximate solution by an orthonormal matrix. Typically, the approximate solution is replaced by the factorQ from its QR decomposition (computed, for example, by the modified Gram-Schmidt method). However, the optimal replacement—the one that is closest in the Frobenius norm—is given by the orthonormal polar factor. Quadratically convergent iteration schemes can be used to compute this factor. In particular, there is a matrix multiplication based iteration that is ideally suited to modern computer architectures. Hence, we argue that perturbing towards the orthonormal polar factor is an attractive choice, and we consider performing a fixed number of iterations. Using the optimality property we show that the perturbations improve the departure from orthonormality without significantly degrading the finite-time global error bound for the ODE solution. Our analysis allows for adaptive time-stepping, where a local error control process is driven by a user-supplied tolerance. Finally, using a recent result of Sun, we show how the global error bound carries through to the case where the orthonormal QR factor is used instead of the orthonormal polar factor. This work was supported by Engineering and Physical Sciences Research Council grants GR/H94634 and GR/K80228.     

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.     

Local error control inSDIRK-methods

This paper describes some problems that are encountered in the implementation of a class of Singly Diagonally Implicit Runge-Kutta (SDIRK) methods. The contribution to the local error from the local truncation error and the residual error from the algebraic systems involved are analysed. A section describes a special interpolation formula. This is used as a prediction stage in the iterative solution of the algebraic equations. A strategy for computing a starting stepsize is presented. The techniques are applied to numerical examples.     

On the convergence of Runge–Kutta methods for stiff non linear differential equations

Summary. This paper studies the convergence properties of general Runge–Kutta methods when applied to the numerical solution of a special class of stiff non linear initial value problems. It is proved that under weaker assumptions on the coefficients of a Runge–Kutta method than in the standard theory of B-convergence, it is possible to ensure the convergence of the method for stiff non linear systems belonging to the above mentioned class. Thus, it is shown that some methods which are not algebraically stable, like the Lobatto IIIA or A-stable SIRK methods, are convergent for the class of stiff problems under consideration. Finally, some results on the existence and uniqueness of the Runge–Kutta solution are also presented. Received November 18, 1996 / Revised version received October 6, 1997     

Interior estimates for the wavelet Galerkin method

Summary. In this paper we derive an interior estimate for the Galerkin method with wavelet-type basis. Such an estimate follows from interior Galerkin equations which are common to a class of methods used in the solution of elliptic boundary value problems. We show that the error in an interior domain can be estimated with the best order of accuracy possible, provided the solution is sufficiently regular in a slightly larger domain, and that an estimate of the same order exists for the error in a weaker norm (measuring the effects from outside the domain ). Examples of the application of such an estimate are given for different problems. Received May 17, 1995 / Revised version received April 26, 1996     

Error growth analysis via stability regions for discretizations of initial value problems

This paper deals with numerical methods for the solution of linear initial value problems. Two main theorems are presented on the stability of these methods. Both theorems give conditions guaranteeing a mild error growth, for one-step methods characterized by a rational function ϕ(z). The conditions are related to the stability regionS={z:z∈ℂ with |ϕ(z)|≤1}, and can be viewed as variants to the resolvent condition occurring in the reputed Kreiss matrix theorem. Stability estimates are presented in terms of the number of time stepsn and the dimensions of the space. The first theorem gives a stability estimate which implies that errors in the numerical process cannot grow faster than linearly withs orn. It improves previous results in the literature where various restrictions were imposed onS and ϕ(z), including ϕ′(z)≠0 forz∈σS andS be bounded. The new theorem is not subject to any of these restrictions. The second theorem gives a sharper stability result under additional assumptions regarding the differential equation. This result implies that errors cannot grow faster thann ^β, with fixed β<1. The theory is illustrated in the numerical solution of an initial-boundary value problem for a partial differential equation, where the error growth is measured in the maximum norm.     

Efficient sixth order methods for nonlinear oscillation problems

A class of two-step (hybrid) methods is considered for solving pure oscillation second order initial value problems. The nonlinear system, which results on applying methods of this type to a nonlinear differential system, may be solved using a modified Newton iteration scheme. From this class the author has derived methods which are fourth order accurate,P-stable, require only two (new) function evaluations per iteration and have a true real perfect square iteration matrix. Now, we propose an extension to sixth order,P-stable methods which require only three (new) function evaluations per iteration and for which the iteration matrix is a true realperfect cube. This implies that at most one real matrix must be factorised at each step. These methods have been implemented in a new variable step, local error controlling code.     

Convergence of linear multistep and one-leg methods for stiff nonlinear initial value problems

To prove convergence of numerical methods for stiff initial value problems, stability is needed but also estimates for the local errors which are not affected by stiffness. In this paper global error bounds are derived for one-leg and linear multistep methods applied to classes of arbitrarily stiff, nonlinear initial value problems. It will be shown that under suitable stability assumptions the multistep methods are convergent for stiff problems with the same order of convergence as for nonstiff problems, provided that the stepsize variation is sufficiently regular.     

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.     

Approximating Runge-Kutta matrices by triangular matrices

The implementation of implicit Runge-Kutta methods requires the solution of large systems of non-linear equations. Normally these equations are solved by a modified Newton process, which can be very expensive for problems of high dimension. The recently proposed triangularly implicit iteration methods for ODE-IVP solvers [5] substitute the Runge-Kutta matrixA in the Newton process for a triangular matrixT that approximatesA, hereby making the method suitable for parallel implementation. The matrixT is constructed according to a simple procedure, such that the stiff error components in the numerical solution are strongly damped. In this paper we prove for a large class of Runge-Kutta methods that this procedure can be carried out and that the diagnoal entries ofT are positive. This means that the linear systems that are to be solved have a non-singular matrix. The research reported in this paper was supported by STW (Dutch Foundation for Technical Sciences).     

Efficient fourth orderP-stable formulae

In this paper, a family of fourth orderP-stable methods for solving second order initial value problems is considered. When applied to a nonlinear differential system, all the methods in the family give rise to a nonlinear system which may be solved using a modified Newton method. The classical methods of this type involve at least three (new) function evaluations per iteration (that is, they are 3-stage methods) and most involve using complex arithmetic in factorising their iteration matrix. We derive methods which require only two (new) function evaluations per iteration and for which the iteration matrix is a true real perfect square. This implies that real arithmetic will be used and that at most one real matrix must be factorised at each step. Also we consider various computational aspects such as local error estimation and a strategy for changing the step size.