About this Journal Submit a Manuscript Table of Contents
Abstract and Applied Analysis
Volume 2013 (2013), Article ID 910624, 11 pages
http://dx.doi.org/10.1155/2013/910624
Research Article

The Use of Phase Lag and Amplification Error Derivatives for the Construction of a Modified Runge-Kutta-Nyström Method

1Department of Business Administration, Faculty of Management and Economy, Technological Educational Institute of Patras, 26 334 Patras, Greece
2Department of Mathematics, College of Sciences, King Saud University, P.O. Box 2455, Riyadh 11451, Saudi Arabia
3Laboratory of Computational Sciences, Department of Computer Science and Technology, Faculty of Sciences and Technology, University of Peloponnese, 22 100 Tripolis, Greece

Received 14 November 2012; Revised 16 March 2013; Accepted 3 April 2013

Academic Editor: Juan Carlos Cortés López

Copyright © 2013 D. F. Papadopoulos and T. E. Simos. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

A new modified Runge-Kutta-Nyström method of fourth algebraic order is developed. The new modified RKN method is based on the fitting of the coefficients, due to the nullification not only of the phase lag and of the amplification error, but also of their derivatives. Numerical results indicate that the new modified method is much more efficient than other methods derived for solving numerically the Schrödinger equation.

1. Introduction

Among the most commonly used methods in the numerical integration of second order differential equations, in which the first derivative terms are omitted, are Runge-Kutta (-Nyström) methods. These methods have been used widely due to their simplicity and their accuracy and mainly because they are one-step methods and thus they require no additional starting values.

In the last decades many researchers developed optimized Runge-Kutta (-Nyström) methods, based mostly on the exponential fitting and the phase lag properties [114].

In the last years, a new methodology has been developed, which is based on the phase lag derivatives. Researchers that have been using the mentioned methodology achieve higher accuracy at their methods [6, 8, 10, 11, 13, 14].

In the present paper a new modified Runge-Kutta-Nyström method is constructed. The new method contains four additional variable coefficients (in comparison with the classical RKN method), which depend on , where is the dominant frequency of the problem and is the step length of integration. In order to evaluate the above coefficients, the new method combines the nullification of phase lag and amplification factor with the nullification of their derivatives.

The new modified RKN method that has been obtained will be used for the numerical solution of second order differential equations and more specifically for the numerical integration of the radial Schrödinger equation.

2. Modified Runge-Kutta-Nyström Method

The -stage modified RKN method for the equation is of the form where

For the classical method , . In the present paper and based on the requirement of the development of the new method, values , are variable and depend on (which is the product of the frequency and the step size ). In Section 4 we will present a development of a four-stage modified Runge-Kutta-Nyström method of algebraic order four.

3. Phase Lag and Amplification Error Analysis of the Modified Runge-Kutta-Nyström Methods

Consider the problem with exact solution

By applying a numerical method for the solution of (4), we obtain a numerical solution of the form Comparing the theoretical (5) and the numerical solution (6), we obtain the following two definitions (see [2] for details).

Definition 1. The difference is called dissipation error.

Definition 2. The difference is called dispersion error or phase lag.

Based on the above definitions, it is easy for one to see that we have two new errors, the dissipation and the dispersion errors. These errors can be expressed via Taylor series around an initial value of the frequency. Vanishing the appropriate derivatives of the Taylor series expansion (first, second, etc.) we optimize the specific error.

In order to develop the new method we use the test equation (4). By applying the modified RKN method (2) and (3) to the test equation (4) we obtain the numerical solution where , , , are polynomials in , completely determined by the parameters of the method (2) and (3).

The eigenvalues of the amplification matrix are the roots of the characteristic equation

In phase analysis one compares the phases of with the phases of the roots of the characteristic equation (10). The following definition is originally formulated by van der Houwen and Sommeijer [1].

Definition 3 (phase lag). Apply the RKN method (2) and (3) to the test equation (4). Then we define the phase lag . 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 be dissipative of order .

From Definition 3, it follows that where and are the determinant and the trace of the amplification matrix , respectively. Also it is already known from [1], that we can write the functions and in the following form:

If at a point , , then the Runge-Kutta-Nyström method has zero dissipation at this point.

We can also put forward an alternative definition for the case of infinite order of phase lag.

Definition 4 (phase lag of order infinity). To obtain phase lag of order infinity the relation must hold.

From Definition 4 we have the following theorem.

Theorem 5. If one has phase lag of order infinity and at a point , , then for more details see [7].

Lemma 6. For the construction of a method with nullification of phase lag, amplification error, and their derivatives, one must satisfy the conditions , , , .

4. Derivation of the New Modified RKN Method

The new method that we are going to develop in this section, is a four-stage explicit Runge-Kutta-Nyström method with the FSAL technique (first stage as last), so the method actually uses three stages at each step for the function evaluations. From (2) and (3), the four-stage explicit modified RKN method can be written in the following form: where

By substituting the coefficients that have been used by the DEP algorithm in [15], (14) takes the following form: where

By applying numerical method (16) to the test equation (4), we compute the polynomials , , , in terms of the modified Runge-Kutta-Nyström parameters. From these polynomials we obtain the expressions of and . Then, according to Lemma 6 we solve the system of the four equations(, , , and ) and thus we obtain the expressions of the coefficients , , which are fully dependent on the product of the step length and the frequency :

For small values of the following Taylor series expansions are used:

5. Algebraic Order

In this section we study the algebraic order of the new modified RKN method. We require that and . By using the chain rule and , the Taylor expansion of the exact solution is

By expanding in Taylor series the corresponding numerical solution of a classical explicit Runge-Kutta-Nyström method, where , we have

By equating (20) and (21) with respect to , we obtain the algebraic order conditions for a fourth algebraic order Runge-Kutta-Nyström method. So the equations obtained for are.

Order 2:

order 3:

order 4:

By following the same procedure for , we obtain the corresponding algebraic order conditions, which are.

Order 1:

order 2:

order 3:

order 4:

For the classical RKN method [15] that we use in the present paper, the above conditions are satisfied, as this method is of fourth algebraic order.

In order to verify the algebraic order of the the new modified explicit Runge-Kutta-Nyström method, first we will produce the algebraic order conditions. To do that, we apply the above methodology for the new method. At first, we want to extract the conditions for and so we expand the modified method (14) in Taylor series. The result from the computations is given below:

By equating (20) and (29) with respect to , we obtain the algebraic order conditions of a fourth algebraic order modified Runge-Kutta-Nyström method of the form (14). The equations obtained for are.

Order 2:

order 3:

order 4:

By following the same procedure for the approximate solution of , we obtain the algebraic order conditions for , which are.

Order 1:

order 2:

order 3:

order 4:

As it is proved, in order to construct an explicit modified Runge-Kutta-Nyström method of fourth algebraic order, (30)–(36) must be satisfied.

For the proposed modified RKN method (16) developed in Section 4, we observe that the Taylor expansions for parameters that were obtained are dependent on . What is important for the new modified method is to maintain the fourth algebraic order as the step size takes small values. As already mentioned, is the product of the frequency and the step size ; thus for (30)–(36) should verify the fourth algebraic order of the new modified method.

For it is easy to prove that , . From the last relation, it is extracted that , , and thus (30)–(36) reduced to (22)–(28) which are satisfied for the constant coefficients of DEP algorithm [15].

So it is proved that the new RKN method is of fourth algebraic order. Moreover the local truncation error in is given from the following equation: accordingly the local truncation error in its first derivative is given from the following equation:

At this point, we have already compute the Taylor expansions of(a)the exact solution and the corresponding numerical solution ,(b)the first derivative of the exact solution and the first derivative of the corresponding numerical solution.

The LTE verifies the fourth algebraic order of the new modified method. From the above procedure the local truncation error in is

Respectively, the local truncation error in is

From (30)–(36) and for we prove that the new method is of fourth algebraic order. Moreover the forms (39) and (40) show that the new modified RKN method is of fourth algebraic order, because all the coefficients up to the power of vanish.

6. Numerical Results

6.1. Schrödinger Equation

The one-dimensional or radial Schrödinger equation has the form

We call the term the centrifugal potential, and the function the electric potential. In (41), is a real number denoting the energy, and is a quantum number. The function denotes the effective potential, where and so . The boundary condition , together with a second boundary condition, for large values of , is determined by the physical considerations.

6.1.1. Woods-Saxon Potential with Positive Energies

For the purpose of our numerical illustration we will take the domain of integration as , using the Woods-Saxon potential:

In the case of positive energies , the potential dies away faster than the centrifugal potential , so for a large number for , Schrödinger equation effectively reduces to

(43) has two linearly independent solutions, and , where and are the spherical Bessel and Neumann functions, respectively. When , the solution of (41) takes the following asymptotic form: where is the scattering phase shift that may be calculated from the below formula: where and and both exist in the asymptotic region.

For positive energies and for , we calculate the phase shift and then we compare it with the accurate value which is . The boundary conditions for this eigenvalue problem are and for large .

We use the following eigenenergies:

6.1.2. Woods-Saxon Potential with Negative Energies

In the case of negative energies , we consider the eigenvalue problem with boundary conditions:

In order to solve this problem numerically, by a chosen eigenvalue, we integrate forward from the point and backward from the point , matching up the solution at some internal point in the range of integration.

For this problem we use the following eigenenergies:

6.1.3. Lennard-Jones Potential

In order to investigate the case of , we examine the Lennard-Jones potential which is given by the equation where , for , and step length of integration . We compare the phase shifts with the values found in [16] and present the decimal digits that the approximate solutions that agree with the reference solution.

6.2. Numerical Results

In this section four Runge-Kutta-Nyström methods are compared (including the new method). These methods have four algebraic orders with four stages; also all of them are using the FSAL properties. The methods used in the comparison have been denoted by(i)MRKNDPAF4: the new fourth-order RKN method with four stages (three effective stages with FSAL property), derived in Section 3,(ii)RKNPAF4: the fourth-order RKN method with four stages (three effective stages with FSAL property), phase lag, and amplification error of order infinity of Papadopoulos et al. [7],(iii)DEPRKN4: the high-order method of pair RKN method of Dormand et al. [15],(iv)EFRKN4: the fourth-order exponential fitted RKN method with four stages (three effective stages with FSAL property) of Franco [4].

One way to measure the efficiency of the method is to compute the accuracy in the decimal digits, that is, , when comparing the phase shift to the actual value versus the computational effort measured by the . The range of integration steps that have been used is for the Woods-Saxon potential with positive energies and for the Woods-Saxon potential with negative energies.

In the case of Lennard-Jones potential, we compare the phase shifts with the values found in [16] and present the decimal digits that the approximate solutions that agree with the reference solution. The results are given in Tables 1 and 2. The step length of integration is .

tab1
Table 1: Lennard-Jones potential with .
tab2
Table 2: Lennard-Jones potential with .

The frequency is given by the suggestion of Ixaru and Rizea [17]:

The numerical results were obtained by using the high-level language MATLAB. In Figures 1, 2, 3, 4, 5, 6, 7 and 8 we display the efficiency curves, that is, the accuracy versus the computational cost measured by the number of function evaluations required by each method.

910624.fig.001
Figure 1: Efficiency for the Schrödinger equation using E = 53.588872.
910624.fig.002
Figure 2: Efficiency for the Schrödinger equation using E = 163.215341.
910624.fig.003
Figure 3: Efficiency for the Schrödinger equation using E = 341.495874.
910624.fig.004
Figure 4: Efficiency for the Schrödinger equation using E = 989.701916.
910624.fig.005
Figure 5: Efficiency for the Schrödinger equation using E = − 49.457788728.
910624.fig.006
Figure 6: Efficiency for the Schrödinger equation using E = − 38.122785096.
910624.fig.007
Figure 7: Efficiency for the Schrödinger equation using E = − 22.588602257.
910624.fig.008
Figure 8: Efficiency for the Schrödinger equation using E = − 3.908232481.

Numerical results indicate that the new method derived in Section 4 is very efficient in solving numerically the Schrödinger equation and more accurate than the other methods at all the eigenvalues.

More specifically in the case of Woods-Saxon potential with positive energies, the new method remains more accurate than the RKNPAF4 and EFRKN4 methods up to two decimals at all eigenvalues. Also our method is more accurate than the classical DEPRKN4 method, by two decimals for the eigenvalue , by three decimals for the eigenvalue , and by four decimals for the next two eigenvalues .

In the case of Woods-Saxon potential with negative energies, the new method has almost the same accuracy with the RKNPAF4 and EFRKN4 methods, but it remains (even a little) more accurate at the majority of the negative eigenvalues. In comparison with the DEPRKN4 method, the new method is much more accurate and specifically by one digit for the eigenvalue and by two digits for the rest of the negative eigenvalues.

At last, for the Lennard-Jones potential the new method is more accurate than all the other methods for eigenvalues and .

7. Conclusions

The modified RKN method, developed in this paper, is much more efficient than the classical one, in any case. The new method remained more efficient for all the eigenvalues and in some cases was more accurate than the other methods up to two decimals. Moreover we observe that the accuracy difference between the new method and the other methods increased as the eigenvalue increased (for details about the original proof see [17]).

Disclosure

T. E. Simos is an active member of the European Academy of Sciences and the European Academy of Sciences and Arts and a corresponding member of the European Academy of Arts, Sciences and Humanities.

References

  1. P. J. van der Houwen and B. P. Sommeijer, “Explicit Runge-Kutta-Nyström methods with reduced phase errors for computing oscillating solutions,” SIAM Journal on Numerical Analysis, vol. 24, no. 3, pp. 595–617, 1987. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  2. P. J. van der Houwen and B. P. Sommeijer, “Predictor-corrector methods for periodic second-order initial-value problems,” IMA Journal of Numerical Analysis, vol. 7, no. 4, pp. 407–422, 1987. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  3. H. Van de Vyver, “An embedded phase-fitted modified Runge-Kutta method for the numerical integration of the radial Schrödinger equation,” Physics Letters A, vol. 352, no. 4, pp. 278–285, 2006.
  4. J. M. Franco, “Exponentially fitted explicit Runge-Kutta-Nyström methods,” Journal of Computational and Applied Mathematics, vol. 167, no. 1, pp. 1–19, 2004. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  5. T. E. Simos and J. V. Aguiar, “A modified Runge-Kutta method with phase-lag of order infinity for the numerical solution of the Schrödinger equation and related problems,” Computers and Chemistry, vol. 25, no. 3, pp. 275–281, 2001.
  6. Z. A. Anastassi, D. S. Vlachos, and T. E. Simos, “A family of Runge-Kutta methods with zero phase-lag and derivatives for the numerical solution of the Schrödinger equation and related problems,” Journal of Mathematical Chemistry, vol. 46, no. 4, pp. 1158–1171, 2009. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  7. D. F. Papadopoulos, Z. A. Anastassi, and T. E. Simos, “A modified phase-fitted and amplification-fitted Runge-Kutta-Nyström method for the numerical solution of the radial Schrödinger equation,” Journal of Molecular Modeling, vol. 16, no. 8, pp. 1339–1346, 2010.
  8. D. F. Papadopoulos and T. E. Simos, “A new methodology for the construction of optimized Runge-Kutta-Nyström methods,” International Journal of Modern Physics C, vol. 22, no. 6, pp. 623–634, 2011. View at Publisher · View at Google Scholar · View at MathSciNet
  9. Z. Kalogiratou, T. Monovasilis, and T. E. Simos, “Computation of the eigenvalues of the Schrödinger equation by exponentially-fitted Runge-Kutta-Nyström methods,” Computer Physics Communications, vol. 180, no. 2, pp. 167–176, 2009. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  10. I. Alolyan and T. E. Simos, “High algebraic order methods with vanished phase-lag and its first derivative for the numerical solution of the Schrödinger equation,” Journal of Mathematical Chemistry, vol. 48, no. 4, pp. 925–958, 2010. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  11. I. Alolyan and T. E. Simos, “A family of ten-step methods with vanished phase-lag and its first derivative for the numerical solution of the Schrödinger equation,” Journal of Mathematical Chemistry, vol. 49, no. 9, pp. 1843–1888, 2011. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  12. A. Konguetsof, “Two-step high order hybrid explicit method for the numerical solution of the Schrödinger equation,” Journal of Mathematical Chemistry, vol. 48, no. 2, pp. 224–252, 2010. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  13. A. Konguetsof, “A hybrid method with phase-lag and derivatives equal to zero for the numerical integration of the Schrödinger equation,” Journal of Mathematical Chemistry, vol. 49, no. 7, pp. 1330–1356, 2011. View at Publisher · View at Google Scholar · View at MathSciNet
  14. A. A. Kosti, Z. A. Anastassi, and T. E. Simos, “An optimized explicit Runge-Kutta method with increased phase-lag order for the numerical solution of the Schrödinger equation and related problems,” Journal of Mathematical Chemistry, vol. 47, no. 1, pp. 315–330, 2010. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  15. J. R. Dormand, M. E. A. El-Mikkawy, and P. J. Prince, “Families of Runge-Kutta-Nyström formulae,” IMA Journal of Numerical Analysis, vol. 7, no. 2, pp. 235–250, 1987. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  16. H. Kobeissi, K. Fakhreddine, and M. Kobeissi, “On a canonical functions approach to the elastic scattering phase-shift problem,” International Journal of Quantum Chemistry, vol. 40, no. 1, pp. 11–21, 1991.
  17. L. G. Ixaru and M. Rizea, “A Numerov-like scheme for the numerical solution of the Schrödinger equation in the deep continuum spectrum of energies,” Computer Physics Communications, vol. 19, no. 1, pp. 23–27, 1980.