NAG FL Interface
s14adf (polygamma_​deriv)

Settings help

FL Name Style:


FL Specification Language:


1 Purpose

s14adf returns a sequence of values of scaled derivatives of the psi function ψ(x) (also known as the digamma function).

2 Specification

Fortran Interface
Subroutine s14adf ( x, n, m, ans, ifail)
Integer, Intent (In) :: n, m
Integer, Intent (Inout) :: ifail
Real (Kind=nag_wp), Intent (In) :: x
Real (Kind=nag_wp), Intent (Out) :: ans(m)
C Header Interface
#include <nag.h>
void  s14adf_ (const double *x, const Integer *n, const Integer *m, double ans[], Integer *ifail)
The routine may be called by the names s14adf or nagf_specfun_polygamma_deriv.

3 Description

s14adf computes m values of the function
w(k,x)=(-1)k+1ψ (k) (x) k! ,  
for x>0, k=n, n+1,,n+m-1, where ψ is the psi function
ψ(x)=ddx lnΓ(x)=Γ(x) Γ(x) ,  
and ψ (k) denotes the kth derivative of ψ.
The routine is derived from the routine PSIFN in Amos (1983). The basic method of evaluation of w(k,x) is the asymptotic series
w(k,x)ε(k,x)+12xk+1 +1xkj=1B2j(2j+k-1)! (2j)!k!x2j  
for large x greater than a machine-dependent value xmin, followed by backward recurrence using
w(k,x)=w(k,x+1)+x-k-1  
for smaller values of x, where ε(k,x)=-lnx when k=0, ε(k,x)= 1kxk when k>0, and B2j, j=1,2,, are the Bernoulli numbers.
When k is large, the above procedure may be inefficient, and the expansion
w(k,x)=j=11(x+j)k+1,  
which converges rapidly for large k, is used instead.

4 References

NIST Digital Library of Mathematical Functions
Amos D E (1983) Algorithm 610: A portable FORTRAN subroutine for derivatives of the psi function ACM Trans. Math. Software 9 494–502

5 Arguments

1: x Real (Kind=nag_wp) Input
On entry: the argument x of the function.
Constraint: x>0.0.
2: n Integer Input
On entry: the index of the first member n of the sequence of functions.
Constraint: n0.
3: m Integer Input
On entry: the number of members m required in the sequence w(k,x), for k=n,,n+m-1.
Constraint: m1.
4: ans(m) Real (Kind=nag_wp) array Output
On exit: the first m elements of ans contain the required values w(k,x), for k=n,,n+m-1.
5: ifail Integer Input/Output
On entry: ifail must be set to 0, −1 or 1 to set behaviour on detection of an error; these values have no effect when no error is detected.
A value of 0 causes the printing of an error message and program execution will be halted; otherwise program execution continues. A value of −1 means that an error message is printed while a value of 1 means that it is not.
If halting is not appropriate, the value −1 or 1 is recommended. If message printing is undesirable, then the value 1 is recommended. Otherwise, the value 0 is recommended. When the value -1 or 1 is used it is essential to test the value of ifail on exit.
On exit: ifail=0 unless the routine detects an error or a warning has been flagged (see Section 6).

6 Error Indicators and Warnings

If on entry ifail=0 or −1, explanatory error messages are output on the current error message unit (as defined by x04aaf).
Errors or warnings detected by the routine:
ifail=1
On entry, x=value.
Constraint: x>0.0.
ifail=2
On entry, n=value.
Constraint: n0.
ifail=3
On entry, m=value.
Constraint: m1.
ifail=4
Computation abandoned due to the likelihood of underflow.
ifail=5
Computation abandoned due to the likelihood of overflow.
ifail=6
There is not enough internal workspace to continue computation. m is probably too large.
ifail=-99
An unexpected error has been triggered by this routine. Please contact NAG.
See Section 7 in the Introduction to the NAG Library FL Interface for further information.
ifail=-399
Your licence key may have expired or may not have been installed correctly.
See Section 8 in the Introduction to the NAG Library FL Interface for further information.
ifail=-999
Dynamic memory allocation failed.
See Section 9 in the Introduction to the NAG Library FL Interface for further information.

7 Accuracy

All constants in s14adf are given to approximately 18 digits of precision. Calling the number of digits of precision in the floating-point arithmetic being used t, then clearly the maximum number of correct digits in the results obtained is limited by p=min(t,18). Empirical tests of s14adf, taking values of x in the range 0.0<x<50.0, and n in the range 1n50, have shown that the maximum relative error is a loss of approximately two decimal places of precision. Tests with n=0, i.e., testing the function -ψ(x), have shown somewhat better accuracy, except at points close to the zero of ψ(x), x1.461632, where only absolute accuracy can be obtained.

8 Parallelism and Performance

Background information to multithreading can be found in the Multithreading documentation.
s14adf is not threaded in any implementation.

9 Further Comments

The time taken for a call of s14adf is approximately proportional to m, plus a constant. In general, it is much cheaper to call s14adf with m greater than 1 to evaluate the function w(k,x), for k=n,,n+m-1, rather than to make m separate calls of s14adf.

10 Example

This example reads values of the argument x from a file, evaluates the function at each value of x and prints the results.

10.1 Program Text

Program Text (s14adfe.f90)

10.2 Program Data

Program Data (s14adfe.d)

10.3 Program Results

Program Results (s14adfe.r)