key: cord-0058890-id0pt7o9
authors: Conte, Dajana; D’Ambrosio, Raffaele; Giordano, Giuseppe; Ixaru, Liviu Gr.; Paternoster, Beatrice
title: User-Friendly Expressions of the Coefficients of Some Exponentially Fitted Methods
date: 2020-08-24
journal: Computational Science and Its Applications - ICCSA 2020
DOI: 10.1007/978-3-030-58799-4_4
sha: 51f3dc4ab38ec3168f9e3f4d9749b7594af12abd
doc_id: 58890
cord_uid: id0pt7o9

The purpose of this work consists in reformulating the coefficients of some exponentially-fitted (EF) methods with the aim of avoiding numerical cancellations and loss of precision. Usually the coefficients of an EF method are expressed in terms of [Formula: see text] , where [Formula: see text] is the frequency and h is the step size. Often, these coefficients exhibit a 0/0 indeterminate form when [Formula: see text] . To avoid this feature we will use two sets of functions, called C and S, which have been introduced by Ixaru in [61]. We show that the reformulation of the coefficients in terms of these functions leads to a complete removal of the indeterminacy and thus the convergence of the corresponding EF method is restored. Numerical results will be shown to highlight these properties.

Exponential fitting is a mathematical procedure to generate numerical methods for different problems with a pronounced oscillatory or hyperbolic behaviour, usually occurring in interpolation, numerical differentiation and quadrature [30, 31, 33, 35, 60, 64, 65, 70, 78] , numerical solution of first order ordinary differential equations [3, 4, 27, 40, 43, 49, 60, 68, 72, 73, 75, 76, 79] , second order differential equations [54, 63] , integral equations [15, 16] , fractional differential equations [1] , partial differential equations [14, 45, 46, 50, 52] . This procedure has been introduced in [62] . Its central idea consists in determining the coefficients of a numerical method by asking that the method is exact for the following set of functions, which is called a fitting space:

where μ may be real or imaginary. The coefficients are functions of the parameter ν = ωh, where ω is the frequency of the oscillatory or hyperbolic functions and h is the step size. The values of μ to be used in (1) are the imaginary μ = iω and real μ = ω, respectively. Often, these coefficients exhibit the indeterminate form 0/0 when ν → 0 such that, in order to restore the convergence of the corresponding numerical methods when ν is small (in practice this depends on how small is the step size h), it is necessary to make use of the Taylor series of the coefficients. Expressed in different words, an accurate computation of the EF coefficients requires the knowledge of four different formulas (an analytic formula valid for big ν and a power series for small ν, for each of the trigonometrical or hyperbolic fitting). In the paper [61] a method was described to replace the four formulas by a single one. The coefficients have been expressed in terms of two sets of particular functions, called C(Z) and S(Z), where Z = ±ν 2 , for real and imaginary μ. A similar method has been introduced in the paper [30] , in which the coefficients are expressed in terms of η m (Z) functions. The work is organized as follows. In Sect. 2 we recall the two sets of C and S functions, and in Sect. 3 the general procedure for the conversion of the coefficient in terms of C and S functions is briefly presented. In Sect. 4 we reformulate the coefficients for the methods in [67, 69] . In Sect. 5 numerical experiments are presented to show how the converted coefficients restore the convergence of the method.

The original ν = μh is replaced by the new Z = ν 2 which is negative if ν is imaginary and positive when ν is real. Thus, Z < 0 and Z > 0 cover the trigonometric and hyperbolic case, respectively. To define the sets of functions C and S we rely on the family of functions η m (Z) functions, m = −1, 0, 1, ..., . . .:

and, for Z = 0,

while for Z = 0,

The two sets of C and S functions are defined as follows: C −1 (Z) and S −1 (Z) are given by the first two η functions,

while the next ones are derived by recurrence for Z = 0,

and by the following values at Z = 0,

, for any n =, 0, 1, 2, . . . .

An important property is the reverse relations:

for n = −1, 0, 1, . . .. For an accurate computation of these functions, it is necessary to introduce their series expansions:

Note: We acknowledge with thanks a recent private communication by Prof. Ander Murua that sets C and S are directly related to the Stumpff functions c n (Z)(n = 0, 1, 2, 3, ..):

For the Stumpff functions see [81] and references therein.

Now we describe the procedure introduced in [61] for the conversation of the coefficients. Let σ(z), a generic coefficient derived by EF technique:

where N (z) and D (z) contains trigonometrical or hyperbolic functions and tend to 0 when z → 0. Let us denote generically by F any of numerator N or denominator D, and treated separately these two functions. The procedure has two steps.

In the first step, F (z) is expressed in terms of Z in the following way:

where k = 0, 1 and

The second step consists in factorizing Z in F (Z) as many times as possible until the form

where m ≥ 0 and F * (Z) = 0.

To be able to factorize Z in F (Z), we have to first evaluate F 0 (Z) at Z = 0. If F 0 (0) = 0, then no Z factorization is possible and the procedure is stopped.

take m = 1 and the following algorithmic applies:

, increase m by 1, and go to (1); 3. If F m (0) = 0, the procedure is stopped and

where,

In the most cases, F m is a linear combination of terms containing functions from sets C and S:

where a j , b j , K are constant, and α j , β j 0. By replacing (8) in (15), the following expression is obtained:

If F m (0) = 0, then the expression of Δ m (Z) is thus obtained:

After applying the procedure described for the numerator and the denominator of the coefficient, the user-friendly reformulation of the coefficient is obtained:

where N * (0) = 0 and D * (0) = 0.

In this Section we reformulate the coefficients of two relevant classes of EF numerical methods [67, 69] in terms of the C and S functions.

The first class of methods, developed by Simos et al. in [67] regards the numerical solution of the second-order Initial Value Problems (IVPs) of the form:

The scheme examined in [67] is of the form:

h 2 d 0ỹ n+1 +d 1 y n +d 2 y n−1 +d 3 y n−2 +d 2 y n−3 +d 1 y n−4 +d 0 y n−5 (20) whereỹ n+1 is determined by solving y n+1 + c 0 y n + c 1 y n−1 + c 2 y n−2 + c 1 y n−3 + c 0 y n−4 + y n−5 = h 2 c 0 y n +c 1 y n−1 +c 2 y n−2 +c 1 y n−3 +c 0 y n−4 .

and usingỹ n+1 = f (x n+1 ,ỹ n+1 ) in (20) . The classical version has the constant coefficients:

see [56] .

In [67] the exponential fitting procedure is applied to produce the following ν dependent expressions for thec coefficients :

− 2(4965191ν 2 + 13671432)ν cos(4ν) + 60500160ν cos(5ν) 

The other coefficients are untouched. They remain the same as in (23). Theoretically we must have lim

for i = 1, 2, 3, but a direct examination of the EF expressions given in (24) (25) (26) , shows that these have an indeterminate form 0/0 for ν → 0 such that, in a numerical approach, a blow up of each coefficient will be obtained when h is decreased (we remind that ν = ωh). This is removed by applying the procedure described in the previous Section. Indeed, now we have: 

The new expressions are also quotients of two ν dependent functions but, as expected, they do no longer exhibit indeterminacy when ν → 0. Also worth mentioning is that the power of Z in (13) for all these three coefficients is m = 4. The use of η functions has the same effect. In fact, by applying the procedure described in [30] , the coefficients expressed by these functions are (see also [23] ):

The full expressions can be obtained using the Mathematica modules in [30] . Both ways of deriving single formulae, instead of four, are then acceptable, and the expected theoretical behavior

is preserved.

We consider the numerical method developed by Ndukum et al. [69] to solve the first-order IVP:

The scheme used is a k-step numerical method of the form:

In the paper [69] the authors presented the coefficients of the method corresponding to k = 1, 2, 3, 4, 5. In the paper [61] the case k = 2 has been considered.

In the following we consider the case k = 3. By applying the exponential fitting procedure, the coefficients are [69] :

5ν − 11ν cos ν + 7ν cos 2ν − ν cos 3ν + 4 sin ν + 2ν 2 sin ν − 2 sin 2ν 7ν cos ν − 17ν cos 2ν + 13ν cos 3ν − 3ν cos 4ν + 4 sin ν + 2ν 2 sin 2ν − 2 sin 2ν (37) α ef 1 = −12ν + 23ν cos ν − 9ν cos 2ν − 3ν cos 3ν + ν cos 4ν − 2 sin ν − 6ν 2 sin ν − 2 sin 2ν + 2 sin 3ν 7ν cos ν − 17ν cos 2ν + 13ν cos 3ν − 3ν cos 4ν + 4 sin ν + 2ν 2 sin 2ν − 2 sin 2ν (38) α ef 2 = 7nu − 5ν cos ν − 15ν cos 2ν + 17ν cos 3ν − 4ν cos 4ν + 2 sin ν + 6ν 2 sin ν + 2 sin 2ν − 2 sin 3ν 7ν cos ν − 17ν cos 2ν + 13ν cos 3ν − 3ν cos 4ν + 4 sin ν + 2ν 2 sin 2ν − 2 sin 2ν 

which all exhibit the 0/0 indeterminacy. The coefficients modified by means of the the procedure described in the previous Section are: 

We observe that in all five coefficient the power of Z in (13) is m = 3.

As for the coefficients expressed in terms of η m (Z) functions, these are: 

Similar to the previous case, the full expression of the coefficients can be obtained using the Mathematica modules in [30] .

In this Section we show the graphs of the behavior of the coefficients and how our reformulation restores the convergence of the corresponding method.

On the left column of Fig. 1 (27) is not verified but (34) holds true. On the right column of the same figure we give additional details. On it we present the deviations of the coefficients computed by the first two approaches with respect to those expressed in terms of the η m (Z) functions. It is seen that the two reformulations (28)- (30) and (31) The same data are presented on Figs. 2 and 3 for the coefficients of Example 2. Again, the coefficients obtained in the frame of the original EF approach of [69] are oscillating and inaccurate when h → 0, in contrast with those in the Also instructive is that, in contrast with Example 1, the blow up of the original EF estimates occurs at values of h much smaller than before, and this is due to the power m in (13) which for the previous case was m = 4 while it is m = 3 by now. This allows concluding that the need of reformulations in the spirit of the approach presented in this paper is more and more stringent when m is increased.

An important issue is that of checking in what extent the accuracy in the evaluation of the coefficients affects the accuracy of the results when solving numerically a differential equation. To illustrate this we consider the following problem:

⎧ ⎪ ⎨ ⎪ ⎩ y = −100y(t) + 99 sin(t) y(0) = 1 y (0) = 11 (52) with t ∈ [0, 20π], whose analytic solution is y(t) = cos(10t) + sin(10t) + sin(t).

In Fig. 4 , we compare the method (20) 

We have shown that the reformulation of the expressions of the coefficients of EF-based numerical methods for differential equations in terms of the sets of functions C and S has two main advantages:

1. it allows reducing the original set of four expressions for each coefficient to a single one with universal use: 2. it removes completely the potential inaccuracy of the numerical solution when the step size h is small.

The procedure is then recommended as a reliable alternative to that based on η functions. Further developments of this research will be oriented to the reformulation, through C and S functions, of existing methods for ordinary differential equations [2, 17, 20, 25, 26, 28, [37] [38] [39] 41, 42, 44, 48, 51, 53, 56, 77, 80] , integral equations [5] [6] [7] [8] 10, 11, 24, 29, 32, 34, 55, 71] , stochastic problems [9, 12, 13, 18, 19, 29, 47] , fractional equations [12, 13, 21, 22, 36] , partial differential equations [57] [58] [59] 66, 74] .

Numerical solution of time fractional diffusion systems

Partitioned general linear methods for separable Hamiltonian problems

Explicit Runge-Kutta methods for initial value problems with oscillating solutions

Exponentially fitted fifthorder two step peer explicit methods

An efficient and fast parallel method for Volterra integral equations of Abel type

High performance numerical methods for Volterra equations with weakly singular kernels

Construction and implementation of two-step continuous methods for Volterra integral equations

Multistep collocation methods for Volterra integrodifferential equations

Stability issues for selected stochastic evolutionary problems: a review

Collocation methods for Volterra integral and integro-differential equations: a review

A family of multistep collocation methods for Volterra integro-differential equations

Two-step collocation methods for fractional differential equations

A spectral method for stochastic fractional differential equations

Exponentially fitted IMEX methods for advection-diffusion problems

High order exponentially fitted methods for Volterra integral equations with periodic solution

Exponential fitting direct quadrature methods for Volterra integral equations

Nearly conservative multivalue methods with extended bounded parasitism

A-stability preserving perturbation of Runge-Kutta methods for stochastic differential equations

Long-term analysis of stochastic theta-methods for damped stochastic oscillators

P-stable general Nystrom methods for y" = f(x, y)

Domain decomposition methods for a class of integropartial differential equations

Optimal Schwarz waveform relaxation for fractional diffusion-wave equations

Regularized exponentially fitted methods for oscillatory problems

Natural volterra Runge-Kutta methods

A pratical approach for the derivation of algebraically stable two-step Runge-Kutta methods

Numerical search for algebraically stable two-step almost collocation methods

Adapted explicit twostep peer methods

GPU acceleration of waveform relaxation methods for large differential systems

On the stability of ϑ-methods for stochastic Volterra integral equations

Some new uses of the ηm(Z) functions

Exponentially-fitted Gauss-Laguerre quadrature rule for integrals over an unbounded interval

A family of multistep collocation methods for Volterra integral equations

Modified Gauss-Laguerre exponential fitting based formulae

Parallel methods for weakly singular Volterra integral equations on GPUs

An exponentially fitted quadrature rule over unbounded intervals

New fractional Lanczos vector polynomials and their application to system of Abel-Volterra integral equations and fractional differential equations

Numerical integration of Hamiltonian problems by G-symplectic methods

Order conditions of general Nyström methods

General Nyström methods in Nordsieck form: error analysis

Exponentially fitted two-step Runge-Kutta methods: construction and parameter selection

Long-term stability of multi-value methods for ordinary differential equations

G-symplecticity implies conjugatesymplecticity of the underlying one-step method

Construction of the EF-based Runge-Kutta methods revisited

Search for highly stable two-step Runge-Kutta methods for ODEs

Adapted numerical methods for advection-reaction-diffusion problems generating periodic wavefronts

Parameter estimation in IMEXtrigonometrically fitted methods for the numerical solution of reaction-diffusion problems

Numerical preservation of longterm dynamics by stochastic two-step methods

Adapted numerical modelling of the Belousov-Zhabotinsky reaction

Exponentially fitted singly diagonally implicit Runge-Kutta methods

Numerical solution of a diffusion problem by exponentially fitted finite difference methods

A general framework for numerical methods solving second order differential problems

Numerical solution of reaction-diffusion systems of λ-ω type by trigonometrically fitted methods

Multivalue collocation methods free from order reduction

Revised exponentially fitted Runge-Kutta-Nyström methods

On the simultaneous approximation of a Hilbert transform and its derivatives on the real semiaxis

Explicit hybrid six-step, sixth order, fully symmetric methods for solving y = f (x, y)

The smoothed particle hydrodynamics method via residual iteration

A reaction-diffusion system of λ-ω type Part II: numerical analysis

Spiral waves for λ-ω systems

Runge-Kutta method with equation dependent coefficients

Exponential and trigonometrical fittings: user-friendly expressions for the coefficients

Exponential Fitting

A conditionally P-stable fourth-order exponentialfitting method for y = f (x, y)

A Gauss quadrature rule for oscillatory integrands

Extended quadrature rules for oscillatory integrands

Plane wave solutions to reaction-diffusion equations

Phase-fitted, six-step methods for solving x = f (t, x)

Functionally fitted explicit two step peer methods

On a family of trigonometrically fitted extended backward differentiation formulas for stiff and oscillatory initial value problems

Mean convergence of an extended Lagrange interpolation process on [0, +∞)

Nyström methods for Fredholm integral equations using equispaced points

A functional fitting Runge-Kutta method with variable coefficients

Present state-of-the-art in exponential fitting. a contribution dedicated to Liviu Ixaru on his 70-th anniversary

Absolute stability of wavetrains can explain spatiotemporal dynamics in reaction-diffusion systems of lambdaomega type

A fourth algebraic order exponentially-fitted Runge-Kutta method for the numerical solution of the Schrödinger equation

An exponentially-fitted Runge-Kutta method for the numerical integration of initial-value problems with periodic or oscillating solutions

Exponentially fitted explicit Runge-Kutta methods

Exponentially fitted quadrature rules of Gauss type for oscillatory integrands

Deferred correction with mono-implicit Runge-Kutta methods for first-order IVPs, numerical methods for differential equations

Explicit two-step peer methods