Abstract

A new diagonally implicit Runge-Kutta-Nyström (DIRKN) method is constructed for solving second order differential equations with oscillatory solutions. The method is originally based on existing DIRKN method derived by Senu et al. which is three-stage and fourth algebraic order. The new derived method has a variable coefficient with phase-lag of order infinity. The numerical experiments are carried out and the results show the efficiency and accuracy of the new method in comparison with the other DIRKN methods in the literature.

1. Introduction

In this paper, we are dealing with the initial value problems (IVPs) related to second order ordinary differential equations (ODEs) in the form: where does not appear explicitly and know that their solutions are periodic. These problems arise in many fields of applied sciences such as astronomy, quantum mechanics, physical chemistry, structural mechanics, and electronics.

For solving this kind of problems, one of the numerical methods is Runge-Kutta-Nyström (RKN) method. RKN method is frequently used due to its computational advantage (see [1]). RKN method consists of two main classes; they are explicit method and implicit method. Generally, it is quite difficult to handle oscillatory problems of the form (1) with classical RKN methods. The term phase-lag was first introduced by Brusa and Nigro [2]. Phase-lag is the angle between the analytical solution and the numerical solution. For solving oscillatory problems, phase-lag is an important property. In [3], van der Houwen and Sommeijer implemented the phase-lag theory to Runge-Kutta (RK) and RKN methods. They presented a few explicit RK and RKN methods with reduced phase errors. Based on the minimal phase-lag theory, Senu et al. constructed a zero-dissipative RKN method in [4]. In [5], Simos derived a Runge-Kutta-Fehlberg method with phase-lag of order infinity. In [6], Papadopoulos et al. extended the idea of phase-lag of order infinity to RKN method and presented a phase-fitted RKN method. Based on idea of phase-lag of order infinity, Papadopoulos and Simos introduced a new methodology to construct optimized RKN methods in [7]. Since then, Kosti et al. constructed two optimized RKN methods with fifth algebraic order in [8, 9]. Notice that all methods mentioned above are explicit methods.

In [10], Sommeijer presented two DIRKN methods with nonempty interval of periodicity. van der Houwen and Sommeijer extended the phase-lag theory to DIRKN methods and presented a few DIRKN methods with high phase-lag order for solving oscillatory problems in [11]. Sharp et al. also presented a few two-stage and three-stage DIRKN methods with high phase-lag order in [12]. Senu et al. extended the DIRKN methods to phase-lag order up to order eight and higher dissipative order in [13, 14]. However, there is no DIRKN method with phase-lag of order infinity being done yet. This motivates us to develop the DIRKN method with phase-lag of order infinity. In this paper, we will construct a three-stage phase-fitted DIRKN method which is based on a three-stage method of algebraic order four derived by Senu et al. [15].

2. Phase Properties of RKN Method

The general form of a -stage implicit RKN method for (1) is given by where

The DIRKN method above can be expressed nicely in a Butcher table as shown as follows:

For diagonally implicit RKN methods, the diagonal elements are equal. We denote the diagonal elements as so that .

The phase-lag error of method (2) is investigated by using the homogeneous test equation in the following:

Applying method (2) to the test equation (5) yields where is the stability matrix of the RKN method and ,   ,   , and  are polynomials in terms of and totally determined by the parameters of method (2). The characteristic equation of can be written as which is the stability polynomial of the RKN method.

It is given that the exact solution of (5) is where Substituting (9) into (8) yields Then, we assume that the eigenvalues of are ,   and they will be called as the amplification factors of the RKN method. The consequent eigenvectors are ,   , where ,   . The numerical solution of (5) is where If and are complex conjugate, then and . By substituting both into (10), we have

Hence, we have the exact solution (10) and the numerical solution (13) of (5) in the similar form. From (10) and (13) we have the following definition.

Definition 1 (phase-lag and amplification error; see [3]). Apply the RKN method (2) to the test equation (5). Then we define the phase-lag as . If , then the RKN method is said to have phase-lag order . In addition, the quantity is called amplification error. If , then the RKN method is said to have dissipation order .

Then, we denote that For DIRKN method, let us denote and in the following form: From Definition 1 it follows that From (12), it leads us to the definition of phase-lag of order infinity.

Definition 2 (phase-lag of order infinity; see [6]). To obtain phase-lag of order infinity the following relation must hold

In addition, when at a point , the method is said to have zero amplification error (zero-dissipative). Hence, we have

Let us denote as the roots of (7); then we have the following definitions.

Definition 3 (absolute stability interval; see [16]). An interval is called the interval of absolute stability of the method (2) if, for all ,   .

Definition 4 (periodicity interval; see [10]). An interval is called the interval of periodicity of method (2) if, for all , are complex conjugate and of modulus one.

3. Construction of the New Method

In this section, we will present the construction of a zero-dissipative DIRKN method with phase-lag of order infinity. The method is originally based on a three-stage DIRKN method with algebraic order four (see Senu et al. [15]) as shown as follows:

Since we want to derive a new method with phase-lag order of order infinity, there must be a variable coefficient. Therefore, we set as free parameter first but let the rest of the coefficients remain the same. Then, we have to compute the stability matrix which is determined by the coefficients in the Butcher table. From matrix , we can compute and in form of (11):

Now, we have and in terms of and . From Definition 2, condition (13) must be satisfied in order to have phase-lag of order infinity. Therefore, we substitute and from (15) into (13) and solve the equation for which yields where

For small value of , we use the Taylor series expansion of as follows: where

Hence, a new method is derived and we denote it by PFDIRKN4. This method has one variable coefficient that depends on the product of the step-length and the frequency . For each specific product of , it helps to nullify the phase-lag error. Hence, the accuracy of the method always depends on the nature of the problem and the properties of the method. For solving oscillatory problems, reducing the phase-lag error of the method is far more important than decreasing its algebraic error. The behavior of value for different value of is illustrated in Figure 1. For small value of , tends to zero; therefore the new method is zero-dissipative too. The behavior of and values of PFDIRKN4 for different value of is illustrated in Figure 2. Table 1 shows a comparison of the properties of the method derived.

4. Problems Tested and Numerical Results

In this section, we will apply the new method to some second order differential equation problems. The following DIRKN methods are selected for the numerical comparisons.(i)PFDIRKN4: the new derived fourth order zero-dissipative phase-fitted DIRKN method.(ii)DIRKN4Senu: the three-stage fourth order DIRKN method derived by Senu et al. [15].(iii)DIRKN4Senu2: the three-stage fourth order DIRKN method with phase-lag order six derived by Senu et al. [15].(iv)DIRKN4Franco: the A-stable fourth order DIRKN method derived by Franco and Gómez [16].(v)DIRKN4Raed: the fourth order DIRKN method derived by Al-Khasawneh et al. [17].(vi)DIRKN4Sharp: the three-stage fourth order with phase-lag of order six DIRKN method by Sharp et al. [12].(vii)DIRKN4Som: the fourth order with phase-lag of order four by Sommeijer [10].

The accuracy criteria taken are calculating the of the maximum global absolute error: where . We test the problems in the following for each method with the step-size ,   .

Problem 1 (homogeneous). Consider Exact solution: .Estimated frequency: .

Problem 2 (inhomogeneous). Consider This case is using .Exact solution: .Estimated frequency: .Source: van der Houwen and Sommeijer [3].

Problem 3 (an almost periodic orbital problem). Consider Exact solution: ,   .Estimated frequency: .Source: Stiefel and Bettis [18].

Problem 4 (two-body problem). Consider Exact solution: ,   .Estimated frequency: .Source: Dormand et al. [19].

Problem 5 (inhomogeneous system). Consider Exact solution: ,   , and parameters and are 20 and 0.1, respectively.Estimated frequency: .Source: Lambert and Watson [20].

Problem 6 (Duffing problem). Consider where and .Exact solution: , where , , , , and .Estimated frequency: .Source: van de Vyver [21].

Problem 7 (linear Strehmel-Weiner problem). Consider Exact solution: Estimated frequency: .Source: Cong [22].

The numerical results are plotted in Figures 3, 4, 5, 6, 7, 8, and 9. Figures 39 display the efficiency curves where the accuracy versus the computational cost measured by the time used by each method in a same computation machine. Figures 39 show the efficiency curves of the methods that consist of PFDIRKN4, DIRKN4Senu, DIRKN4Senu2, DIRKN4Franco, DIRKN4Raed, DIRKN4Sharp, and DIRKN4Som. These efficiency curves display a clear comparison among those methods.

From Figures 3 to 9, we found that the new method PFDIRKN4 is the most accurate and efficient for solving Problems 17, followed by DIRKN4Senu2, DIRKN4Senu, DIRKN4Som, DIRKN4Sharp, DIRKN4Franco, and DIRKN4Raed. From Figures 3 to 9, we can conclude that high phase-lag order of the method is very important for solving oscillating problems. Since both methods PFDIRKN4 and DIRKN4Senu2 are far more accurate than other methods which are having higher phase-lag order. In an addition, we can observe that the new method is having the same computational cost with the corresponding original method but higher accuracy.

5. Conclusion

In this paper, we have derived a new DIRKN method for solving second order IVPs which are oscillatory in nature. The new method is derived based on a fourth algebraic order DIRKN method [15] with zero dissipation. The new method has a variable coefficient that will help to nullify the phase-lag error. Numerical results show that the new method is more accurate and efficient for solving second order differential equations with oscillating solutions in nature.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.