A Zero-Dissipative Phase-Fitted Fourth Order Diagonally Implicit Runge-Kutta-Nyström Method for Solving Oscillatory Problems

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.


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 phaselag 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].

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  11 =  22 = ⋅ ⋅ ⋅ =   = .
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  2 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  1 ,  2 and they will be called as the amplification factors of the RKN method.The consequent eigenvectors are [1,  1 ]  , [1,  2 ]  , where   =   /(  −   ),  = 1, 2. The numerical solution of ( 5) is where If  1 and  2 are complex conjugate, then  1,2 = || exp(±) and  1,2 = || exp(±).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.
Then, we denote that For DIRKN method, let us denote  and  in the following form: From Definition 1 it follows that ) ,          = √  ( 2 ).( 16) 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 () = 0 at a point , the method is said to have zero amplification error (zero-dissipative).Hence, we have Let us denote  1,2 as the roots of ( 7); then we have the following definitions.

Construction of the New Method
In this section, we will present the construction of a zerodissipative 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  31 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  2 and  31 .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  31 which yields where For small value of , we use the Taylor series expansion of  31 as follows: where Hence, a new method is derived and we denote it by PFDIRKN4.This method has one variable coefficient  31 that depends on the product of the step-length ℎ and the frequency V.For each specific product of Vℎ, 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  31 value for different value of  is illustrated in Figure 1.For small value of , () tends to zero; therefore the new method is zerodissipative too.The behavior of ‖ (+1) ‖ 2 and ‖ (+1) ‖ 2 values of PFDIRKN4 for different value of  is illustrated in Figure 2. Table 1 shows a comparison of the properties of the method derived.

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.
(iii) DIRKN4Senu2: the three-stage fourth order DIRKN method with phase-lag order six derived by Senu et al. [15].
(vi) DIRKN4Sharp: the three-stage fourth order with phase-lag of order six DIRKN method by Sharp et al. [12].
The numerical results are plotted in Figures 3, 4, 5, 6, 7, 8, and 9. Figures 3-9 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 3-9 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 5-11, 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.

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.

Table 1 :
Summary of the properties of the methods.