Mathematical Problems in Engineering

Volume 2012 (2012), Article ID 867948, 12 pages

doi:10.1155/2012/867948

Research Article

An Optimized Runge-Kutta Method for the Numerical Solution of the Radial Schrödinger Equation

Qinghe Ming, Yanping Yang, and Yonglei Fang

School of Mathematics and Statistics, Zaozhuang University, Zaozhuang 277160, China

Received 20 July 2012; Accepted 30 August 2012

Academic Editor: Gradimir Milovanovic

Copyright © 2012 Qinghe Ming et al. 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.

An optimized explicit modified Runge-Kutta (RK) method for the numerical integration of the radial Schrödinger equation is presented in this paper. This method has frequency-depending coefficients with vanishing dispersion, dissipation, and the first derivative of dispersion. Stability and phase analysis of the new method are examined. The numerical results in the integration of the radial Schrödinger equation with the Woods-Saxon potential are reported to show the high efficiency of the new method.

1. Introduction

In this paper, we are concerned with the numerical integration of the one-dimensional Schrödinger equation of the form

where the real number is the energy and the function is the effective potential satisfying as . Two boundary conditions are associated with this equation: one is , and the other imposed at large is determined by physical considerations. The form of this second boundary condition depends crucially on the sign of the energy . Such problems are frequently encountered in a variety of scientific fields and engineering applications [1–9]. Concerning the oscillatory character of the solution to the Schrödinger equation (1.1), there have appeared a lot of numerical integrators of adapted type, a pronounced class of which is based on important properties such as the phase lag and the amplification (see [10–18]). These are actually two different kinds of truncation errors. The first is the angle between the analytical solution and the numerical solution, and the second is the distance from a standard cyclic solution. If a good frequency is estimated in advance, then it is a good choice to construct numerical methods with zero dispersion or/and zero dissipation. These techniques are called phase fitted or/and zero dissipation. Related work can be founded in [19–21]. For Runge-Kutta methods, Simos and Aguiar [18] constructed a modified Runge-Kutta method for the numerical integration of the Schrödinger equation by phase fitting based on the fifth-order RK method. Recently, Van de Vyver [16] gave an embedded pair of modified RK methods by nullifying the phase-lags of the fifth-order method and the fourth-order method. And in [22], Tsitouras and Simos constructed phase-fitted and zero dissipation fifth-order Runge-Kutta method for the numerical solution of oscillatory problems.

In this paper, inspired by the ideas in [23–28], we construct a new kind of modified fifth-order Runge-Kutta method by nullifying the dispersion, the dissipation, and the first derivative of the dispersion. In Section 2, the preliminaries of the phase properties of explicit modified Runge-Kutta methods are introduced. In Section 3, the coefficients of a new kind of optimized modified RK method are obtained. Section 4 examines the stability and phase properties of the new method. In Section 5, the numerical experiments are reported.

2. Preliminaries

We begin by considering the numerical integration of the initial value problem (IVP) of first-order differential equations in the following form:

whose solution shares an oscillatory character. We follow the convention to assume that the frequency is known to be in advance or can be accurately estimated. An -stage-modified explicit Runge-Kutta (RK) method has the following scheme:

where the coefficients ,  ,  ,   are constants, is the step size, and the parameters are even functions of . It is convenient to express the modified RK method (2.2) by the Butcher tableau as follows:

or simply by (). The extra-frequency-depending parameters ,  ,   are introduced to tune the traditional RK method to the special oscillatory structure of the problem. We assume that ,   so that as , the modified RK method (2.2) reduces to a traditional RK method. An alternative approach adopted by, for example, exponential/trigonometric fitting techniques, is to let some of the coefficients ,  ,  ,   be functions of (see [16, 18, 29]).

Applying the modified RK method (2.2) to the test equation as follows:

yields

A comparison of the numerical solution with the exact solution leads to the notions of phase-lag and dissipation error defined as follows.

Definition 2.1. The following two quantities are called the phase lag (or dispersion) and the amplification factor error (or dissipation error), respectively:

The method is said to be dispersive of order and dissipative of order if

If and , the method is called phase fitted (zero dispersive) and amplification-fitted (zero dissipative), respectively.

For modified RK method (2.2), we have

where

are polynomials in , which are completely defined by the Runge-Kutta coefficients ,  ,  , and . Therefore, we have

Based on the fifth algebraic order six-stage Dormand and Prince Runge-Kutta method, Simos and Aguiar [18] obtained an explicit modified RK method with one parameter (taking the orthers for ) determined by nullifying the quantity . In [22], Tsitouras and Simos presented an optimized Runge-Kutta method by nullifying the dispersion and the dissipation. In this paper, we construct a new optimized Runge-Kutta method by nullifying the dispersion, the dissipation, and the first derivative of the dispersion.

3. Construction of the New Method

In this section, we are concerned with the following Runge-Kutta method given by the Butcher tableau as follows:

If we choose , the classical Runge-Kutta method with order fifth derived by Dormand and Prince [30] is recovered.

In order to construct the new embedded RK pair, we set free and keep the rest of the coefficients. Motivated by the ideas in [23–28], we obtain the dispersion, the dissipation, and the first derivative of the dispersion of this method, which depend on as follows:

where

Now, solving (3.2), we get values in terms of . Instead of giving the very complicated expressions for , for the purpose of practical computation, we present their Taylor expansions as follows:

In order to check the algebraic order of the newly obtained modified RK method, we note that the order conditions listed in [31] for traditional RK methods are not sufficient for the modified RK method (2.2). Writing

we obtain the following additional conditions for the modified RK method (2.2) to be of up to order five (see [16]):

(i)order 3 requires:

(ii)order 4 requires in addition:

(iii)order 5 requires in addition:

By simple calculation, it is verified that the new method is of algebraic order fifth. We denote the new method as MODRK5PLDPLAM.

4. Analysis of Stability and Phase Properties

In this section, we are interested in the stability and phase properties of the new method. Lambert and Watson’s stability theory [32] was reformulated by Coleman and Ixaru [33] for the periodicity of exponentially fitted symmetric methods for . Van de Vyver [34] adapted this theory to RK methods. Following Van de Vyver’s approach, we consider the test equation as follows:

Applying the modified RK method (2.2) to test (4.1) yields the difference equation

where

with the identity matrix.

Definition 4.1 (see [34]). For the modified RK method (2.2) with stability function , the region in the - plane

is called the region of imaginary stability. And any closed curve defined by is a stability boundary of the method.

In Figure 1 we plot the region of imaginary stability for the method MODRK5PLDPLAM.

Figure 1: Regions of imaginary stability of the MODDPHARK5 method.

Definition 4.2 (see [34]). For the modified RK method (2.2) with stability function , the quantities

are called the phase lag (dispersion) and amplification factor error (dissipation), respectively. If

the method is said to be of phase-lag order and dissipation order , respectively, where the and are called the phase-lag constant and dissipation constant, respectively.

We note that, by definition, when (), it must be true that and . In general, since the fitting frequency is just an estimate of the true frequency. Therefore the order of and in Definition 4.2 measure to what extent a modified RK method is accurate in phase and dissipation. Denoting the ratio , we obtain the following expressions for the phase lag and the dissipation error of the new method MODRK5PLDPLAM:

Thus, the method MODRK5PLDPLAM has a phase lag of order six and a dissipation of order five.

5. Numerical Experiments

In this section, we test the numerical performance of the new fifth-order method in the integration of the radial Schrödinger equation with the well-known Woods-Saxon potential, respectively. We compare the new method with some existing highly efficient methods in the literature.

The methods we choose for comparison are as follows:

(i)PHARK5S: the phase-fitted fifth-order RK method given by Simos in [17],(ii)MODPHARK5S: the modified phase-fitted fifth-order RK method given by Simos and Aguiar in [18],(iii)MODPHARK5V: the higher-order method of the modified phase-fitted embedded RK5 (2.4) pair given by Van de Vyver in [16],(iv)ARK5: an adapted fifth-order RK method given by Fang et al. in [35],(v)PHADISRK5S: the phase-fitted and zero dissipation fifth-order RK method given by Tsitouras and Simos in [22],(vi)MODRK5PLDPLAM: the phase-fitted fifth-order method derived in this paper.

We consider the numerical integration of the Schrödinger equation (1.1) with the well-known Woods-Saxon potential

where . The problem is solved in the interval . Following [16, 36–38], we choose the fitting frequency

In the numerical experiment we consider the resonance problem (), the numerical results are compared with the analytical solution of the Woods-Saxon potential, rounded to six decimal places. In Figures 2, 3, 4, and 5, we plot the error versus (with the integration step-size ) for , 163.215341, 341.495874, and 989.701916, respectively.

Figure 2: Efficiency curves for .

Figure 3: Efficiency curves for .

Figure 4: Efficiency curves for .

Figure 5: Efficiency curves for .

6. Conclusions and Discussions

Based on the classical fifth RK method of Dormand and Prince [30], a new optimized explicit modified RK method with modifying parameters is obtained by nullifying the dispersion, the dissipation, and the first derivative of the dispersion. The numerical results stated in Figures 2–5 illustrate the higher efficiency of the new method compared to some highly efficient methods in the recent literature [16–18, 22, 35].

Acknowledgments

The authors are deeply grateful to the anonymous referees for their constructive comments and valuable suggestions. This research is partially supported by NSFC (no. 11101357), the foundation of Shandong Outstanding Young Scientists Award Project (no. BS2010SF031), the foundation of Scientific Research Project of Shandong Universities (no. J11LG69), and NSF of Shandong Province, China (no. ZR2011AL006).

References

A. C. Allison, “The numerical solution of coupled differential equations arising from the Schrödinger equation,” Journal of Computational Physics, vol. 6, pp. 378–391, 1970.

J. M. Blatt, “Practical points concerning the solution of the Schrödinger equation,” Journal of Computational Physics, vol. 1, pp. 378–391, 1967.

J. W. Cooley, “An improved eigenvalue corrector formula for solving the Schrödinger equation for central fields,” Mathematics of Computation, vol. 15, pp. 363–374, 1961.

G. Avdelas, T. E. Simos, and J. Vigo-Aguiar, “An embedded exponentially-fitted Runge-Kutta method for the numerical solution of the Schrödinger equation and related peri-odic initial-value problems,” Computer Physics Communications, vol. 131, pp. 52–67, 2000.

J. Vigo-Aguiar and T. E. Simos, “Review of multistep methods for the numerical solution of the radial Schrödinger equation,” International Journal of Quantum Chemistry, vol. 103, pp. 278–290, 2005.

J. Vigo-Aguiar, J. Martín-Vaquero, and R. Criado, “On the stability of exponential fitting BDF algorithms,” Journal of Computational and Applied Mathematics, vol. 175, no. 1, pp. 183–194, 2005.

T. E. Simos and J. Vigo-Aguiar, “A new modified Runge-Kutta-Nystrom method with phase-lag of order infinity for the numerical solution of the Schrödinger equation and related problems,” International Journal of Modern Physics C, vol. 11, pp. 1195–1208, 2000.

T. E. Simos, “Exponentially and trigonometrically fitted methods for the solution of the Schrödinger equation,” Acta Applicandae Mathematicae, vol. 110, no. 3, pp. 1331–1352, 2010.

Z. Kalogiratou, Th. Monovasilis, and T. E. Simos, “Symplectic integrators for the numerical solution of the Schrödinger equation,” Journal of Computational and Applied Mathematics, vol. 158, no. 1, pp. 83–92, 2003.

T. E. Simos, I. T. Famelis, and C. Tsitouras, “Zero dissipative, explicit Numerov-type methods for second order IVPs with oscillating solutions,” Numerical Algorithms, vol. 34, no. 1, pp. 27–40, 2003.

T. E. Simos, “Dissipative trigonometrically-fitted methods for linear second-order IVPs with oscillating solution,” Applied Mathematics Letters, vol. 17, no. 5, pp. 601–607, 2004.

K. Tselios and T. E. Simos, “Runge-Kutta methods with minimal dispersion and dissipation for problems arising from computational acoustics,” Journal of Computational and Applied Mathematics, vol. 175, no. 1, pp. 173–181, 2005.

Z. A. Anastassi and T. E. Simos, “An optimized Runge-Kutta method for the solution of orbital problems,” Journal of Computational and Applied Mathematics, vol. 175, no. 1, pp. 1–9, 2005.

S. Stavroyiannis and T. E. Simos, “Optimization as a function of the phase-lag order of nonlinear explicit two-step P-stable method for linear periodic IVPs,” Applied Numerical Mathematics, vol. 59, no. 10, pp. 2467–2474, 2009.

T. E. Simos, “A family of fifth algebraic order trigonometricallyffitted Runge-Kutta methods for the numerical solution of the Schrödinger equation,” Computational Materials Science, vol. 34, pp. 342–354, 2005.

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, pp. 278–285, 2006.

T. E. Simos, “An embedded Runge-Kutta method with phase-lag of order infinity for the numerical solution of the Schrödinger equation,” International Journal of Modern Physics C, vol. 11, no. 6, pp. 1115–1133, 2000.

T. E. Simos and J. V. Aguiar, “A modified phase-fitted Runge-Kutta method for the numerical solution of the Schrödinger equation,” Journal of Mathematical Chemistry, vol. 30, no. 1, pp. 121–131, 2001.

A. D. Raptis and T. E. Simos, “A four-step phase-fitted method for the numerical integration of second order initial value problems,” BIT. Numerical Mathematics, vol. 31, no. 1, pp. 160–168, 1991.

B. Paternoster, “A phase-fitted collocation-based Runge-Kutta-Nyström method,” Applied Numerical Mathematics, vol. 35, no. 4, pp. 339–355, 2000.

T. E. Simos, “A two-step method with phase-lag of order infinity for the numerical integration of second order periodic initial-value problems,” International Journal of Computer Mathematics, vol. 39, pp. 135–140, 1991.

Ch. Tsitouras and T. E. Simos, “Optimized Runge-Kutta pairs for problems with oscillating solutions,” Journal of Computational and Applied Mathematics, vol. 147, no. 2, pp. 397–409, 2002.

A. A. Kosti, Z. A. Anastassi, and T. E. Simos, “Construction of an optimized explicit Runge-Kutta-Nyström method for the numerical solution of oscillatory initial value problems,” Computers & Mathematics with Applications, vol. 61, no. 11, pp. 3381–3390, 2011.

A. A. Kosti, Z. A. Anastassi, and T. E. Simos, “An optimized explicit Runge-Kutta-Nyström method for the numerical solution of orbital and related periodical initial value problems,” Computer Physics Communications, vol. 183, no. 3, pp. 470–479, 2012.

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.

I. Alolyan and T. E. Simos, “A new hybrid two-step method with vanished phase-lag and its first and second derivatives for the numerical solution of the Schrödinger equation and related problems,” Journal of Mathematical Chemistry, vol. 50, no. 7, pp. 1861–1881, 2012.

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.

I. Alolyan and T. E. Simos, “A family of eight-step methods with vanished phase-lag and its derivatives for the numerical integration of the Schrödinger equation,” Journal of Mathematical Chemistry, vol. 49, no. 3, pp. 711–764, 2011.

G. Vanden Berghe, H. De Meyer, M. Van Daele, and T. Van Hecke, “Exponentially-fitted explicit Runge-Kutta methods,” Computer Physics Communications, vol. 123, no. 1–3, pp. 7–15, 1999.

J. R. Dormand and P. J. Prince, “A family of embedded Runge-Kutta formulae,” Journal of Computational and Applied Mathematics, vol. 6, no. 1, pp. 19–26, 1980.

E. Hairer, S. P. Nørsett, and G. Wanner, Solving Ordinary Differential Equations I, Nonstiff Problems, Springer, Berlin, Germany, 2nd edition, 1993.

J. D. Lambert and I. A. Watson, “Symmetric multistep methods for periodic initial value problems,” Journal of the Institute of Mathematics and its Applications, vol. 18, no. 2, pp. 189–202, 1976.

J. P. Coleman and L. Gr. Ixaru, “P-stability and exponential-fitting methods for ,” IMA Journal of Numerical Analysis, vol. 16, no. 2, pp. 179–199, 1996.

H. Van de Vyver, “Stability and phase-lag analysis of explicit Runge-Kutta methods with variable coefficients for oscillatory problems,” Computer Physics Communications, vol. 173, no. 3, pp. 115–130, 2005.

Y. L. Fang, Y. Z. Song, and X. Y. Wu, “New embedded pairs of explicit Runge-Kutta methods with FSAL properties adapted to the numerical integration of oscillatory problems,” Physics Letters A, vol. 372, no. 44, pp. 6551–6559, 2008.

H. Van de Vyver, “Comparison of some special optimized fourth-order RungeCKutta methods for the numerical solution of the Schrödinger equation,” Computer Physics Communications, vol. 166, no. 2, Article ID 109C122, 2005.

H. Van de Vyver, “Modified explicit Runge-Kutta methods for the numerical solution of the Schrödinger equation,” Applied Mathematics and Computation, vol. 171, no. 2, pp. 1025–1036, 2005.

L. Gr. 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 Communication, vol. 19, pp. 23–27, 1980.