gonum.org/v1/gonum@v0.14.0/mathext/internal/amos/amoslib/zs1s2.f (about)

     1        SUBROUTINE ZS1S2(ZRR, ZRI, S1R, S1I, S2R, S2I, NZ, ASCLE, ALIM,
     2       * IUF)
     3  C***BEGIN PROLOGUE  ZS1S2
     4  C***REFER TO  ZBESK,ZAIRY
     5  C
     6  C     ZS1S2 TESTS FOR A POSSIBLE UNDERFLOW RESULTING FROM THE
     7  C     ADDITION OF THE I AND K FUNCTIONS IN THE ANALYTIC CON-
     8  C     TINUATION FORMULA WHERE S1=K FUNCTION AND S2=I FUNCTION.
     9  C     ON KODE=1 THE I AND K FUNCTIONS ARE DIFFERENT ORDERS OF
    10  C     MAGNITUDE, BUT FOR KODE=2 THEY CAN BE OF THE SAME ORDER
    11  C     OF MAGNITUDE AND THE MAXIMUM MUST BE AT LEAST ONE
    12  C     PRECISION ABOVE THE UNDERFLOW LIMIT.
    13  C
    14  C***ROUTINES CALLED  ZABS,ZEXP,ZLOG
    15  C***END PROLOGUE  ZS1S2
    16  C     COMPLEX CZERO,C1,S1,S1D,S2,ZR
    17        DOUBLE PRECISION AA, ALIM, ALN, ASCLE, AS1, AS2, C1I, C1R, S1DI,
    18       * S1DR, S1I, S1R, S2I, S2R, ZEROI, ZEROR, ZRI, ZRR, ZABS
    19        INTEGER IUF, IDUM, NZ
    20        DATA ZEROR,ZEROI  / 0.0D0 , 0.0D0 /
    21        NZ = 0
    22        AS1 = ZABS(CMPLX(S1R,S1I,kind=KIND(1.0D0)))
    23        AS2 = ZABS(CMPLX(S2R,S2I,kind=KIND(1.0D0)))
    24        IF (S1R.EQ.0.0D0 .AND. S1I.EQ.0.0D0) GO TO 10
    25        IF (AS1.EQ.0.0D0) GO TO 10
    26        ALN = -ZRR - ZRR + DLOG(AS1)
    27        S1DR = S1R
    28        S1DI = S1I
    29        S1R = ZEROR
    30        S1I = ZEROI
    31        AS1 = ZEROR
    32        IF (ALN.LT.(-ALIM)) GO TO 10
    33        CALL ZLOG(S1DR, S1DI, C1R, C1I, IDUM)
    34        C1R = C1R - ZRR - ZRR
    35        C1I = C1I - ZRI - ZRI
    36        CALL ZEXP(C1R, C1I, S1R, S1I)
    37        AS1 = ZABS(CMPLX(S1R,S1I,kind=KIND(1.0D0)))
    38        IUF = IUF + 1
    39     10 CONTINUE
    40        AA = DMAX1(AS1,AS2)
    41        IF (AA.GT.ASCLE) THEN
    42          RETURN
    43        END IF
    44        S1R = ZEROR
    45        S1I = ZEROI
    46        S2R = ZEROR
    47        S2I = ZEROI
    48        NZ = 1
    49        IUF = 0
    50        RETURN
    51        END