Package netlib.slatec

Class Drd

java.lang.Object
netlib.slatec.Drd

public class Drd
extends java.lang.Object
F2J-Version of Slatec, DRD.F, Revision 920501.

 Changes since Revision 920501:

 - 2012-11-21, 2013-04-30, Reiser:
    - Made LOLIM, C1,...,C4 static constants. Note that UPLIM is not a
      constant but depends on ERRTOL.
    - Removed initialization on first call (if (first) ...).
 - 2012-12-16, Reiser:
    - corrected unintended casting to float (due to mixed mode arithmetic bug
      in F2J and the usage of single notation in original DRD) that yielded
      tuplim = infinity:
        tuplim =
           (double)(((float)Math.pow(D1mach.d1mach(1), ((1.0E0f / 3.0E0f)))));
 - 2013-04-30, Reiser:
    - Introduced runtime-settable ERRTOL.
    - Default compile-time ERRTOL settings: Set AUTO_ERRTOL=false for tighter
      default ERRTOL. Set AUTO_ERRTOL=true for original SLATEC error bounds.

 -----------------------------------------------------------------------------

  Produced by f2java.  f2java is part of the Fortran-
  -to-Java project at the University of Tennessee Netlib
  numerical software repository.

  Original authorship for the BLAS and LAPACK numerical
  routines may be found in the Fortran source, available at
  http://www.netlib.org.

  Fortran input file: drd.f
  f2java version: 0.8.1

DECK DRD
 C***BEGIN PROLOGUE  DRD
 C***PURPOSE  Compute the incomplete or complete elliptic integral of
 C            the 2nd kind. For X and Y nonnegative, X+Y and Z positive,
 C            DRD(X,Y,Z) = Integral from zero to infinity of
 C                                -1/2     -1/2     -3/2
 C                      (3/2)(t+X)    (t+Y)    (t+Z)    dt.
 C            If X or Y is zero, the integral is complete.
 C***LIBRARY   SLATEC
 C***CATEGORY  C14
 C***TYPE      DOUBLE PRECISION (RD-S, DRD-D)
 C***KEYWORDS  COMPLETE ELLIPTIC INTEGRAL, DUPLICATION THEOREM,
 C             INCOMPLETE ELLIPTIC INTEGRAL, INTEGRAL OF THE SECOND KIND,
 C             TAYLOR SERIES
 C***AUTHOR  Carlson, B. C.
 C             Ames Laboratory-DOE
 C             Iowa State University
 C             Ames, IA  50011
 C           Notis, E. M.
 C             Ames Laboratory-DOE
 C             Iowa State University
 C             Ames, IA  50011
 C           Pexton, R. L.
 C             Lawrence Livermore National Laboratory
 C             Livermore, CA  94550
 C***DESCRIPTION
 C
 C   1.     DRD
 C          Evaluate an INCOMPLETE (or COMPLETE) ELLIPTIC INTEGRAL
 C          of the second kind
 C          Standard FORTRAN function routine
 C          Double precision version
 C          The routine calculates an approximation result to
 C          DRD(X,Y,Z) = Integral from zero to infinity of
 C                              -1/2     -1/2     -3/2
 C                    (3/2)(t+X)    (t+Y)    (t+Z)    dt,
 C          where X and Y are nonnegative, X + Y is positive, and Z is
 C          positive.  If X or Y is zero, the integral is COMPLETE.
 C          The duplication theorem is iterated until the variables are
 C          nearly equal, and the function is then expanded in Taylor
 C          series to fifth order.
 C
 C   2.     Calling Sequence
 C
 C          DRD( X, Y, Z, IER )
 C
 C          Parameters On Entry
 C          Values assigned by the calling routine
 C
 C          X      - Double precision, nonnegative variable
 C
 C          Y      - Double precision, nonnegative variable
 C
 C                   X + Y is positive
 C
 C          Z      - Double precision, positive variable
 C
 C
 C
 C          On Return    (values assigned by the DRD routine)
 C
 C          DRD     - Double precision approximation to the integral
 C
 C
 C          IER    - Integer
 C
 C                   IER = 0 Normal and reliable termination of the
 C                           routine. It is assumed that the requested
 C                           accuracy has been achieved.
 C
 C                   IER >  0 Abnormal termination of the routine
 C
 C
 C          X, Y, Z are unaltered.
 C
 C   3.    Error Messages
 C
 C         Value of IER assigned by the DRD routine
 C
 C                  Value assigned         Error message printed
 C                  IER = 1                MIN(X,Y) .LT. 0.0D0
 C                      = 2                MIN(X + Y, Z ) .LT. LOLIM
 C                      = 3                MAX(X,Y,Z) .GT. UPLIM
 C
 C
 C   4.     Control Parameters
 C
 C                  Values of LOLIM, UPLIM, and ERRTOL are set by the
 C                  routine.
 C
 C          LOLIM and UPLIM determine the valid range of X, Y, and Z
 C
 C          LOLIM  - Lower limit of valid arguments
 C
 C                    Not less  than 2 / (machine maximum) ** (2/3).
 C
 C          UPLIM  - Upper limit of valid arguments
 C
 C                 Not greater than (0.1D0 * ERRTOL / machine
 C                 minimum) ** (2/3), where ERRTOL is described below.
 C                 In the following table it is assumed that ERRTOL will
 C                 never be chosen smaller than 1.0D-5.
 C
 C
 C                    Acceptable values for:   LOLIM      UPLIM
 C                    IBM 360/370 SERIES   :   6.0D-51     1.0D+48
 C                    CDC 6000/7000 SERIES :   5.0D-215    2.0D+191
 C                    UNIVAC 1100 SERIES   :   1.0D-205    2.0D+201
 C                    CRAY                 :   3.0D-1644   1.69D+1640
 C                    VAX 11 SERIES        :   1.0D-25     4.5D+21
 C
 C
 C          ERRTOL determines the accuracy of the answer
 C
 C                 The value assigned by the routine will result
 C                 in solution precision within 1-2 decimals of
 C                 "machine precision".
 C
 C          ERRTOL    Relative error due to truncation is less than
 C                    3 * ERRTOL ** 6 / (1-ERRTOL) ** 3/2.
 C
 C
 C
 C        The accuracy of the computed approximation to the integral
 C        can be controlled by choosing the value of ERRTOL.
 C        Truncation of a Taylor series after terms of fifth order
 C        introduces an error less than the amount shown in the
 C        second column of the following table for each value of
 C        ERRTOL in the first column.  In addition to the truncation
 C        error there will be round-off error, but in practice the
 C        total error from both sources is usually less than the
 C        amount given in the table.
 C
 C
 C
 C
 C          Sample choices:  ERRTOL   Relative truncation
 C                                    error less than
 C                           1.0D-3    4.0D-18
 C                           3.0D-3    3.0D-15
 C                           1.0D-2    4.0D-12
 C                           3.0D-2    3.0D-9
 C                           1.0D-1    4.0D-6
 C
 C
 C                    Decreasing ERRTOL by a factor of 10 yields six more
 C                    decimal digits of accuracy at the expense of one or
 C                    two more iterations of the duplication theorem.
 C
 C *Long Description:
 C
 C   DRD Special Comments
 C
 C
 C
 C          Check: DRD(X,Y,Z) + DRD(Y,Z,X) + DRD(Z,X,Y)
 C          = 3 / SQRT(X * Y * Z), where X, Y, and Z are positive.
 C
 C
 C          On Input:
 C
 C          X, Y, and Z are the variables in the integral DRD(X,Y,Z).
 C
 C
 C          On Output:
 C
 C
 C          X, Y, Z are unaltered.
 C
 C
 C
 C          ********************************************************
 C
 C          WARNING: Changes in the program may improve speed at the
 C                   expense of robustness.
 C
 C
 C
 C    -------------------------------------------------------------------
 C
 C
 C   Special double precision functions via DRD and DRF
 C
 C
 C                  Legendre form of ELLIPTIC INTEGRAL of 2nd kind
 C
 C                  -----------------------------------------
 C
 C
 C                                             2         2   2
 C                  E(PHI,K) = SIN(PHI) DRF(COS (PHI),1-K SIN (PHI),1) -
 C
 C                     2      3             2         2   2
 C                  -(K/3) SIN (PHI) DRD(COS (PHI),1-K SIN (PHI),1)
 C
 C
 C                                  2        2            2
 C                  E(K) = DRF(0,1-K ,1) - (K/3) DRD(0,1-K ,1)
 C
 C                         PI/2     2   2      1/2
 C                       = INT  (1-K SIN (PHI) )  D PHI
 C                          0
 C
 C                  Bulirsch form of ELLIPTIC INTEGRAL of 2nd kind
 C
 C                  -----------------------------------------
 C
 C                                               2 2    2
 C                  EL2(X,KC,A,B) = AX DRF(1,1+KC X ,1+X ) +
 C
 C                                              3          2 2    2
 C                                 +(1/3)(B-A) X DRD(1,1+KC X ,1+X )
 C
 C
 C
 C
 C                  Legendre form of alternative ELLIPTIC INTEGRAL
 C                  of 2nd kind
 C
 C                  -----------------------------------------
 C
 C
 C
 C                            Q     2       2   2  -1/2
 C                  D(Q,K) = INT SIN P  (1-K SIN P)     DP
 C                            0
 C
 C
 C
 C                                     3          2     2   2
 C                  D(Q,K) = (1/3) (SIN Q) DRD(COS Q,1-K SIN Q,1)
 C
 C
 C
 C
 C                  Lemniscate constant  B
 C
 C                  -----------------------------------------
 C
 C
 C
 C
 C                       1    2    4 -1/2
 C                  B = INT  S (1-S )    DS
 C                       0
 C
 C
 C                  B = (1/3) DRD (0,2,1)
 C
 C
 C                  Heuman's LAMBDA function
 C
 C                  -----------------------------------------
 C
 C
 C
 C                  (PI/2) LAMBDA0(A,B) =
 C
 C                                    2                2
 C                 = SIN(B) (DRF(0,COS (A),1)-(1/3) SIN (A) *
 C
 C                            2               2         2       2
 C                  *DRD(0,COS (A),1)) DRF(COS (B),1-COS (A) SIN (B),1)
 C
 C                            2       3             2
 C                  -(1/3) COS (A) SIN (B) DRF(0,COS (A),1) *
 C
 C                           2         2       2
 C                   *DRD(COS (B),1-COS (A) SIN (B),1)
 C
 C
 C
 C                  Jacobi ZETA function
 C
 C                  -----------------------------------------
 C
 C                             2                 2       2   2
 C                  Z(B,K) = (K/3) SIN(B) DRF(COS (B),1-K SIN (B),1)
 C
 C
 C                                       2             2
 C                             *DRD(0,1-K ,1)/DRF(0,1-K ,1)
 C
 C                               2       3           2       2   2
 C                            -(K /3) SIN (B) DRD(COS (B),1-K SIN (B),1)
 C
 C
 C ---------------------------------------------------------------------
 C
 C***REFERENCES  B. C. Carlson and E. M. Notis, Algorithms for incomplete
 C                 elliptic integrals, ACM Transactions on Mathematical
 C                 Software 7, 3 (September 1981), pp. 398-403.
 C               B. C. Carlson, Computing elliptic integrals by
 C                 duplication, Numerische Mathematik 33, (1979),
 C                 pp. 1-16.
 C               B. C. Carlson, Elliptic integrals of the first kind,
 C                 SIAM Journal of Mathematical Analysis 8, (1977),
 C                 pp. 231-242.
 C***ROUTINES CALLED  D1MACH, XERMSG
 C***REVISION HISTORY  (YYMMDD)
 C   790801  DATE WRITTEN
 C   890531  Changed all specific intrinsics to generic.  (WRB)
 C   890531  REVISION DATE from Version 3.2
 C   891214  Prologue converted to Version 4.0 format.  (BAB)
 C   900315  CALLs to XERROR changed to CALLs to XERMSG.  (THJ)
 C   900326  Removed duplicate information from DESCRIPTION section.
 C           (WRB)
 C   900510  Modify calls to XERMSG to put in standard form.  (RWC)
 C   920501  Reformatted the REFERENCES section.  (WRB)
 C***END PROLOGUE  DRD

 
Author:
Stefan Reiser (s.reiser@tu-bs.de)
  • Field Summary

    Fields 
    Modifier and Type Field Description
    static double LOLIM  
  • Method Summary

    Modifier and Type Method Description
    static double drd​(double x, double y, double z, intW ier)  
    static double getErrtol()  
    static double getUplim()  
    static void setErrtol​(double errtol_)
    Sets ERRTOL (and as a side effect UPLIM).

    Methods inherited from class java.lang.Object

    clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait
  • Field Details

    • LOLIM

      public static final double LOLIM
  • Method Details