Stability of Runge–Kutta–Nyström methods

AbstractIn this paper, a general and detailed study of linear stability of Runge–Kutta–Nyström (RKN) methods is given. In the case that arbitrarily stiff problems are integrated, we establish a condition that RKN methods must satisfy so that a uniform bound for stability can be achieved. This condition is not satisfied by any method in the literature. Therefore, a stable method is constructed and some numerical comparisons are made

Exploring efficient: numerical methods for differential equations

Numerical analysis is a way to do higher mathematical problems on a computer, a technique widely used by scientists and engineers to solve their problems. A major advantage of numerical analysis is that a numerical answer can be obtained even when a problem has no “analytical” solution. Results from numerical analysis are an approximation, which can be made as accurate as desired. The analysis of errors in numerical methods is a critically important part of the study of numerical analysis. Hence, we will see in this research that computation of the error is a must as it is a way to measure the efficiency of the numerical methods developed. Numerical methods require highly tedious and repetitive computations that can only be done using the computer. Hence in this research, it is shown that computer programs must be written for the implementation of numerical methods. In the early part of related research the computer language used was Fortran. Subsequently more and more computer programs used the C programming language. Additionally, now computations can also be carried out using softwares like MATLAB, MATHEMATICA and MAPLE. Many physical problems that arise from ordinary differential equations (ODEs) have magnitudes of eigenvalues which vary greatly, and such systems are commonly known as stiff systems. Stiff systems usually consist of a transient solution, that is, a solution which varies rapidly at the beginning of the integration. This phase is referred to as the transient phase and during this phase, accuracy rather than stability restricts the stepsize of the numerical methods used. Thus the generally the structure of the solutions suggests application of specific methods for non-stiff equations in the transient phase and specific methods for stiff equations during the steady-state phase in a manner whereby computational costs can be reduced. Consequently, in this research we developed embedded Runge-Kutta methods for solving stiff differential equations so that variable stepsize codes can be used in its implementation. We have also included intervalwise partitioning, whereby the system is considered as non-stiff first, and solved using the method with simple iterations, and once stiffness is detected, the system is solved using the same method, but with Newton iterations. By using variable stepsize code and intervalwise partitioning, we have been able to reduce the computational costs. With the aim of increasing the computational efficiency of the Runge-Kutta methods, we have also developed methods of higher order with less number of stages or function evaluations. The method used is an extension of the classical Runge-Kutta method and the approximation at the current point is based on the information at the current internal stage as well as the previous internal stage. This is the idea underlying the construction of Improved Runge-Kutta methods, so that the resulting method will give better accuracy. Usually higher order ordinary differential equations are solved by converting them into a system of first order ODEs and using numerical methods suitable for first order ODEs. However it is more efficient, in terms of accuracy, number of function evaluations as well as computational time, if the higher order ODEs can be solved directly (without being converted to a system of first order ODEs), using numerical methods. In this research we developed numerical methods, particularly Runge-Kutta type methods, which can directly solve special third order and fourth order ODEs. Special second order ODE is an ODE which does not depend on the first derivative. The solution from this type of ODE often exhibits a pronounced oscillatory character. It is well known that it is difficult to obtain accurate numerical results if the ODEs are oscillatory in nature. In order to address this problem a lot of research has been focused on developing methods which have high algebraic order, reduced phase-lag or dispersion and reduced dissipation. Phaselag is the angle between the true and approximate solution, while dissipation is the difference between the approximate solution and the standard cyclic solution. If a method has high algebraic order, high order of dispersion and dissipation, then the numerical solutions obtained will be very accurate. Hence in this research we have developed numerical methods, specifically hybrid methods which have all the above mentioned properties. If the solutions are oscillatory in nature, it means that the solutions will have components which are trigonometric functions, that is, sine and cosine functions. In order to get accurate numerical solutions we thus phase-fitted the methods using trigonometric functions. In this research, it is proven that trigonometrically-fitting the hybrid methods and applying them to solve oscillatory delay differential equations result in better numerical results. These are the highlights of my research journey, though a lot of work has also been done in developing numerical methods which are multistep in nature, for solving higher order ODEs, as well as implementation of methods developed for solving fuzzy differential equations and partial differential equations, which are not covered here

Two point block multistep methods with trigonometric−fitting for solving oscillatory problems

In this paper, we present the absolute stability of the existing 2-point implicit block multistep step methods of step number k = 3 and k = 5 and solving special second order ordinary differential equations (ODEs). The methods are then trigonometrically fitted so that they are suitable for solving highly oscillatory problems arising from the special second order ODEs. Their explicit counterparts are also trigonometrically fitted so that in the implementation the methods can act as a predictor-corrector pairs. The numerical results based on the integration over a large interval are given to show the performance of the proposed methods. From the numerical results we can conclude that the new trigonometrically-fitted methods are superior in terms of accuracy and execution time, compared to the existing methods in the scientific literature when used for solving problems which are oscillatory in nature

Numerical algorithm of block method for general second order ODEs using variable step size

This paper outlines an alternative algorithm for solving general second order ordinary differential equations (ODEs). Normally, the numerical method was designed for solving higher order ODEs by converting it into an n-dimensional first order equations with implementation of constant step length. Nevertheless, this involved a lot of computational complexity which led to consumption a lot of time. Consequently, a direct block multistep method with utilization of variable step size strategy is proposed. This method was developed for computing the solution at four points simultaneously and the derivation based on numerical integration as well as using interpolation approach. The convergence of the proposed method is justified under suitable conditions of stability and consistency. Five numerical examples are considered and some comparisons are made with the existing methods for demonstrating the validity and reliability of the proposed algorithm

NySALT: Nystr\"{o}m-type inference-based schemes adaptive to large time-stepping

Large time-stepping is important for efficient long-time simulations of deterministic and stochastic Hamiltonian dynamical systems. Conventional structure-preserving integrators, while being successful for generic systems, have limited tolerance to time step size due to stability and accuracy constraints. We propose to use data to innovate classical integrators so that they can be adaptive to large time-stepping and are tailored to each specific system. In particular, we introduce NySALT, Nystr\"{o}m-type inference-based schemes adaptive to large time-stepping. The NySALT has optimal parameters for each time step learnt from data by minimizing the one-step prediction error. Thus, it is tailored for each time step size and the specific system to achieve optimal performance and tolerate large time-stepping in an adaptive fashion. We prove and numerically verify the convergence of the estimators as data size increases. Furthermore, analysis and numerical tests on the deterministic and stochastic Fermi-Pasta-Ulam (FPU) models show that NySALT enlarges the maximal admissible step size of linear stability, and quadruples the time step size of the St\"{o}rmer--Verlet and the BAOAB when maintaining similar levels of accuracy.Comment: 26 pages, 7 figure

エネルギー関数を持つ発展方程式に対する幾何学的数値計算法

学位の種別: 課程博士審査委員会委員 : (主査)東京大学教授 松尾 宇泰, 東京大学教授 中島 研吾, 東京大学准教授 鈴木 秀幸, 東京大学准教授 長尾 大道, 東京大学准教授 齋藤 宣一University of Tokyo(東京大学

Error Analysis of Exponential Integrators for Nonlinear Wave-Type Equations

This thesis is concerned with the time integration of certain classes of nonlinear evolution equations in Hilbert spaces by exponential integrators. We aim to prove error bounds which can be established by including only quantities given by a wellposedness result. In the first part, we consider semilinear wave equations and introduce a class of first- and second-order exponential schemes. A standard error analysis is not possible due to the lack of regularity. We have to employ appropriate filter functions as well as the integration by parts and summation by parts formulas in order to obtain optimal error bounds. In the second part, we propose two exponential integrators of first and second order applied to a class of quasilinear wave-type equations. By a detailed investigation of the differentiability of the right-hand side we derive error bounds in different norms. In the framework we can treat quasilinear Maxwell’s equations in full space and on a smooth domain as well as a class of quasilinear wave equations. In both parts, we include numerical examples to confirm our theoretical findings

Towards New High-Order Operator Splitting Time-Integration Methods

Operator splitting (OS) methods represent a powerful strategy to solve an extensive range of mathematical models in the form of differential equations. They have a long and celebrated history, having been successfully used for well over half a century to provide efficient numerical solutions to challenging problems. In fact, OS methods are often the only viable way to solve many problems in practice. The simplest, and perhaps, most well-known OS methods are Lie--Trotter--Godunov and the Strang--Marchuk methods. They compute a numerical solution that is first-, and second-order accurate in time, respectively. OS methods can be derived by imposing order conditions using the Campbell--Baker--Hausdorff formula. It follows that, by setting the appropriate order conditions, it is possible to derive OS methods of any desired order. An important observation regarding OS methods with order higher than two is that, according to the Sheng--Suzuki theorem, at least one of their defining coefficients must be negative. Therefore, the time integration with OS methods of order higher than two has not been considered suitable to solve deterministic parabolic problems, because the necessary backward time integration would cause instabilities. Throughout this thesis, we focus our attention on high-order (i.e., order higher than two) OS methods. We successfully assess the convergence properties of some higher-order OS methods on diffusion-reaction problems describing cardiac electrophysiology and on an advection-diffusion-reaction problem describing chemical combustion. Furthermore, we compare the efficiency performance of higher-order methods to second-order methods. For all the cases considered, we confirm an improved efficiency performance compared to methods of lower order. Next, we observe how, when using OS and Runge--Kutta type methods to advance the time integration, we can construct a unique extended Butcher tableau with a similar structure to the ones describing Generalized Additive Runge--Kutta (GARK) methods. We define a combination of methods to be OS-GARK methods. We apply linear stability analysis to OS-GARK methods; this allows us to conveniently analyze the stability properties of any combination of OS and Runge--Kutta methods. Doing so, we are able to perform an eigenvalue analysis to understand and improve numerically unstable solutions