naginterfaces.library.sum.invlaplace_​crump

naginterfaces.library.sum.invlaplace_crump(f, t, relerr, alphab, tfac=0.8, mxterm=100, data=None)

invlaplace_crump estimates values of the inverse Laplace transform of a given function using a Fourier series approximation. A bound on the exponential order of the inverse is required.

For full information please refer to the NAG Library document for c06la

https://www.nag.com/numeric/nl/nagdoc_29.3/flhtml/c06/c06laf.html

Parameters

fcallable fp = f(p, data=None)

$f$ must evaluate the function $F\left(p\right)$ for a given value of $p$.

Parameters

pcomplex

The argument $p$.

dataarbitrary, optional, modifiable in place

User-communication data for callback functions.

Returns

fpcomplex

The value $F\left(p\right)$.

tfloat, array-like, shape $\left(n\right)$

Each $t[\text{j}-1]$ must specify a point at which the inverse Laplace transform is required, for $\text{j}=1,2,\dots ,n$.

relerrfloat

The required relative error in the values of the inverse Laplace transform. If the absolute value of the inverse is less than $relerr$, absolute accuracy is used instead. $relerr$ must be in the range $0.0\le relerr<1.0$. If $relerr$ is set too small or to $0.0$, the function uses a value sufficiently larger than machine precision.

alphabfloat

${\alpha }_b$, an upper bound for $\alpha$ (see Notes). Usually, ${\alpha }_b$ should be specified equal to, or slightly larger than, the value of $\alpha$. If ${\alpha }_b<\alpha$ then the prescribed accuracy may not be achieved or completely incorrect results may be obtained. If ${\alpha }_b$ is too large invlaplace_crump will be inefficient and convergence may not be achieved.

Note: it is as important to specify ${\alpha }_b$ correctly as it is to specify the correct function for inversion.

tfacfloat, optional

$t_{fac}$, a factor to be used in calculating the parameter $\tau$. Larger values (e.g., $5.0$) may be specified for difficult problems, but these may require very large values of $mxterm$.

mxtermint, optional

The maximum number of (complex) terms to be used in the evaluation of the Fourier series.

dataarbitrary, optional

User-communication data for callback functions.

Returns

valinvfloat, ndarray, shape $\left(n\right)$

An estimate of the value of the inverse Laplace transform at $t=t[\text{j}-1]$, for $\text{j}=1,2,\dots ,n$.

errestfloat, ndarray, shape $\left(n\right)$

An estimate of the error in $valinv\left[j\right]$. This is usually an estimate of relative error but, if $valinv[j-1]<relerr$, $errest\left[j\right]$ estimates the absolute error. $errest[j-1]$ is unreliable when $valinv[j-1]$ is small but slightly greater than $relerr$.

ntermsint

The number of (complex) terms actually used.

naint

The number of values of $a$ used by the function. See Further Comments.

alowfloat

The smallest value of $a$ used in the algorithm. This may be used for checking the value of $alphab-$ see Further Comments.

ahighfloat

The largest value of $a$ used in the algorithm. This may be used for checking the value of $alphab-$ see Further Comments.

nfevalint

The number of calls to $f$ made by the function.

Raises

NagValueError

(errno $1$)

On entry, $tfac=⟨value⟩$.

Constraint: $tfac>0.5$.

(errno $1$)

On entry, $relerr=⟨value⟩$.

Constraint: $0.0\le relerr<1.0$.

(errno $1$)

On entry, $mxterm=⟨value⟩$.

Constraint: $mxterm\ge 1$.

(errno $1$)

On entry, $n=⟨value⟩$.

Constraint: $n\ge 1$.

(errno $2$)

On entry, the elements of $t$ are not strictly increasing.

(errno $2$)

On entry, $t\left[0\right]=⟨value⟩$.

Constraint: $0.0\le t\left[0\right]<t\left[1\right]<\cdots <t[n-1]$.

(errno $3$)

On entry, $t[n-1]$ is too large: $t[n-1]=⟨value⟩$. If necessary, scale the problem as described in Further Comments.

(errno $4$)

Required accuracy cannot be obtained. Try increasing $tfac$, $alphab$, or both. $tfac=⟨value⟩$ and $alphab=⟨value⟩$.

(errno $5$)

Convergence failure in epsilon algorithm: $mxterm=⟨value⟩$. Some values of $valinv[j-1]$ may be calculated to the desired accuracy; this may be determined by examining the values of $errest[j-1]$. Try reducing the range of $t$ or increasing $mxterm$. If $errno$ = 5 still results, try reducing $tfac$.

Warns

NagAlgorithmicWarning

(errno $6$)

All values of $valinv$ have been calculated, but are not of the required accuracy; the values of $errest[j-1]$ should be examined carefully. Try reducing the range of $t$, or increasing $tfac$, $alphab$ or both.

Notes

No equivalent traditional C interface for this routine exists in the NAG Library.

Given a function $F\left(p\right)$ defined for complex values of $p$, invlaplace_crump estimates values of its inverse Laplace transform by Crump’s method (see Crump (1976)). (For a definition of the Laplace transform and its inverse, see the C06 Introduction.)

Crump’s method applies the epsilon algorithm (see Wynn (1956)) to the summation in Durbin’s Fourier series approximation (see Durbin (1974))

$$f\left(t_j\right)\simeq \frac{e^{a{t}_j}}{\tau }[\frac{1}{2}F\left(a\right)-\sum _{k=1}^{\infty }\{Re\left(F(a+\frac{k\pi i}{\tau })\right)\cos \left(\frac{k\pi t_j}{\tau }\right)-Im\left(F(a+\frac{k\pi i}{\tau })\right)\sin \left(\frac{k\pi t_j}{\tau }\right)\}]\text{,}$$

for $j=1,2,\dots ,n$, by choosing $a$ such that a prescribed relative error should be achieved. The method is modified slightly if $t=0.0$ so that an estimate of $f\left(0.0\right)$ can be obtained when it has a finite value. $\tau$ is calculated as $t_{fac}\times max(0.01,t_j)$, where $t_{fac}>0.5$. You specify $t_{fac}$ and ${\alpha }_b$, an upper bound on the exponential order $\alpha$ of the inverse function $f\left(t\right)$. $\alpha$ has two alternative interpretations:

1. $\alpha$ is the smallest number such that

   $$\left|f\left(t\right)\right|\le m\times exp\left(\alpha t\right)$$

   for large $t$,

2. $\alpha$ is the real part of the singularity of $F\left(p\right)$ with largest real part.

The method depends critically on the value of $\alpha$. See Further Comments for further details. The function calculates at least two different values of the argument $a$, such that $a>{\alpha }_b$, in an attempt to achieve the requested relative error and provide error estimates. The values of $t_{\text{j}}$, for $\text{j}=1,2,\dots ,n$, must be supplied in monotonically increasing order. The function calculates the values of the inverse function $f\left(t_j\right)$ in decreasing order of $j$.

References

Crump, K S, 1976, Numerical inversion of Laplace transforms using a Fourier series approximation, J. Assoc. Comput. Mach. (23), 89–96

Durbin, F, 1974, Numerical inversion of Laplace transforms: An efficient improvement to Dubner and Abate’s method, Comput. J. (17), 371–376

Wynn, P, 1956, On a device for computing the $e_m\left(S_n\right)$ transformation, Math. Tables Aids Comput. (10), 91–96