Production of the Reduction Formula of Seventh Order Runge-Kutta Method with Step Size Control of an Ordinary Differential Equation

Georgios D. Trikkaliotis, Maria Ch. Gousidou-Koutita

Department of Mathematics, Aristotle University of Thessaloniki, Thessaloniki, Greece.

DOI: 10.4236/am.2022.134023                

Abstract

The purpose of the present work is to construct a nonlinear equation system (85 × 53) using Butcher’s Table and then by solving this system to find the values of all set parameters and finally the reduction formula of the Runge-Kutta (7,9) method (7^th order and 9 stages) for the solution of an Ordinary Differential Equation (ODE). Since the system of high order conditions required to be solved is too complicated, we introduce a subsystem from the original system where all coefficients are found with respect to 9 free parameters. These free parameters, as well as some others in addition, are adjusted in such a way to furnish more efficient R-K methods. We use the MATLAB software to solve several of the created subsystems for the comparison of our results which have been solved analytically.

Keywords

Initial Value Problem, Runge-Kutta Methods, Ordinary Differential Equations

1. Fundamental Principles

From the research carried out in publications of works related to the creation of R-K methods for solving Ordinary Differential Equations it was observed that the description of the creation of these methods was done in a general way and for all classes of normal differential equations. The approach to create this method was not simple but was partly complex. So it was decided to create and propose a process for creating a R-K method that will be simple, understandable and applicable.

A system of ordinary differential equations of the form

$y^{\prime }=f\left(x,y\right),y\left(x_0\right)=y_0$ (1)

with $x_0\in ℝ$, $y,y^{\prime }\in {ℝ}^m$ and $f:ℝ\times {ℝ}^m\to {ℝ}^m$, is called Initial Value Problem (IVP).

Runge-Kutta methods are commonly used numerical methods for addressing (1). They usually presented in a so-called Butcher table (Table 1) [2] [3]:

The table contains on the 1^st column the coefficients c_i, the matrix A with the coefficients of a_ij, which appear in the Formulae of K_i, and w_i the coefficients in Formula of y_i+1.

In this type of table, we have $w^{\text{T}},c\in {ℝ}^m$ while $A\in {ℝ}^{m\times m}$. Then, the method shares m stages and in case that c₁ = 0 and A is strictly lower triangular, it is evaluated explicitly.

The solution of a differential equation is a continuous curve y(x) that passes through the point (x₀,y₀) and satisfies $y^{\prime }\left(x\right)=f\left(x,y\right)$. Numerical solution of a differential equation is a distinct set of values of y(x) which is an approach to the continuous solution of the y(x) curve.

Carl David Tolmé Runge [4] and Martin Wilhelm Kutta [5] introduced the methods bearing their names almost in the turning of the 19th century. Runge and Kutta observed that the derivation of high-order derivatives that appear in the Taylor method can be avoided. In this method we place the problem with indeterminate parameters and make the result at the highest order using calculations of f(x,y) inside (x_n,y_n) and (x_n+1,y_n+1)intervals. The derivatives in the Taylor form are replaced by calculating f(x,y)at a number of points inside (x_n,y_n) and (x_n+1,y_n+1) intervals.

Runge was the first to present a 2nd order R-K method by combining a sequence of Euler formulas [4]. Some years later, Kutta managed to construct a 4 stages 4th order method [5]. Nyström showed a method (5,6) of 5^th order and 6 stages [6]. Fehlberg [7], Shanks [8] and Lawson [9] showed 5^th order methods of 6 stages too. 6^th order methods have been presented by Butcher [2], Fehlberg [7], Shanks [8] and Lawson [9]. Huta’s 6^th order method of 8 stages is the most popular [10]. Higher order R-K methods have been presented by Shanks [8], Felhberg [7], Feagin [11] [12] Hairer [13], Butcher [14] [15], Curtis [16], Famelis [17], Papakostas [17] [18], Tsitouras [17] [18] [19] and others.

Some problems that could be solved in this paper:

· We want with analytic way to derive the RK(7,9) method and we introduce our method for that.

· We give first arbitrary variables with values of the existing table of RK(7,9) method in order to compare our method of solving the non-linear system.

Table 1. The so-called butcher table.

· We suggest and some others arbitrary variables which lead to desired Tables of the RK(7,9) method, because not all the arbitrary values lead to desired results of the method.

Firstly, we present the Introduction of the paper with historical references, and then in section 2 we give an analytic approach of the Runge-Kutta method (7,9) 7^th order with 9 stages method. Finally, we give the conclusions of our work of a certain set of values of the parameters of the method.

2. Presentation of the Runge-Kutta 7^th Order 9 Stages Method

The reduction formula of R-K methods for an ordinary differential equation is given by the relation $y_{n+1}=y_n+{\displaystyle {\sum }_{i=1}^{\nu }{w}_i{K}_i}$ (2), with w_i acting as coefficients of weight, ν the number of steps and $K_i=hf\left(x_n+c_{i}h,y_n+{\displaystyle {\sum }_{j=1}^{i-1}{\alpha }_{ij}{K}_j}\right),c_1=0,i=\text{2},\text{3},\cdots ,\nu$ (3) withh the step of the method. The parameters w_i, c_i and α_ij must be specified. In every R-K method the relations ${\displaystyle {\sum }_{i=1}^{\nu }{w}_i}=1$ (4) και $c_i={\displaystyle {\sum }_{j=1}^{i-1}{\alpha }_{ij}},i=\text{2},\text{3},\cdots ,\nu$ (5) must be valid.

Runge-Kutta (7,9) method is a method of 7^th order and 9 stages and we use the coefficients obtained for $r=1,2,\cdots ,7$ (6), where r is the order of ODE, from Butcher’s Table[1] from whom the equations of the nonlinear 85 × 53 system result.

The values of w_i,c_i and α_ij will be found by the solving this system as well as the K_i and the reduction formula for the solution of the differential equation.

The equations of the system are numbered from (8), (9), ···, (92) and introducing the abbreviation:

$P_{\kappa \lambda }={\alpha }_{\kappa }{}_2{c}_2^{\lambda }+{\alpha }_{\kappa }{}_3{c}_3^{\lambda }+\cdots +{\alpha }_{\kappa }{}_{\kappa -1}{c}_{\kappa -1}^{\lambda }\text{with}\kappa =3,4,\cdots ,9\text{and}\lambda =1,2,3,4,5$ (7)

[7] the following system is obtained:

${\displaystyle {\sum }_{\kappa =1}^9{w}_{\kappa }=1}$ (8)

${\displaystyle {\sum }_{\kappa =2}^9{w}_{\kappa }{c}_{\kappa }=\frac{1}{2}}$ (9)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 1}=\frac{1}{6}}$ (10)

${\displaystyle {\sum }_{\kappa =2}^9{w}_{\kappa }{c}_{\kappa }^2=\frac{1}{3}}$ (11)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}\right)=\frac{1}{24}}$ (12)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 2}=\frac{1}{12}}$ (13)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }{P}_{\kappa 1}=\frac{1}{8}}$ (14)

${\displaystyle {\sum }_{\kappa =2}^9{w}_{\kappa }{c}_{\kappa }^3=\frac{1}{4}}$ (15)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{120}}$ (16)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 2}\right)=\frac{1}{60}}$ (17)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{c}_{\lambda }{P}_{\lambda 1}\right)=\frac{1}{40}}$ (18)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 3}=\frac{1}{20}}$ (19)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}\right)=\frac{1}{30}}$ (20)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }{P}_{\kappa 2}=\frac{1}{15}}$ (21)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 1}^2=\frac{1}{20}}$ (22)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }^2{P}_{\kappa 1}=\frac{1}{10}}$ (23)

${\displaystyle {\sum }_{\kappa =2}^9{w}_{\kappa }{c}_{\kappa }^4=\frac{1}{5}}$ (24)

${\displaystyle {\sum }_{\kappa =6}^9{w}_{\kappa }\left\{{\displaystyle {\sum }_{\lambda =5}^{\kappa -1}{\alpha }_{\kappa \lambda }\left[{\displaystyle {\sum }_{\mu =4}^{\lambda -1}{\alpha }_{\lambda \mu }\left({\displaystyle {\sum }_{\nu =3}^{\mu -1}{\alpha }_{\mu \nu }{P}_{\nu 1}}\right)}\right]}\right\}=\frac{1}{720}}$ (25)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 2}}\right)}\right]=\frac{1}{360}}$ (26)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{c}_{\mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{240}}$ (27)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}\right)=\frac{1}{120}}$ (28)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }{c}_{\lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{180}}$ (29)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{c}_{\lambda }{P}_{\lambda 2}\right)=\frac{1}{90}}$ (30)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}^2\right)=\frac{1}{120}}$ (31)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{c}_{\lambda }^2{P}_{\lambda 1}\right)=\frac{1}{60}}$ (32)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 4}=\frac{1}{30}}$ (33)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }{c}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{144}}$ (34)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 2}\right)=\frac{1}{72}}$ (35)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }{c}_{\lambda }}{P}_{\lambda 1}\right)=\frac{1}{48}}$ (36)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }{P}_{\kappa 3}=\frac{1}{24}}$ (37)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{P}_{\kappa 1}\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}\right)=\frac{1}{72}}$ (38)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 1}{P}_{\kappa 2}=\frac{1}{36}}$ (39)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }^2\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}\right)=\frac{1}{36}}$ (40)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }^2{P}_{\kappa 2}=\frac{1}{18}}$ (41)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }{P}_{\kappa 1}^2=\frac{1}{24}}$ (42)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }^3{P}_{\kappa 1}=\frac{1}{12}}$ (43)

${\displaystyle {\sum }_{\kappa =2}^9{w}_{\kappa }{c}_{\kappa }^5=\frac{1}{6}}$ (44)

${\displaystyle {\sum }_{\kappa =7}^9{w}_{\kappa }\langle {\displaystyle {\sum }_{\lambda =6}^{\kappa -1}{\alpha }_{\kappa \lambda }}\left\{{\displaystyle {\sum }_{\mu =5}^{\lambda -1}{\alpha }_{\lambda \mu }\left[{\displaystyle {\sum }_{\nu =4}^{\mu -1}{\alpha }_{\mu \nu }\left({\displaystyle {\sum }_{\rho =3}^{\nu -1}{\alpha }_{\nu \rho }{P}_{\rho 1}}\right)}\right]}\right\}\rangle =\frac{1}{5040}}$ (45)

${\displaystyle {\sum }_{\kappa =6}^9{w}_{\kappa }\left\{{\displaystyle {\sum }_{\lambda =5}^{\kappa -1}{\alpha }_{\kappa \lambda }\left[{\displaystyle {\sum }_{\mu =4}^{\lambda -1}{\alpha }_{\lambda \mu }\left({\displaystyle {\sum }_{\nu =3}^{\mu -1}{\alpha }_{\mu \nu }{P}_{\nu 2}}\right)}\right]}\right\}=\frac{1}{2520}}$ (46)

${\displaystyle {\sum }_{\kappa =6}^9{w}_{\kappa }\left\{{\displaystyle {\sum }_{\lambda =5}^{\kappa -1}{\alpha }_{\kappa \lambda }\left[{\displaystyle {\sum }_{\mu =4}^{\lambda -1}{\alpha }_{\lambda \mu }\left({\displaystyle {\sum }_{\nu =3}^{\mu -1}{\alpha }_{\mu \nu }{c}_{\nu }{P}_{\nu 1}}\right)}\right]}\right\}=\frac{1}{1680}}$ (47)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 3}}\right)}\right]=\frac{1}{840}}$ (48)

${\displaystyle {\sum }_{\kappa =6}^9{w}_{\kappa }\left\{{\displaystyle {\sum }_{\lambda =5}^{\kappa -1}{\alpha }_{\kappa \lambda }\left[{\displaystyle {\sum }_{\mu =4}^{\lambda -1}{\alpha }_{\lambda \mu }{c}_{\mu }\left({\displaystyle {\sum }_{\nu =3}^{\mu -1}{\alpha }_{\mu \nu }{P}_{\nu 1}}\right)}\right]}\right\}=\frac{1}{1260}}$ (49)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{c}_{\mu }{P}_{\mu 2}}\right)}\right]=\frac{1}{630}}$ (50)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 1}^2}\right)}\right]=\frac{1}{840}}$ (51)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{c}_{\mu }^2{P}_{\mu 1}}\right)}\right]=\frac{1}{420}}$ (52)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 4}\right)=\frac{1}{210}}$ (53)

${\displaystyle {\sum }_{\kappa =6}^9{w}_{\kappa }\left\{{\displaystyle {\sum }_{\lambda =5}^{\kappa -1}{\alpha }_{\kappa \lambda }{c}_{\lambda }\left[{\displaystyle {\sum }_{\mu =4}^{\lambda -1}{\alpha }_{\lambda \mu }\left({\displaystyle {\sum }_{\nu =3}^{\mu -1}{\alpha }_{\mu \nu }{P}_{\nu 1}}\right)}\right]}\right\}=\frac{1}{1008}}$ (54)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }{c}_{\lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 2}}\right)}\right]=\frac{1}{504}}$ (55)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }{c}_{\lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{c}_{\mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{336}}$ (56)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{c}_{\lambda }{P}_{\lambda 3}\right)=\frac{1}{168}}$ (57)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }{P}_{\lambda 1}\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{504}}$ (58)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}{P}_{\lambda 3}\right)=\frac{1}{252}}$ (59)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }{c}_{\lambda }^2\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{252}}$ (60)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{c}_{\lambda }^2{P}_{\lambda 2}\right)=\frac{1}{126}}$ (61)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{c}_{\lambda }{P}_{\lambda 1}^2\right)=\frac{1}{168}}$ (62)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{c}_{\lambda }^3{P}_{\lambda 1}\right)=\frac{1}{84}}$ (63)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 5}=\frac{1}{42}}$ (64)

${\displaystyle {\sum }_{\kappa =6}^9{w}_{\kappa }{c}_{\kappa }\left\{{\displaystyle {\sum }_{\lambda =5}^{\kappa -1}{\alpha }_{\kappa \lambda }\left[{\displaystyle {\sum }_{\mu =4}^{\lambda -1}{\alpha }_{\lambda \mu }\left({\displaystyle {\sum }_{\nu =3}^{\mu -1}{\alpha }_{\mu \nu }{P}_{\nu 1}}\right)}\right]}\right\}=\frac{1}{840}}$ (65)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }{c}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 2}}\right)}\right]=\frac{1}{420}}$ (66)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }{c}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{c}_{\mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{280}}$ (67)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 3}\right)=\frac{1}{140}}$ (68)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }{c}_{\kappa }\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }{c}_{\lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{210}}$ (69)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{c}_{\lambda }{P}_{\lambda 2}\right)=\frac{1}{105}}$ (70)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}^2\right)=\frac{1}{140}}$ (71)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{c}_{\lambda }^2{P}_{\lambda 1}\right)=\frac{1}{70}}$ (72)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }{P}_{\kappa 4}=\frac{1}{35}}$ (73)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }{P}_{\kappa 1}\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{336}}$ (74)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{P}_{\kappa 1}\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 2}\right)=\frac{1}{168}}$ (75)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{P}_{\kappa 1}\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }{c}_{\lambda }}{P}_{\lambda 1}\right)=\frac{1}{112}}$ (76)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 1}{P}_{\kappa 3}=\frac{1}{56}}$ (77)

${\displaystyle {\sum }_{\kappa =5}^9{w}_{\kappa }{c}_{\kappa }^2\left[{\displaystyle {\sum }_{\lambda =4}^{\kappa -1}{\alpha }_{\kappa \lambda }\left({\displaystyle {\sum }_{\mu =3}^{\lambda -1}{\alpha }_{\lambda \mu }{P}_{\mu 1}}\right)}\right]=\frac{1}{168}}$ (78)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_k^2\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 2}\right)=\frac{1}{84}}$ (79)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_k^2\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }{c}_{\lambda }}{P}_{\lambda 1}\right)=\frac{1}{56}}$ (80)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }^2{P}_{\kappa 3}=\frac{1}{28}}$ (81)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}\right)}^2=\frac{1}{252}}$ (82)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{P}_{\kappa 2}\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}\right)=\frac{1}{126}}$ (83)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 2}^2=\frac{1}{63}}$ (84)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }{P}_{\kappa 1}\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}\right)=\frac{1}{84}}$ (85)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }{P}_{\kappa 1}{P}_{\kappa 2}=\frac{1}{42}}$ (86)

${\displaystyle {\sum }_{\kappa =4}^9{w}_{\kappa }{c}_{\kappa }^3\left({\displaystyle {\sum }_{\lambda =3}^{\kappa -1}{\alpha }_{\kappa \lambda }}{P}_{\lambda 1}\right)=\frac{1}{42}}$ (87)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }^3{P}_{\kappa 2}=\frac{1}{21}}$ (88)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{P}_{\kappa 1}^3=\frac{1}{56}}$ (89)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }^2{P}_{\kappa 1}^2=\frac{1}{28}}$ (90)

${\displaystyle {\sum }_{\kappa =3}^9{w}_{\kappa }{c}_{\kappa }^4{P}_{\kappa 1}=\frac{1}{14}}$ (91)

${\displaystyle {\sum }_{\kappa =2}^9{w}_{\kappa }{c}_{\kappa }^6=\frac{1}{7}}$ (92)

In the system of (8), (9), (11), (15), (24), (44) and (92) equations we set as c₂ = c₃ = 1/12, c₄ = 1/6, c₅ = 2/6, c₆ = 3/6,c₇ = 4/6, c₈ = 5/6, c₉ = 6/6. Τhe values of c₂, c₃, ··· , c₉ are chosen to be in ascending order and as small and different from each other as possible. We set in addition $w_2=w_3=0$ (93) and the resulting solution is:

$w_1=w_9=\frac{41}{840}$ (94)

$w_4=w_8=\frac{216}{840}$ (95)

$w_5=w_7=\frac{27}{840}$ (96)

and

$w_6=\frac{272}{840}$ (97)

Since the above equations become somewhat lengthy, we introduce the following abbreviations: [7]

$P_{41}={\alpha }_{42}{c}_2+{\alpha }_{43}{c}_3=\frac{{\alpha }_{42}+{\alpha }_{43}}{12}=\frac{S_4}{12}$ (98)

$P_{51}={\alpha }_{52}{c}_2+{\alpha }_{53}{c}_3+{\alpha }_{54}{c}_4=\frac{{\alpha }_{52}+{\alpha }_{53}+2{\alpha }_{54}}{12}=\frac{S_5}{12}$ (99)

$P_{61}=\frac{{\alpha }_{62}+{\alpha }_{63}+2{\alpha }_{44}+4{\alpha }_{65}}{12}=\frac{S_6}{12}$ (100)

$P_{71}=\frac{{\alpha }_{72}+{\alpha }_{73}+2{\alpha }_{74}+4{\alpha }_{75}+6{\alpha }_{76}}{12}=\frac{S_7}{12}$ (101)

$P_{81}=\frac{{\alpha }_{82}+{\alpha }_{83}+2{\alpha }_{84}+4{\alpha }_{85}+6{\alpha }_{86}+8{\alpha }_{87}}{12}=\frac{S_8}{12}$ (102)

$P_{91}=\frac{{\alpha }_{92}+{\alpha }_{93}+2{\alpha }_{94}+4{\alpha }_{95}+6{\alpha }_{96}+8{\alpha }_{97}+10{\alpha }_{98}}{12}=\frac{S_9}{12}$ (103)

$P_{42}={\alpha }_{42}{c}_2^2+{\alpha }_{43}{c}_3^2=\frac{{\alpha }_{42}+{\alpha }_{43}}{144}=\frac{S_4}{144}$ (104)

$P_{52}={\alpha }_{52}{c}_2^2+{\alpha }_{53}{c}_3^2+{\alpha }_{54}{c}_4^2=\frac{{\alpha }_{52}+{\alpha }_{53}+4{\alpha }_{54}}{144}=\frac{S_5+2{\alpha }_{54}}{144}$ (105)

$P_{62}={\alpha }_{62}{c}_2^2+{\alpha }_{63}{c}_3^2+{\alpha }_{64}{c}_4^2+{\alpha }_{65}{c}_5^2=\frac{S_6+2{\alpha }_{64}+12{\alpha }_{65}}{144}$ (106)

$P_{72}={\alpha }_{72}{c}_2^2+{\alpha }_{73}{c}_3^2+{\alpha }_{74}{c}_4^2+{\alpha }_{75}{c}_5^2+{\alpha }_{76}{c}_6^2=\frac{S_7+2{\alpha }_{74}+12{\alpha }_{75}+30a_{76}}{144}$ (107)

$P_{82}=\frac{S_8+2{\alpha }_{84}+12{\alpha }_{85}+30a_{86}+56a_{87}}{144}$ (108)

$P_{92}=\frac{S_9+2{\alpha }_{94}+12{\alpha }_{95}+30a_{96}+56a_{97}+90a_{98}}{144}$ (109)

$P_{43}=\frac{a_{42}+a_{43}}{1728}=\frac{S_4}{1728}$ (110)

$P_{53}=\frac{a_{52}+a_{53}+8a_{54}}{1728}=\frac{S_5+6{\alpha }_{54}}{1728}$ (111)

$P_{63}=\frac{S_6+6{\alpha }_{64}+60{\alpha }_{65}}{1728}$ (112)

$P_{73}=\frac{S_7+6{\alpha }_{74}+60{\alpha }_{75}+210a_{76}}{1728}$ (113)

$P_{83}=\frac{S_8+6{\alpha }_{84}+60{\alpha }_{85}+210a_{86}+504a_{87}}{1728}$ (114)

$P_{93}=\frac{S_9+6{\alpha }_{94}+60{\alpha }_{95}+210a_{96}+504a_{97}+990a_{98}}{1728}$ (115)

$P_{44}=\frac{S_4}{20736}$ (116)

$P_{54}=\frac{S_5+14{\alpha }_{54}}{20736}$ (117)

$P_{64}=\frac{S_6+14a_{64}+252{\alpha }_{65}}{20736}$ (118)

$P_{74}=\frac{S_7+14a_{74}+252{\alpha }_{75}+1290a_{76}}{20736}$ (119)

$P_{84}=\frac{S_8+14a_{84}+252{\alpha }_{85}+1290a_{86}+4088{\alpha }_{87}}{20736}$ (120)

$P_{94}=\frac{S_9+14a_{94}+252{\alpha }_{95}+1290a_{96}+4088{\alpha }_{97}+9990{\alpha }_{98}}{20736}$ (121)

$P_{45}=\frac{S_4}{248832}$ (122)

$P_{55}=\frac{S_5+30{\alpha }_{54}}{248832}$ (123)

$P_{65}=\frac{S_6+30{\alpha }_{64}+1020{\alpha }_{65}}{248832}$ (124)

$P_{75}=\frac{S_7+30{\alpha }_{74}+1020{\alpha }_{75}+7770{\alpha }_{76}}{248832}$ (125)

$P_{85}=\frac{S_8+30{\alpha }_{84}+1020{\alpha }_{85}+7770{\alpha }_{86}+32760{\alpha }_{87}}{248832}$ (126)

$P_{95}=\frac{S_9+30{\alpha }_{94}+1020{\alpha }_{95}+7770{\alpha }_{96}+32760{\alpha }_{97}+99990{\alpha }_{98}}{248832}$ (127)

Then we substitute the defined abbreviations in the original system, as well as the found values of c₂, c₃, c₄, c₅, c₆, c₇, c₈, c₉, w₂, w₃, w₄, w₅, w₆, w₇, w₈, w₉, and as a result the system is simplified.

In the system of (10), (14), (23), (43) and (91) we express S₅, S₆, S₇, S₈ and S₉ as a function of S₄ and by substituting them in (22) we find that: $S_4=\frac{1}{6}$ (128) and $S_5=\frac{4}{6}$ (129), $S_6=\frac{9}{6}$ (130), $S_7=\frac{16}{6}$ (131), $S_8=\frac{25}{6}$ (132), $S_9=\frac{36}{6}=6$ (133).

To continue we set ${\alpha }_{42}=0$ (134). From the abbreviation $S_4={\alpha }_{42}+{\alpha }_{\text{43}}$ and from the relation ${\alpha }_{41}+{\alpha }_{42}+{\alpha }_{43}=c_{\text{4}}$ we obtain that ${\alpha }_{43}=\frac{2}{12}$ (135) and ${\alpha }_{41}=0$ (136).

In the system of equations (13), (18), (21), (32), (33), (36), (37), (39), (41), (62), (63), (64), (72), (73), (76), (77), (80), (81), (86) and (88) we substitute the values found above, omitting the equations which are a linear combination of equations of the system and also considering α₉₄, α₉₅, α₉₆, α₉₇, α₉₈ as parameters, the 10 × 15 linear system (A) is obtained:

$27{\alpha }_{54}=-822+41\left({\alpha }_{94}+6{\alpha }_{95}+15{\alpha }_{96}+28{\alpha }_{97}+45{\alpha }_{98}\right)$ (137)

$272{\alpha }_{64}=102-41\left(4{\alpha }_{94}+18{\alpha }_{95}+30{\alpha }_{96}+28{\alpha }_{97}\right)$ (138)

$272{\alpha }_{65}=653-41\left({\alpha }_{95}+5{\alpha }_{96}+14{\alpha }_{97}+30{\alpha }_{98}\right)$ (139)

$27{\alpha }_{74}=102+41\left(6{\alpha }_{94}+18{\alpha }_{95}+15{\alpha }_{96}\right)$ (140)

$27{\alpha }_{75}=-208+41\left(3{\alpha }_{95}+10{\alpha }_{96}+14{\alpha }_{97}\right)$ SYSTEM (A) (141)

$270{\alpha }_{76}=-2560+41\left(10{\alpha }_{96}+56{\alpha }_{97}+180{\alpha }_{98}\right)$ (142)

$216{\alpha }_{84}=-822-41\left(4{\alpha }_{94}+6{\alpha }_{95}\right)$ (143)

$216{\alpha }_{85}=653-41\left(3{\alpha }_{95}+5{\alpha }_{96}\right)$ (144)

$216\left(5{\alpha }_{86}\right)=-1280-41\left(10{\alpha }_{96}+28{\alpha }_{97}\right)$ (145)

$216\left(7{\alpha }_{87}\right)=1683-41\left(7{\alpha }_{97}+45{\alpha }_{98}\right)$ (146)

From (20) and (85) equations we obtain:

$\begin{array}{l}27\left({\alpha }_{54}\right)+272\left[4\left({\alpha }_{64}+4{\alpha }_{65}\right)\right]+27\left[10\left({\alpha }_{74}+4{\alpha }_{75}+9{\alpha }_{76}\right)\right]\\ +216\left[20\left({\alpha }_{84}+4{\alpha }_{85}+9{\alpha }_{86}+16{\alpha }_{87}\right)\right]\\ +41\left[35\left({\alpha }_{94}+4{\alpha }_{95}+9{\alpha }_{96}+16{\alpha }_{97}+25{\alpha }_{98}\right)\right]=49824\end{array}$ (147)

which along with the equations of system (A) and after setting: ${\alpha }_{65}={\alpha }_{75}={\alpha }_{84}={\alpha }_{86}=0$ (148) results that:

${\alpha }_{54}=\frac{160}{9}$ (149)

${\alpha }_{64}=-\frac{297}{68}$ (150)

${\alpha }_{65}=0$ (151)

${\alpha }_{74}=\frac{157}{3}$ (152)

${\alpha }_{75}=0$ (153)

${\alpha }_{76}=\frac{157}{45}$ (154)

${\alpha }_{84}=0$ (155)

${\alpha }_{85}=-\frac{65}{72}$ (156)

${\alpha }_{86}=0$ (157)

${\alpha }_{87}=\frac{29}{90}$ (158)

${\alpha }_{94}=-\frac{1211}{326}$ (159)

${\alpha }_{95}=-\frac{205}{237}$ (160)

${\alpha }_{96}=\frac{419}{90}$ (161)

${\alpha }_{96}=-\frac{25}{9}$ (162)

${\alpha }_{98}=\frac{27}{25}$ (163)

We found above α₄₃ = 1/6 and setting ${\alpha }_{32}=\frac{11}{12}$ (164) ( ${\alpha }_{31}=-\frac{10}{12}$ (165)), from the system of (30), (34), (58), (60) equations, implies that: ${\alpha }_{53}=\frac{4}{9}$ (166) ${\alpha }_{63}=\frac{18}{5}$ (167) ${\alpha }_{73}=-\frac{23}{2}$ (168) ${\alpha }_{83}=-\frac{3}{14}$ (169) and from Equation (12) result ${\alpha }_{93}=-\frac{119}{6}$ (170).

From the abbreviations: S₅ = 4/6, S₆ = 9/6, S₇ = 16/6, S₈ = 25/6 and S₉ = 36/6 = 6 results that: ${\alpha }_{52}=-\frac{106}{3}$ (171) ${\alpha }_{62}=\frac{199}{30}$ (172) ${\alpha }_{72}=-\frac{683}{6}$ (173) ${\alpha }_{82}=\frac{379}{70}$ (174) and ${\alpha }_{92}=\frac{283}{14}$ (175).

From relations ${\displaystyle {\sum }_{j=1}^{i-1}{\alpha }_{ij}}=c_i$ , $i=2,3,4,5,6,7,8,9$ we obtain: ${\alpha }_{51}=\frac{157}{9}$ (176) ${\alpha }_{61}=-\frac{161}{30}$ (177) ${\alpha }_{71}=\frac{3158}{45}$ (178) ${\alpha }_{81}=-\frac{53}{14}$ (179) and ${\alpha }_{91}=\frac{56}{25}$ (180).

According to the so-called Butcher’s Table(Table2) the (7,9) R-K method is given as below:

therefore

$K_1=hf\left(x_n,y_n\right)$ (181)

$K_2=hf\left(x_n+\frac{h}{12},y_n+\frac{K_1}{12}\right)$ (182)

$K_3=hf\left(x_n+\frac{h}{12},y_n+\frac{-10K_1+11K_2}{12}\right)$ (183)

$K_4=hf\left(x_n+\frac{2h}{12},y_n+\frac{2K_3}{12}\right)$ (184)

Table 2. For choices values of arbitrary constants. c₂ = c₃ = 1/12, c₄ = 1/6, c₅ = 2/6, c₆ = 3/6, c₇ = 4/6, c₈ = 5/6, c₉ = 6/6.

$K_5=hf\left(x_n+\frac{4h}{12},y_n+\frac{157K_1-318K_2+4K_3+160K_4}{9}\right)$ (185)

$K_6=hf\left(x_n+\frac{6h}{12},y_n+\frac{-322K_1+199K_2+108K_3-131K_5}{30}\right)$ (186)

$K_7=hf\left(x_n+\frac{8h}{12},y_n+\frac{3158K_1}{45}-\frac{638K_2}{6}-\frac{23K_3}{2}+\frac{157K_4}{3}+\frac{157K_6}{45}\right)$ (187)

$K_8=hf\left(x_n+\frac{10h}{12},y_n-\frac{53K_1}{14}+\frac{38K_2}{7}-\frac{3K_3}{14}-\frac{65K_5}{72}+\frac{29K_7}{90}\right)$ (188)

$\begin{array}{l}{K}_9=hf(x_n+h,y_n+\frac{56K_1}{25}+\frac{283k_2}{14}-\frac{119K_3}{6}-\frac{26K_4}{7}-\frac{13K_5}{15}\\ +\frac{149K_6}{32}-\frac{25K_7}{9}+\frac{27K_8}{25})\end{array}$ (189)

and the reduction formula for the solution of the Differential Equation is:

$y_{n+1}=y_n+\frac{41K_1+216K_4+27K_5+272K_6+27K_7+216K_8+41K_9}{840}$ (190)

3. Conclusion

This paper is concerned with training the coefficients of a 7^th order and 9 stages Runge-Kutta method for addressing initial value problems. As the presented method is 9 stages, we use a set of 9 free parameters. After optimizing the free parameters (coefficients), we concluded to a certain set of values of them. This set of values was found to outperform other representatives in a wide range of relevant problems.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Butcher, J.C. (1963) Coefficients for the Study of Runge-Kutta Integration Processes. Journal of the Australian Mathematical Society, 3, 185-201. https://doi.org/10.1017/S1446788700027932
[2]	Butcher, J.C. (1964) On Runge-Kutta Processes of High Order. Journal of the Australian Mathematical Society, 4, 179-194. https://doi.org/10.1017/S1446788700023387
[3]	Butcher, J.C. (2003) Numerical Methods for Ordinary Differential Equations. John Wiley & Sons, Chichester. https://doi.org/10.1002/0470868279
[4]	Runge, C. (1895) Ueber die numerische Auflöung von Differentialgleichungen. Mathematische Annalen, 46, 167-178. https://doi.org/10.1007/BF01446807
[5]	Kutta, W. (1901) Beitrag zur naherungsweisen Integration von Differentialgleichungen. Zeitschrift für angewandte Mathematik und Physik, 46, 435-453.
[6]	Nyström, E.J. (1925) Ueber die numerische Integration von Differentialgleichungen. Acta Societatis Scientiarum Fennicae, 50, 1-54.
[7]	Fehlberg, E. (1968) Classical Fifth-, Sixth-, Seventh-, and Eight-Order Runge-Kutta Formulas with Stepsize Control. NASA Technical Report No. R-287, National Aeronautics and Space Administration, Washington DC, 4-13.
[8]	Shanks, M. (1966) Solutions of Differential Equations by Evaluation of Functions. Mathematics of Computation, 20, 21-38. https://doi.org/10.1090/S0025-5718-1966-0187406-1
[9]	Gousidou-Koutita, M. (2009) Numerical Methods with Applications of Ordinary and Partial Differential Equations. University Lectures 2008/9, Aristotle University of Thessaloniki, Thessaloniki, Greece.
[10]	Hŭta, A. (1956) Une amélioration de la méthode de Runge-Kutta-Nyström pour la résolution numérique des équations différentielles du premièr ordre. Acta Facultatis Berum Naturalium Universitatis Comenianae, 1, 201-224.
[11]	Feagin, T. (2009) High-Order Explicit Runge-Kutta Methods Using m-Symmetry. University of Houston, Houston.
[12]	Feagin, T. (2006) A Tenth-Order Runge-Kutta Method with Error Estimate. University of Houston, Houston.
[13]	Hairer, E. (1978) A Runge-Kutta Method of Order 10. IMA Journal of Applied Mathematics, 21, 47-59. https://doi.org/10.1093/imamat/21.1.47
[14]	Butcher, J.C. (1964) Implicit Runge-Kutta Processes. Mathematics of Computation, 18, 50-64. https://doi.org/10.1090/S0025-5718-1964-0159424-9
[15]	Butcher, J.C. (1965) On the Attainable Order of Runge-Kutta Methods. Mathematics of Computation, 19, 408-417. https://doi.org/10.1090/S0025-5718-1965-0179943-X
[16]	Curtis, A.R. (1990) An Eighth Order Runge-Kutta Process with Eleven Function Evaluations per Step. Numerische Mathematik, 16, 268-277. https://doi.org/10.1007/BF02219778
[17]	Famelis, I.T.H., Papakostas, S.N. and Tsitouras, C. (2004) Symbolic Derivation of Runge-Kutta Order Conditions. Journal of Symbolic Computation, 37, 311-327. https://doi.org/10.1016/j.jsc.2003.07.001
[18]	Papakostas, S.N. and Tsitouras, C. (1999) High Phase-Lag Order Runge-Kutta and Nystrom Pairs. SIAM Journal on Scientific Computing, 21, 747-763. https://doi.org/10.1137/S1064827597315509
[19]	Tsitouras, C. (2001) Optimized Explicit Runge-Kutta Pair of Orders 9(8). Applied Numerical Mathematics, 38, 121-134. https://doi.org/10.1016/S0168-9274(01)00025-3