International Journal of Mathematics and Mathematical SciencesVolume 2011 (2011), Article ID 328197, 17 pagesdoi:10.1155/2011/328197

Research Article

A Variable Step-Size Exponentially Fitted Explicit Hybrid Method for Solving Oscillatory Problems

F. Samat, F. Ismail, and M. B. Suleiman

Department of Mathematics, Faculty of Science, Universiti Putra Malaysia (UPM), Selangor, 43400 Serdang, Malaysia

Received 24 May 2011; Revised 9 September 2011; Accepted 9 September 2011

Academic Editor: Naseer Shahzad

Copyright © 2011 F. Samat 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 exponentially fitted explicit hybrid method for solving oscillatory problems is obtained. This method has four stages. The first three stages of the method integrate exactly differential systems whose solutions can be expressed as linear combinations of , while the last stage of this method integrates exactly systems whose solutions are linear combinations of . This method is implemented in variable step-size code basing on an embedding approach. The stability analysis is given. Numerical experiments that have been carried out show the efficiency of our method.

1. Introduction

There has been great interest in the research of new methods for numerically solving the second-order initial value problems of the form whose solution exhibits specific oscillatory behavior. Such problems arises in celestial mechanics, in quantum mechanical scattering problems, and elsewhere, and they can be solved by using general purpose methods or by using codes specially adapted to the structure or to the solution of the problem. In the case of adapted hybrid methods, particular algorithms have been proposed by several authors (e.g., see [1–3]) to solve these classes of problems.

This paper concerns derivation of a new method for the numerical integration of (1.1) with oscillatory solution through the usage of an exponential-fitting technique proposed by Vanden Berghe and Van Daele [4]. The coefficients of the new method depend on the product of frequency and step size, thus the method can only be applied when a good estimate of the dominant frequency of the solution is known in advance. The exponential-fitting technique is generally described in the context of linear multistep methods as follows.

Consider linear -step methods of the form where , , and do not both vanish. Associate with the method is the linear operator where . The linear operator is said to integrate exactly the function if . We now impose the requirement that the operator integrates exactly the following functions: with . Then, we solve the algebraic equations for and , . In particular, the step-by-step exponential-fitting procedure is described as follows.

Step 1. Find the maximal such that integrates exactly the set of power functions and that the resulting coefficients have the same value as those of a classical method. The classical method is the method with constant coefficients.

Step 2. Consider this set of functions with . Using obtained from Step 1, we find suitable and which correspond to a finite set consisting of power functions and exponentials. Then, solve equations obtained by assuming that integrates exactly elements of this set.

Considering each formula stage of a hybrid method as a linear multistep method with nonequidistant grid points; we derive our exponentially fitted hybrid method. This approach has also been employed by Simos and Vigo-Aguiar [1] and D’Ambrosio et al. [5] in their derivation of exponential-fitting methods.

2. Embedded Pairs of Hybrid Methods

For the numerical solution of (1.1), we consider the explicit hybrid method which has been established by Franco [6]: where and represent and , respectively. The method requires function evaluations or stages at each step of integration. The explicit hybrid method can be represented by the Butcher tableau: Adopting the concept of embedded pair of Runge-Kutta-Nystrom (RKN) methods (see [7, 8]), we define an embedded pair of hybrid methods to be based on the hybrid method () of order and another hybrid method () of order , represented by the following tableau: Embedded pairs of explicit hybrid methods are used in variable step-size algorithm because they provide cheap error estimation. A local error estimation is determined by the formula where and are solutions obtained using the higher-order and the lower-order formula, respectively. The LTE is used to control the step-size of which the procedure is given by

(i) if , then ,

(ii) if , then ,

(iii) if , then and repeat the step,

where, from numerical experiments, div is chosen to be 2¹⁷. We do not allow step-size change after each step because it would contribute to unnecessary rounding errors. If the step is acceptable (i.e., tol/div < LTE < div · tol and LTE ≤ tol/div), then we adopt the widely used of performing local extrapolation: although the LTE is the error estimation for the lower-order formula, the solutions obtained by using the higher-order formula are actually accepted at each point.

3. Construction of Hybrid Methods

3.1. The Sixth-Order Hybrid Method of Franco

In this section, we consider the sixth-order hybrid method which is zero dissipative and dispersive of order six derived in [6]. The interval of periodicity of the method is (0, 2.75). In order for the method to be implemented in variable step-size code as described in Section 2, we derive a three-stage fourth-order explicit hybrid method. Order conditions for the fourth-order hybrid methods with as given in [9] are listed as follows.

Order conditions: The two-step explicit hybrid methods (2.1) are a subclass of the two-step hybrid methods considered by Coleman [9] by taking ,  , and   . Substituting , ,  , and    together with coefficients and of the sixth-order formula into (3.1), we get the following values of coefficients for the three-stage fourth-order method: The embedded 6(4) pair of hybrid methods is denoted as EHM6(4) method and represented by the tableau:

3.2. Variable Step-Size Exponentially Fitted Explicit Hybrid Method

We consider the four-stage explicit hybrid method given by the following tableau: Note that the -values and some of the -values are taken from the sixth-order hybrid method of Franco as described briefly in Section 3.1. The formula of a four-stage explicit hybrid method is Associate with (3.6), (3.7), and (3.8) the linear operator , we have Note that the operator integrates exactly 1 and for any and . The coefficients of the variable step-size exponentially fitted hybrid method are obtained by solving equations arising from the choice of integers , , and described as follows.

Step 1. Consider the linear operator that associated with (3.6). Solving equations obtained by assuming that and imposing to integrate exactly and , we obtain and which is the same value as the coefficients of EHM6(4) method. This means that . Choosing ,   gives us . Setting and imposing to integrate exactly leads the following equations: Solving these equations, we get coefficients and .

Step 2. Consider the linear operator that associated with (3.7). Solving equations obtained by assuming that ,  ,   and imposing to integrate exactly and , we obtain and which are of the same value as the coefficients of EHM6(4) method. This means that . Choosing , gives us . Setting ,  ,  , and , we get the following equations: Solving these equations, we get and .

Step 3. Consider the linear operator that associated with (3.8). Solving equations obtained by assuming that ,  ,  ,   and imposing to integrate exactly and , we obtain and which are of the same value as the coefficients of EHM6(4) method. This means that . Choosing and gives us . Setting ,  ,  ,  ,  , and leads to the following equations: Solving these equations, we get and .

Step 4. The linear operator associating with (3.9) is given as where . Solving equations obtained by assuming that ,  ,   and imposing to integrate exactly , , , and , we obtain which are of the same value as the coefficients of EHM6(4) method. This means that . Choosing and gives us . Setting ,  ,  , and , the following equations are obtained: Solving these equations, we get , , , , and .

Step 5. Consider Associate with (3.17) is the linear operator : where . Solving equations obtained by assuming that ,   and imposing to integrate exactly , and , we obtain which are of the same value as the coefficients of EHM6(4) method. This means that . Choosing and gives us . Setting , and , the following equations are obtained: Solving these equations we get , , and .

Let . The Taylor series expansions for coefficients of the variable step-size exponentially fitted explicit hybrid method: It is noted that, if then the coefficients of EHM6(4) method will be recovered. This is the important feature of the new method which will be denoted by EEHM6(4) method. In our program code, we use as an imaginary number and thus, can be written as . We first convert all coefficients (which are in term of ) to trigonometric form by using computer algebra package Maple. Then, we substitute into the resulting expressions. For small , the coefficients are subject to heavy cancellations. Since varies throughout the integration, therefore it is convenient to use Taylor series expansions for the coefficients of the method.

4. Stability

Now our interest lies in the stability analysis for an exponentially fitted explicit hybrid method of the form (2.1) in which its coefficients are functions of one frequency and the step-length . Applying this method to the test problem we obtain which is written in vector form where , , and . Solving for in (4.2) and then substitute it into (4.3) leads to the following recursion: where and .

The characteristic equation associated with (4.4) is The subsequent definitions are introduced to provide the stability concept for exponentially fitted hybrid methods corresponding to the characteristic equation (4.5). Here, we follow the main ideas given by Coleman and Ixaru [3].

Definition 4.1. For exponentially fitted hybrid methods corresponding to the characteristic equation (4.5) and satisfying , a region of stability Ω is a region of the - plane throughout which . Any closed curve defined by is a stability boundary.

Definition 4.2. For exponentially fitted hybrid methods corresponding to the characteristic equation (4.5) and satisfying , a region of stability Ω is a region of the - plane such that, for every , Any closed curve defined by and is a stability boundary.In Figure 1, we illustrate the stability region of the exponentially fitted method which is based on the higher-order formula of EHM6(4) method whereas, in Figure 2, we show the stability region of the exponentially fitted method which is based on the lower-order formula of EHM6(4) method. From Figures 1 and 2, it is observed that the stability regions get smaller as increases.

(a)

(b)

(c)

Figure 1: Stability region of the exponentially fitted method which is based on the higher-order formula of EHM6(4) method. (a) (, ), (b) (, ), (c) (, ).

(a)

(b)

(c)

Figure 2: Stability region of the exponentially fitted method which is based on the lower-order formula of EHM6(4) method. (a) (, ), (b) (, ), (c) (, ).

5. Numerical Experiments

The following are abbreviations representing the codes that have been used in the comparisons.(i)EEHM6(4): Variable step-size exponentially fitted explicit hybrid method derived in this paper. This method requires four function evaluations at each step.(ii)EHM6(4): Embedded explicit hybrid method 6(4) pair derived in this paper. This method requires four function evaluations at each step.(iii)FRKN4(3): Modified embedded explicit RKN 4(3) pair proposed by Van de Vyver [10]. This method requires three function evaluations at each step.

We have applied all the above codes to four problems to evaluate their efficiency. For controlling the step-size in FRKN4(3) code, we use the step-size control procedure discussed in [10]. In Figures 3, 4, 5, and 6, the efficiency curves of natural logarithm of maximum global error (MAXGE) versus the computational effort measured by the natural logarithm of number of function evaluations (NFE) are depicted. The maximum global error (MAXGE) is given by this formula where is the true solution and is the computed solution. In the following test problems, the tolerances used are , . The exact solutions are used to obtain all the starting values.

Figure 3: Efficiency curves for Problem 1.

Figure 4: Efficiency curves for Problem 2.

Figure 5: Efficiency curves for Problem 3.

Figure 6: Efficiency curves for Problem 4.

Problem 1 (perturbed system). Source: Franco [11]: with and Solution: For EEHM6(4) and FRKN4(3) codes, we choose for the first component whereas for the second component.

Problem 2 (linear oscillatory problem). Source: Franco [6]: with the theoretical solution: For EEHM6(4) and FRKN4(3) codes, we choose .

Problem 3 (perturbed system). Source: Hans Van de Vyver [10]: with and Solution: For EEHM6(4) and FRKN4(3) codes, we choose .

Problem 4 (the undamped Duffing equation). Source: Van de Vyver [10]: where and . The exact solution computed by the Galerkin method with a precision of the coefficients is given by where , , , , and . For EEHM6(4) and FRKN4(3) codes, we choose .

From Figures 3–5, it is obvious that EEHM6(4) performs more efficiently than EHM6(4) and FRKN4(3) codes for solving Problems 1–3. For Problem 4, EEHM6(4) is as efficient as EHM6(4) code.

6. Conclusions

The embedding approach for hybrid methods has been devised by adopting the embedding concept in RKN methods. By using this embedding approach, we present Franco’s sixth-order explicit hybrid method for variable step-size code, denoted as EHM6(4) method. Based on EHM6(4) method, the new variable step-size exponentially fitted hybrid method has been developed and it is denoted as EEHM6(4) method. The numerical results show that EEHM6(4) method is efficient for the numerical solution of oscillatory problems. Furthermore, EEHM6(4) method is preferable to handle perturbed oscillators.

References

T. E. Simos and J. Vigo-Aguiar, “A dissipative exponentially-fitted method for the numerical solution of the Schrödinger equation and related problems,” Computer Physics Communications, vol. 152, no. 3, pp. 274–294, 2003.

T. E. Simos and P. S. Williams, “P-stable hybrid exponentially-fitted method for the numerical integration of the Schrodinger equation,” Computer Physics Communications, vol. 131, no. 1, pp. 109–119, 2000.

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.

G. Vanden Berghe and M. Van Daele, “Trigonometric polynomial or exponential fitting approach?” Journal of Computational and Applied Mathematics, vol. 233, no. 4, pp. 969–979, 2009.

R. D'Ambrosio, E. Esposito, and B. Paternoster, “Exponentially fitted two-step hybrid methods for ,” Journal of Computational and Applied Mathematics, vol. 235, no. 16, pp. 4888–4897, 2011.

J. M. Franco, “A class of explicit two-step hybrid methods for second-order IVPs,” Journal of Computational and Applied Mathematics, vol. 187, no. 1, pp. 41–57, 2006.

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.

S. Filippi and J. Graf, “New Runge-Kutta-Nyström formula-pairs of order 8(7), 9(8), 10(9) and 11(10) for differential equations of the form ,” Journal of Computational and Applied Mathematics, vol. 14, no. 3, pp. 361–370, 1986.

J. P. Coleman, “Order conditions for a class of two-step methods for ,” IMA Journal of Numerical Analysis, vol. 23, no. 2, pp. 197–220, 2003.

H. Van de Vyver, “A Runge-Kutta-Nyström pair for the numerical integration of perturbed oscillators,” Computer Physics Communications, vol. 167, no. 2, pp. 129–142, 2005.

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.