C346 Elliptic Integrals of First, Second, and Third Kind

Function subprograms RELI1, RELI2, RELI3 and DELI1, DELI2, DELI3 calculate, for real argument x, the elliptic integrals of the first, second and third kind, respectively.

On CDC and Cray computers, the double-precision versions DELI1, DELI2 and DELI3 are not available.

Mainly for reasons of numerical stability, the algorithms are based on the following definitions:

First kind:

$F$1(x,k') = ∫0x{dξ(1+ξ2)(1+k'2ξ2)}(k'2≥0).

Second kind:

$F$2(x,k',a,b) = ∫0x{a+bξ2(1+ξ2)(1+ξ2)(1+k'2ξ2)} dξ(k'2≥0).

Third kind:

$F$3(x,k',p) = ∫0x{1+ξ2(1+pξ2)(1+ξ2)(1+k'2ξ2)} dξ(k'2≥0, px2≠-1).

Note that $F$₁(x,k') = F₂(x,k',1,1) =F₃(x,k',1) . For p < 0, the integral $F$₃ is defined by its principal value.

For the integral of the second kind, a special entry-mode argument is provided which allows $F$₂(x,k',a,b) to be calculated when , i.e. when $k\text{'}$ is imaginary.

Other common definitions of the elliptic integrals and their relations to $F$₁ , $F$₂ , $F$₃ are listed here for convenience ($k2+k\text{'}2=1$ ):

First kind:

$F(k,\phi )\; =\; \int$0φ{dψ1-k2sin2ψ} = F1(tanφ,k')(|φ| ≤π/2, |k| < 1),

$F(y,k)\; =\; \int$0y{dη(1-η2)(1-k2η2)} = F1( y/1-y2,k' )(|y| < 1, |k| < 1).

Second kind:

$E(k,\phi )\; =\; \int$0φ1-k2sin2ψ dψ = F2(tanφ,k',1,k'2)(|φ| ≤π/2, |k| ≤1),

$E(y,k)\; =\; \int$0y{1-k2η21-η2} dη = F2( y/1-y2,k',1,k'2)(|y| < 1, |k| ≤1).

Third kind:

$\Pi (\phi ,h,k)$	$=$	$\int$0φ{dψ(1+hsin2ψ)1-k2sin2ψ}	$=$	$F$3(tanφ,k',h+1)
$\Pi (y,h,k)$	$=$	$\int$0y{dη(1+hη2)(1-η2)(1-k2η2)}	$=$	$F$3( y/1-y2,k',h+1 )

Structure:

FUNCTION subprograms
User Entry Names: RELI1, RELI2, RELI3, DELI1, DELI2, DELI3
Files Referenced: Unit 6
External References: ASINH, DASINH (B102), MTLMTR (N002), ABEND (Z035)

Usage:

In any arithmetic expression, with $AKP=k\text{'}$ ,

RELI1(X,AKP)	DELI1(X,AKP)	$F$1(X,k') ,
RELI2(X,AKP,A,B,MODE)	DELI2(X,AKP,A,B,MODE)	$F$2(X,k'',A,B) ,
RELI3(X,AKP,P)	DELI3(X,AKP,P)	$F$3(X,k',P) ,

where RELI1, RELI2, RELI3 are of type REAL, where DELI1, DELI2, DELI3 are of type DOUBLE PRECISION, and X, AKP, A, B and P have the same type as the function name. MODE is of type INTEGER.

The notation $k\text{'}\text{'}$ indicates that, when calling RELI2 or DRELI2, the parameters AKP and MODE must be set as follows:
If $k\text{'}2>\; 0$ : $MODE\; =\; +1$ and $AKP=\; k\text{'}$ ,
if : $MODE\; =\; -1$ and $AKP=$ Im $k\text{'}\; =\; -ik\text{'}$

(real).

Method:

The evaluation of $F$₁ and $F$₂ is based on the Landen transformation, that of $F$₃ on the Bartky transformation. $F$₂ for is calculated by using a transformation of the arguments. See Ref. 1 and 2 for details.

Accuracy:

The REAL functions (except on CDC and Cray computers) have full single-precision accuracy. The REAL functions on CDC and Cray, and the DOUBLE PRECISION functions on all computers have an accuracy approximately two significant digits less than the machine precision.

Restrictions:

1. RELI2 and DELI2: \ AKP*X**2 < 1 if MODE =-1.
2. RELI2 and DELI2: \ $MODE\; =\; +1$ or $-1$ .
3. RELI3 and DELI3: P*X**2 $\ne -$ 1.

Error handling:

Error C346.1: Restriction 1 is not satisfied.
Error C346.2: Restriction 2 is not satisfied.
Error C346.3: Restriction 3 is not satisfied.
In all cases, the function value is set equal to zero, and a message is written on Unit 6, unless subroutine MTLSET (N002) has been called.

The subprograms are based on the Algol60 procedures el1, el2 in Ref. 1 and el3 in Ref. 2.

References:

1. R. Bulirsch, Numerical calculation of elliptic integrals and elliptic functions, Numer. Math. 7 (1965) 78--90.
2. R. Bulirsch, Numerical calculation of elliptic integrals and elliptic functions III, Numer. Math. 13 (1969) 305--315.

C347