amos

Documentation

Overview

Package amos implements functions from Algorithm 644, Transactions on Mathematical Software

Index

• func ABS(a int) int
• func DGAMLN(z float64) float64
• func MAX(a int, b int) int
• func MIN(a int, b int) int
• func ZABS(ZR float64, ZI float64) float64
• func ZACAI(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, ...) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, ...)
• func ZACON(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, ...) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, ...)
• func ZAIRY(ZR float64, ZI float64, ID int, KODE int) (float64, float64, int, int)
• func ZASYI(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, float64, ...)
• func ZBESH(ZR float64, ZI float64, FNU float64, KODE int, M int, N int, CYR []float64, ...) (float64, float64, float64, int, int, int, []float64, []float64, int, int)
• func ZBESI(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, int)
• func ZBESJ(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, int)
• func ZBESK(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, int)
• func ZBESY(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, []float64, ...)
• func ZBINU(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, float64, ...)
• func ZBIRY(ZR float64, ZI float64, ID int, KODE int) (float64, float64, int)
• func ZBKNU(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, float64, ...)
• func ZBUNI(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, int, int, ...)
• func ZBUNK(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, ...) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, ...)
• func ZDIV(AR float64, AI float64, BR float64, BI float64) (float64, float64)
• func ZEXP(AR float64, AI float64) (float64, float64)
• func ZFLOATS(Z complex128) (float64, float64)
• func ZKSCL(ZRR float64, ZRI float64, FNU float64, N int, YR []float64, YI []float64, ...) (float64, float64, float64, int, []float64, []float64, int, float64, float64, ...)
• func ZLOG(AR float64, AI float64) (float64, float64)
• func ZMLRI(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, float64)
• func ZMLT(AR float64, AI float64, BR float64, BI float64) (float64, float64)
• func ZRATI(ZR float64, ZI float64, FNU float64, N int, CYR []float64, CYI []float64, ...) (float64, float64, float64, int, []float64, []float64, float64)
• func ZS1S2(ZRR float64, ZRI float64, S1R float64, S1I float64, S2R float64, S2I float64, ...) (float64, float64, float64, float64, float64, float64, int, float64, float64, ...)
• func ZSERI(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, float64, ...)
• func ZSHCH(ZR float64, ZI float64) (float64, float64, float64, float64)
• func ZSQRT(AR float64, AI float64) (float64, float64)
• func ZUCHK(YR float64, YI float64, NZ int, ASCLE float64, TOL float64) (float64, float64, int, float64, float64)
• func ZUNHJ(ZR float64, ZI float64, FNU float64, IPMTR int, TOL float64, PHIR float64, ...) (float64, float64, float64, int, float64, float64, float64, float64, float64, ...)
• func ZUNI1(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, int, float64, ...)
• func ZUNI2(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, int, float64, ...)
• func ZUNIK(ZRR float64, ZRI float64, FNU float64, IKFLG int, IPMTR int, TOL float64, ...) (float64, float64, float64, int, int, float64, int, float64, float64, float64, ...)
• func ZUNK1(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, ...) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, ...)
• func ZUNK2(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, ...) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, ...)
• func ZUOIK(ZR float64, ZI float64, FNU float64, KODE int, IKFLG int, N int, YR []float64, ...) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, ...)
• func ZWRSK(ZRR float64, ZRI float64, FNU float64, KODE int, N int, YR []float64, ...) (float64, float64, float64, int, int, []float64, []float64, int, []float64, ...)

Functions

func ABS

func ABS(a int) int

ABS returns the absolute value of a

func DGAMLN

func DGAMLN(z float64) float64

DGAMLN computes the natural log of the gamma function

func MAX

func MAX(a int, b int) int

MAX returns the maximum value of a and b

func MIN

func MIN(a int, b int) int

MIN returns the minimum value of a and b

func ZABS

func ZABS(ZR float64, ZI float64) float64

ZABS computes the absolute value or magnitude of a complex number real and imag parts

func ZACAI

func ZACAI(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, YI []float64, NZ int, RL float64, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, float64, float64, float64)

ZACAI applies the analytic continuation formula

K(FNU,ZN*EXP(MP))=K(FNU,ZN)*EXP(-MP*FNU) - MP*I(FNU,ZN) MP=PI*MR*CMPLX(0.0,1.0)

to continue the K function from the right half to the left half z plane for use with ZAIRY where FNU=1/3 or 2/3 and N=1. ZACAI is the same as ZACON with the parts for larger orders an recurrence removed. A recursive call to ZACON can result if ZACON is called from ZAIRY.

func ZACON

func ZACON(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, YI []float64, NZ int, RL float64,
	FNUL float64, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, float64, float64, float64, float64)

ZACON APPLIES THE ANALYTIC CONTINUATION FORMULA

K(FNU,ZN*EXP(MP))=K(FNU,ZN)*EXP(-MP*FNU) - MP*I(FNU,ZN)

MP=PI*MR*CMPLX(0.0,1.0)

TO CONTINUE THE K FUNCTION FROM THE RIGHT HALF TO THE LEFT HALF Z PLANE

func ZAIRY

func ZAIRY(ZR float64, ZI float64, ID int, KODE int) (float64, float64, int, int)

ZAIRY computes the Airy functions A(z) and DAI(z) for complex z For KODE=1, ZAIRY computes the complex airy function Ai(z) or its derivative dAi(z)/dz on ID=0 or ID=1 respectively. For KODE=2, a scaling option Exp(ζ)*Ai(z) or Exp(ζ)*dAi(z)/dz is provided to remove the exponential decay in -𝛑/3 < arg(z) < 𝛑/3 and the exponential growth in 𝛑/3 < |arg(z)| < 𝛑 where ζ=(2/3)*z^(3/2).

While the airy functions Ai(z) and dAi(z)/dz are analytic in the whole z plane, the corresponding scaled functions defined for KODE=2 have a cut along the negative real axis. Defintions and notation are found in the NBS Handbook of Mathematical Functions.

Input

ZR,ZI  - Z=CMPLX(ZR,ZI)
ID     - order of derivative, ID=0 OR ID=1
KODE   - a parameter to indicate the scaling option
         KODE= 1 returns AI=AI(Z) on ID=0 or AI=DAI(Z)/DZ on ID=1
         KODE= 2 returns AI=EXP(ζ)*AI(Z) on ID=0 or AI=EXP(ζ)*DAI(Z)/DZ on ID=1
                 where ζ=(2/3)*Z^(3/2)

Ouput

AIR,AII- complex answer depending on the choices for ID and KODE
NZ     - underflow indicator
         	NZ= 0   , normal return
         	NZ= 1   , AI=0 due to underflow in -𝛑3 < ARG(Z) < 𝛑/3 on KODE=1
IERR   - error flag
         IERR=0, normal return - computation completed
         IERR=1, input error   - no computation
         IERR=2, overflow      - no computation, real(ζ) too large on KODE=1
         IERR=3, |Z| large     - computation completed losses of signifcance by argument reduction
                                 produce less than half of machine accuracy
         IERR=4, |Z| too large - no computation complete loss of accuracy by argument reduction
         IERR=5, error         - no computation, algorithm termination condition not met

AI and DAI are computed for |Z| > 1.0 from the K Bessel functions by

AI(Z)=C*SQRT(Z)*K(1/3,ζ)
DAI(Z)=-C*Z*K(2/3,ζ)
C=1.0/(𝛑*SQRT(3.0))
ζ=(2/3)*Z^(3/2)

with the power series for |z| <= 1.0

func ZASYI

func ZASYI(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, YI []float64, NZ int,
	RL float64, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, []float64, []float64, int,
	float64, float64, float64, float64)

ZASYI COMPUTES THE I BESSEL FUNCTION FOR REAL(Z) >= 0.0 BY MEANS OF THE ASYMPTOTIC EXPANSION FOR LARGE CABS(Z) IN THE REGION CABS(Z) > MAX(RL,FNU*FNU/2). NZ=0 IS A NORMAL RETURN. NZ < 0 INDICATES AN OVERFLOW ON KODE=1.

func ZBESH

func ZBESH(ZR float64, ZI float64, FNU float64, KODE int, M int, N int, CYR []float64, CYI []float64, NZ int, IERR int) (float64, float64, float64, int, int, int, []float64, []float64, int, int)

ZBESH COMPUTE THE HANKEL FUNCTION OF A COMPLEX ARGUMENT

FOR KODE=1 COMPUTES AN N MEMBER SEQUENCE OF COMPLEX HANKEL (BESSEL) FUNCTIONS CY(J)=H(M,FNU+J-1,Z) FOR KINDS M=1 OR 2, REAL, NONNEGATIVE ORDERS FNU+J-1, J=1,...,N, AND COMPLEX Z != CMPLX(0.0,0.0) IN THE CUT PLANE -PI < ARG(Z) <= PI. ON KODE=2, ZBESH RETURNS THE SCALED HANKEL FUNCTIONS

CY(I)=EXP(-MM*Z*I)*H(M,FNU+J-1,Z) MM=3-2*M, I**2=-1.

WHICH REMOVES THE EXPONENTIAL BEHAVIOR IN BOTH THE UPPER AND LOWER HALF PLANES. DEFINITIONS AND NOTATION ARE FOUND IN THE NBS HANDBOOK OF MATHEMATICAL FUNCTIONS (REF. 1).

INPUT ZR,ZI,FNU ARE DOUBLE PRECISION

ZR,ZI  - Z=CMPLX(ZR,ZI), Z != CMPLX(0.0e0,0.0e0),
         -PT < ARG(Z) <=  PI
FNU    - ORDER OF INITIAL H FUNCTION, FNU >= 0.0e0
KODE   - A PARAMETER TO INDICATE THE SCALING OPTION
         KODE= 1  RETURNS
                  CY(J)=H(M,FNU+J-1,Z),   J=1,...,N
             = 2  RETURNS
                  CY(J)=H(M,FNU+J-1,Z)*EXP(-I*Z*(3-2M))
                       J=1,...,N  ,  I**2=-1
M      - KIND OF HANKEL FUNCTION, M=1 OR 2
N      - NUMBER OF MEMBERS IN THE SEQUENCE, N >= 1

OUTPUT CYR,CYI ARE DOUBLE PRECISION

CYR,CYI- DOUBLE PRECISION VECTORS WHOSE FIRST N COMPONENTS
         CONTAIN REAL AND IMAGINARY PARTS FOR THE SEQUENCE
         CY(J)=H(M,FNU+J-1,Z)  OR
         CY(J)=H(M,FNU+J-1,Z)*EXP(-I*Z*(3-2M))  J=1,...,N
         DEPENDING ON KODE, I**2=-1.
NZ     - NUMBER OF COMPONENTS SET TO ZERO DUE TO UNDERFLOW,
         NZ= 0   , NORMAL RETURN
         NZ > 0 , FIRST NZ COMPONENTS OF CY SET TO ZERO DUE
                   TO UNDERFLOW, CY(J)=CMPLX(0.0e0,0.0e0)
                   J=1,...,NZ WHEN Y > 0.0 AND M=1 OR
                   Y < 0.0 AND M=2. FOR THE COMPLMENTARY
                   HALF PLANES, NZ STATES ONLY THE NUMBER
                   OF UNDERFLOWS.
IERR   - ERROR FLAG
         IERR=0, NORMAL RETURN - COMPUTATION COMPLETED
         IERR=1, INPUT ERROR   - NO COMPUTATION
         IERR=2, OVERFLOW      - NO COMPUTATION, FNU TOO
                 LARGE OR CABS(Z) TOO SMALL OR BOTH
         IERR=3, CABS(Z) OR FNU+N-1 LARGE - COMPUTATION DONE
                 BUT LOSSES OF SIGNIFCANCE BY ARGUMENT
                 REDUCTION PRODUCE LESS THAN HALF OF MACHINE
                 ACCURACY
         IERR=4, CABS(Z) OR FNU+N-1 TOO LARGE - NO COMPUTA-
                 TION BECAUSE OF COMPLETE LOSSES OF SIGNIFI-
                 CANCE BY ARGUMENT REDUCTION
         IERR=5, ERROR              - NO COMPUTATION,
                 ALGORITHM TERMINATION CONDITION NOT MET

func ZBESI

func ZBESI(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, CYI []float64, NZ int, IERR int) (float64, float64, float64, int, int, []float64, []float64, int, int)

ZBESI COMPUTE THE I-BESSEL FUNCTION OF A COMPLEX ARGUMENT

FOR KODE=1 COMPUTES AN N MEMBER SEQUENCE OF COMPLEX BESSEL FUNCTIONS CY(J)=I(FNU+J-1,Z) FOR REAL, NONNEGATIVE ORDERS FNU+J-1, J=1,...,N AND COMPLEX Z IN THE CUT PLANE -PI < ARG(Z) <= PI. ON KODE=2, ZBESI RETURNS THE SCALED FUNCTIONS

CY(J)=EXP(-ABS(X))*I(FNU+J-1,Z) J = 1,...,N , X=REAL(Z)

WITH THE EXPONENTIAL GROWTH REMOVED IN BOTH THE LEFT AND RIGHT HALF PLANES FOR Z TO INFINITY. DEFINITIONS AND NOTATION ARE FOUND IN THE NBS HANDBOOK OF MATHEMATICAL FUNCTIONS (REF. 1).

INPUT ZR,ZI,FNU ARE DOUBLE PRECISION

ZR,ZI  - Z=CMPLX(ZR,ZI),  -PI < ARG(Z) <=  PI
FNU    - ORDER OF INITIAL I FUNCTION, FNU >= 0.0e0
KODE   - A PARAMETER TO INDICATE THE SCALING OPTION
         KODE= 1  RETURNS
                  CY(J)=I(FNU+J-1,Z), J=1,...,N
             = 2  RETURNS
                  CY(J)=I(FNU+J-1,Z)*EXP(-ABS(X)), J=1,...,N
N      - NUMBER OF MEMBERS OF THE SEQUENCE, N >= 1

OUTPUT CYR,CYI ARE DOUBLE PRECISION

CYR,CYI- DOUBLE PRECISION VECTORS WHOSE FIRST N COMPONENTS
         CONTAIN REAL AND IMAGINARY PARTS FOR THE SEQUENCE
         CY(J)=I(FNU+J-1,Z)  OR
         CY(J)=I(FNU+J-1,Z)*EXP(-ABS(X))  J=1,...,N
         DEPENDING ON KODE, X=REAL(Z)
NZ     - NUMBER OF COMPONENTS SET TO ZERO DUE TO UNDERFLOW,
         NZ= 0   , NORMAL RETURN
         NZ > 0 , LAST NZ COMPONENTS OF CY SET TO ZERO
                   TO UNDERFLOW, CY(J)=CMPLX(0.0e0,0.0e0)
                   J = N-NZ+1,...,N
IERR   - ERROR FLAG
         IERR=0, NORMAL RETURN - COMPUTATION COMPLETED
         IERR=1, INPUT ERROR   - NO COMPUTATION
         IERR=2, OVERFLOW      - NO COMPUTATION, REAL(Z) TOO
                 LARGE ON KODE=1
         IERR=3, CABS(Z) OR FNU+N-1 LARGE - COMPUTATION DONE
                 BUT LOSSES OF SIGNIFCANCE BY ARGUMENT
                 REDUCTION PRODUCE LESS THAN HALF OF MACHINE
                 ACCURACY
         IERR=4, CABS(Z) OR FNU+N-1 TOO LARGE - NO COMPUTA-
                 TION BECAUSE OF COMPLETE LOSSES OF SIGNIFI-
                 CANCE BY ARGUMENT REDUCTION
         IERR=5, ERROR              - NO COMPUTATION,
                 ALGORITHM TERMINATION CONDITION NOT MET

func ZBESJ

func ZBESJ(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, CYI []float64, NZ int, IERR int) (float64, float64, float64, int, int, []float64, []float64, int, int)

ZBESJ COMPUTE THE J-BESSEL FUNCTION OF A COMPLEX ARGUMENT

ON KODE=1, CBESJ COMPUTES AN N MEMBER SEQUENCE OF COMPLEX BESSEL FUNCTIONS CY(I)=J(FNU+I-1,Z) FOR REAL, NONNEGATIVE ORDERS FNU+I-1, I=1,...,N AND COMPLEX Z IN THE CUT PLANE -PI < ARG(Z) <= PI. ON KODE=2, CBESJ RETURNS THE SCALED FUNCTIONS

CY(I)=EXP(-ABS(Y))*J(FNU+I-1,Z) I = 1,...,N , Y=AIMAG(Z)

WHICH REMOVE THE EXPONENTIAL GROWTH IN BOTH THE UPPER AND LOWER HALF PLANES FOR Z TO INFINITY. DEFINITIONS AND NOTATION ARE FOUND IN THE NBS HANDBOOK OF MATHEMATICAL FUNCTIONS (REF. 1).

INPUT ZR,ZI,FNU ARE DOUBLE PRECISION

ZR,ZI  - Z=CMPLX(ZR,ZI),  -PI < ARG(Z) <=  PI
FNU    - ORDER OF INITIAL J FUNCTION, FNU >= 0.0e0
KODE   - A PARAMETER TO INDICATE THE SCALING OPTION
         KODE= 1  RETURNS
                  CY(I)=J(FNU+I-1,Z), I=1,...,N
             = 2  RETURNS
                  CY(I)=J(FNU+I-1,Z)EXP(-ABS(Y)), I=1,...,N
N      - NUMBER OF MEMBERS OF THE SEQUENCE, N >= 1

OUTPUT CYR,CYI ARE DOUBLE PRECISION

CYR,CYI- DOUBLE PRECISION VECTORS WHOSE FIRST N COMPONENTS
         CONTAIN REAL AND IMAGINARY PARTS FOR THE SEQUENCE
         CY(I)=J(FNU+I-1,Z)  OR
         CY(I)=J(FNU+I-1,Z)EXP(-ABS(Y))  I=1,...,N
         DEPENDING ON KODE, Y=AIMAG(Z).
NZ     - NUMBER OF COMPONENTS SET TO ZERO DUE TO UNDERFLOW,
         NZ= 0   , NORMAL RETURN
         NZ > 0 , LAST NZ COMPONENTS OF CY SET  ZERO DUE
                   TO UNDERFLOW, CY(I)=CMPLX(0.0e0,0.0e0),
                   I = N-NZ+1,...,N
IERR   - ERROR FLAG
         IERR=0, NORMAL RETURN - COMPUTATION COMPLETED
         IERR=1, INPUT ERROR   - NO COMPUTATION
         IERR=2, OVERFLOW      - NO COMPUTATION, AIMAG(Z)
                 TOO LARGE ON KODE=1
         IERR=3, CABS(Z) OR FNU+N-1 LARGE - COMPUTATION DONE
                 BUT LOSSES OF SIGNIFCANCE BY ARGUMENT
                 REDUCTION PRODUCE LESS THAN HALF OF MACHINE
                 ACCURACY
         IERR=4, CABS(Z) OR FNU+N-1 TOO LARGE - NO COMPUTA-
                 TION BECAUSE OF COMPLETE LOSSES OF SIGNIFI-
                 CANCE BY ARGUMENT REDUCTION
         IERR=5, ERROR              - NO COMPUTATION,
                 ALGORITHM TERMINATION CONDITION NOT MET

func ZBESK

func ZBESK(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, CYI []float64, NZ int, IERR int) (float64, float64, float64, int, int, []float64, []float64, int, int)

ZBESK COMPUTES K-BESSEL FUNCTIONS OF COMPLEX ARGUMENT

ON KODE=1, CBESK COMPUTES AN N MEMBER SEQUENCE OF COMPLEX BESSEL FUNCTIONS CY(J)=K(FNU+J-1,Z) FOR REAL, NONNEGATIVE ORDERS FNU+J-1, J=1,...,N AND COMPLEX Z != CMPLX(0.0,0.0) IN THE CUT PLANE -PI < ARG(Z) <= PI. ON KODE=2, CBESK RETURNS THE SCALED K FUNCTIONS,

CY(J)=EXP(Z)*K(FNU+J-1,Z) , J=1,...,N,

WHICH REMOVE THE EXPONENTIAL BEHAVIOR IN BOTH THE LEFT AND RIGHT HALF PLANES FOR Z TO INFINITY. DEFINITIONS AND NOTATION ARE FOUND IN THE NBS HANDBOOK OF MATHEMATICAL FUNCTIONS (REF. 1).

INPUT ZR,ZI,FNU ARE DOUBLE PRECISION

ZR,ZI  - Z=CMPLX(ZR,ZI), Z != CMPLX(0.0e0,0.0e0),
         -PI < ARG(Z) <=  PI
FNU    - ORDER OF INITIAL K FUNCTION, FNU >= 0.0e0
N      - NUMBER OF MEMBERS OF THE SEQUENCE, N >= 1
KODE   - A PARAMETER TO INDICATE THE SCALING OPTION
         KODE= 1  RETURNS
                  CY(I)=K(FNU+I-1,Z), I=1,...,N
             = 2  RETURNS
                  CY(I)=K(FNU+I-1,Z)*EXP(Z), I=1,...,N

OUTPUT CYR,CYI ARE DOUBLE PRECISION

CYR,CYI- DOUBLE PRECISION VECTORS WHOSE FIRST N COMPONENTS
         CONTAIN REAL AND IMAGINARY PARTS FOR THE SEQUENCE
         CY(I)=K(FNU+I-1,Z), I=1,...,N OR
         CY(I)=K(FNU+I-1,Z)*EXP(Z), I=1,...,N
         DEPENDING ON KODE
NZ     - NUMBER OF COMPONENTS SET TO ZERO DUE TO UNDERFLOW.
         NZ= 0   , NORMAL RETURN
         NZ > 0 , FIRST NZ COMPONENTS OF CY SET TO ZERO DUE
                   TO UNDERFLOW, CY(I)=CMPLX(0.0e0,0.0e0),
                   I=1,...,N WHEN X >= 0.0. WHEN X < 0.0
                   NZ STATES ONLY THE NUMBER OF UNDERFLOWS
                   IN THE SEQUENCE.

IERR   - ERROR FLAG
         IERR=0, NORMAL RETURN - COMPUTATION COMPLETED
         IERR=1, INPUT ERROR   - NO COMPUTATION
         IERR=2, OVERFLOW      - NO COMPUTATION, FNU IS
                 TOO LARGE OR CABS(Z) IS TOO SMALL OR BOTH
         IERR=3, CABS(Z) OR FNU+N-1 LARGE - COMPUTATION DONE
                 BUT LOSSES OF SIGNIFCANCE BY ARGUMENT
                 REDUCTION PRODUCE LESS THAN HALF OF MACHINE
                 ACCURACY
         IERR=4, CABS(Z) OR FNU+N-1 TOO LARGE - NO COMPUTA-
                 TION BECAUSE OF COMPLETE LOSSES OF SIGNIFI-
                 CANCE BY ARGUMENT REDUCTION
         IERR=5, ERROR              - NO COMPUTATION,
                 ALGORITHM TERMINATION CONDITION NOT MET

func ZBESY

func ZBESY(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, CYI []float64, NZ int, CWRKR []float64, CWRKI []float64, IERR int) (float64, float64, float64, int, int, []float64, []float64, int, []float64, []float64, int)

ZBESY COMPUTEs THE Y-BESSEL FUNCTION OF A COMPLEX ARGUMENT

ON KODE=1, CBESY COMPUTES AN N MEMBER SEQUENCE OF COMPLEX BESSEL FUNCTIONS CY(I)=Y(FNU+I-1,Z) FOR REAL, NONNEGATIVE ORDERS FNU+I-1, I=1,...,N AND COMPLEX Z IN THE CUT PLANE -PI < ARG(Z) <= PI. ON KODE=2, CBESY RETURNS THE SCALED FUNCTIONS

CY(I)=EXP(-ABS(Y))*Y(FNU+I-1,Z) I = 1,...,N , Y=AIMAG(Z)

WHICH REMOVE THE EXPONENTIAL GROWTH IN BOTH THE UPPER AND LOWER HALF PLANES FOR Z TO INFINITY. DEFINITIONS AND NOTATION ARE FOUND IN THE NBS HANDBOOK OF MATHEMATICAL FUNCTIONS (REF. 1).

INPUT ZR,ZI,FNU ARE DOUBLE PRECISION

ZR,ZI  - Z=CMPLX(ZR,ZI), Z != CMPLX(0.0e0,0.0e0),
         -PI < ARG(Z) <=  PI
FNU    - ORDER OF INITIAL Y FUNCTION, FNU >= 0.0e0
KODE   - A PARAMETER TO INDICATE THE SCALING OPTION
         KODE= 1  RETURNS
                  CY(I)=Y(FNU+I-1,Z), I=1,...,N
             = 2  RETURNS
                  CY(I)=Y(FNU+I-1,Z)*EXP(-ABS(Y)), I=1,...,N
                  WHERE Y=AIMAG(Z)
N      - NUMBER OF MEMBERS OF THE SEQUENCE, N >= 1
CWRKR, - DOUBLE PRECISION WORK VECTORS OF DIMENSION AT
CWRKI    AT LEAST N

OUTPUT CYR,CYI ARE DOUBLE PRECISION

CYR,CYI- DOUBLE PRECISION VECTORS WHOSE FIRST N COMPONENTS
         CONTAIN REAL AND IMAGINARY PARTS FOR THE SEQUENCE
         CY(I)=Y(FNU+I-1,Z)  OR
         CY(I)=Y(FNU+I-1,Z)*EXP(-ABS(Y))  I=1,...,N
         DEPENDING ON KODE.
NZ     - NZ=0 , A NORMAL RETURN
         NZ > 0 , NZ COMPONENTS OF CY SET TO ZERO DUE TO
         UNDERFLOW (GENERALLY ON KODE=2)
IERR   - ERROR FLAG
         IERR=0, NORMAL RETURN - COMPUTATION COMPLETED
         IERR=1, INPUT ERROR   - NO COMPUTATION
         IERR=2, OVERFLOW      - NO COMPUTATION, FNU IS
                 TOO LARGE OR CABS(Z) IS TOO SMALL OR BOTH
         IERR=3, CABS(Z) OR FNU+N-1 LARGE - COMPUTATION DONE
                 BUT LOSSES OF SIGNIFCANCE BY ARGUMENT
                 REDUCTION PRODUCE LESS THAN HALF OF MACHINE
                 ACCURACY
         IERR=4, CABS(Z) OR FNU+N-1 TOO LARGE - NO COMPUTA-
                 TION BECAUSE OF COMPLETE LOSSES OF SIGNIFI-
                 CANCE BY ARGUMENT REDUCTION
         IERR=5, ERROR              - NO COMPUTATION,
                 ALGORITHM TERMINATION CONDITION NOT MET

func ZBINU

func ZBINU(ZR float64, ZI float64, FNU float64, KODE int, N int, CYR []float64, CYI []float64, NZ int,
	RL float64, FNUL float64, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, []float64, []float64, int, float64, float64, float64, float64, float64)

ZBINU COMPUTES THE I FUNCTION IN THE RIGHT HALF Z PLANE

func ZBIRY

func ZBIRY(ZR float64, ZI float64, ID int, KODE int) (float64, float64, int)

ZBIRY COMPUTES THE COMPLEX AIRY FUNCTION BI(Z) OR ITS DERIVATIVE DBI(Z)/DZ ON ID=0 OR ID=1 RESPECTIVELY FOR KODE=1 FOR KODE=2, A SCALING OPTION CEXP(-AXZTA)*BI(Z) OR CEXP(-AXZTA)* DBI(Z)/DZ IS PROVIDED TO REMOVE THE EXPONENTIAL BEHAVIOR IN BOTH THE LEFT AND RIGHT HALF PLANES WHERE ZTA=(2/3)*Z*CSQRT(Z)=CMPLX(XZTA,YZTA) AND AXZTA=ABS(XZTA). DEFINTIONS AND NOTATION ARE FOUND IN THE NBS HANDBOOK OF MATHEMATICAL FUNCTIONS (REF. 1). INPUTS

ZR,ZI  - Z=CMPLX(ZR,ZI)
ID     - ORDER OF DERIVATIVE, ID=0 OR ID=1
KODE   - A PARAMETER TO INDICATE THE SCALING OPTION
         KODE= 1  return S
                  BI=BI(Z)                 ON ID=0 OR
                  BI=DBI(Z)/DZ             ON ID=1
             = 2  return S
                  BI=CEXP(-AXZTA)*BI(Z)     ON ID=0 OR
                  BI=CEXP(-AXZTA)*DBI(Z)/DZ ON ID=1 WHERE
                  ZTA=(2/3)*Z*CSQRT(Z)=CMPLX(XZTA,YZTA)
                  AND AXZTA=ABS(XZTA)

OUTPUT

BIR,BII- COMPLEX ANSWER DEPENDING ON THE CHOICES FOR ID AND
         KODE
IERR   - ERROR FLAG
         IERR=0, NORMAL return  - COMPUTATION COMPLETED
         IERR=1, INPUT ERROR   - NO COMPUTATION
         IERR=2, OVERFLOW      - NO COMPUTATION, REAL(Z)
                 TOO LARGE ON KODE=1
         IERR=3, CABS(Z) LARGE      - COMPUTATION COMPLETED
                 LOSSES OF SIGNIFCANCE BY ARGUMENT REDUCTION
                 PRODUCE LESS THAN HALF OF MACHINE ACCURACY
         IERR=4, CABS(Z) TOO LARGE  - NO COMPUTATION
                 COMPLETE LOSS OF ACCURACY BY ARGUMENT
                 REDUCTION
         IERR=5, ERROR              - NO COMPUTATION,
                 ALGORITHM TERMINATION CONDITION NOT MET

func ZBKNU

func ZBKNU(ZR float64, ZI float64, FNU float64, KODE int, N int,
	YR []float64, YI []float64, NZ int, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int,
	[]float64, []float64, int, float64, float64, float64)

ZBKNU COMPUTES THE K BESSEL FUNCTION IN THE RIGHT HALF Z PLANE.

func ZBUNI

func ZBUNI(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, YI []float64, NZ int, NUI int, NLAST int,
	FNUL float64, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, []float64, []float64, int, int, int, float64, float64, float64, float64)

ZBUNI COMPUTES THE I BESSEL FUNCTION FOR LARGE CABS(Z) > FNUL AND FNU+N-1 < FNUL. THE ORDER IS INCREASED FROM FNU+N-1 GREATER THAN FNUL BY ADDING NUI AND COMPUTING ACCORDING TO THE UNIFORM ASYMPTOTIC EXPANSION FOR I(FNU,Z) ON IFORM=1 AND THE EXPANSION FOR J(FNU,Z) ON IFORM=2

func ZBUNK

func ZBUNK(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, YI []float64, NZ int,
	TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, int, []float64, []float64, int,
	float64, float64, float64)

ZBUNK COMPUTES THE K BESSEL FUNCTION FOR FNU > FNUL. ACCORDING TO THE UNIFORM ASYMPTOTIC EXPANSION FOR K(FNU,Z) IN ZUNK1 AND THE EXPANSION FOR H(2,FNU,Z) IN ZUNK2

func ZDIV

func ZDIV(AR float64, AI float64, BR float64, BI float64) (float64, float64)

ZDIV returns complex divide of a complex number real and imag parts

func ZEXP

func ZEXP(AR float64, AI float64) (float64, float64)

ZEXP returns complex exponential function of a complex number real and imag parts

func ZFLOATS

func ZFLOATS(Z complex128) (float64, float64)

ZFLOATS splits a complex number into real and imag parts

func ZKSCL

func ZKSCL(ZRR float64, ZRI float64, FNU float64, N int, YR []float64, YI []float64, NZ int,
	RZR float64, RZI float64, ASCLE float64, TOL float64, ELIM float64) (float64, float64, float64, int, []float64, []float64, int,
	float64, float64, float64, float64, float64)

ZKSCL SET K FUNCTIONS TO ZERO ON UNDERFLOW, CONTINUE RECURRENCE ON SCALED FUNCTIONS UNTIL TWO MEMBERS COME ON SCALE, THEN RETURN WITH MIN(NZ+2,N) VALUES SCALED BY 1/TOL.

func ZLOG

func ZLOG(AR float64, AI float64) (float64, float64)

ZLOG returns complex logarithm of a complex number real and imag parts

func ZMLRI

func ZMLRI(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, YI []float64,
	NZ int, TOL float64) (float64, float64, float64, int, int, []float64, []float64, int, float64)

ZMLRI COMPUTES THE I BESSEL FUNCTION FOR RE(Z) >= 0.0 BY THE MILLER ALGORITHM NORMALIZED BY A NEUMANN SERIES.

func ZMLT

func ZMLT(AR float64, AI float64, BR float64, BI float64) (float64, float64)

ZMLT returns complex multipliTION of complex numberSs real and imag parts

func ZRATI

func ZRATI(ZR float64, ZI float64, FNU float64, N int, CYR []float64, CYI []float64, TOL float64) (float64, float64, float64, int, []float64, []float64, float64)

ZRATI COMPUTES RATIOS OF I BESSEL FUNCTIONS BY BACKWARD RECURRENCE. THE STARTING INDEX IS DETERMINED BY FORWARD RECURRENCE AS DESCRIBED IN J. RES. OF NAT. BUR. OF STANDARDS-B, MATHEMATICAL SCIENCES, VOL 77B, P111-114, SEPTEMBER, 1973, BESSEL FUNCTIONS I AND J OF COMPLEX ARGUMENT AND INTEGER ORDER,

BY D. J. SOOKNE.

func ZS1S2

func ZS1S2(ZRR float64, ZRI float64, S1R float64, S1I float64, S2R float64, S2I float64, NZ int, ASCLE float64,
	ALIM float64, IUF int) (float64, float64, float64, float64, float64, float64, int, float64, float64, int)

ZS1S2 tests for a possible underflow resulting from the addition of the I and K functions in the analytic continuation formula where S1=K function and S2=I function. For KODE=1 the I and K functions are different orders of magnitude, but for KODE=2 they can be of the same order of magnitude and the maximum must be at least one precision above the underflow limit.

func ZSERI

func ZSERI(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, YI []float64, NZ int,
	TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, []float64, []float64, int,
	float64, float64, float64)

ZSERI COMPUTES THE I BESSEL FUNCTION FOR REAL(Z) >= 0.0 BY MEANS OF THE POWER SERIES FOR LARGE CABS(Z) IN THE REGION CABS(Z) <= 2*SQRT(FNU+1). NZ=0 IS A NORMAL RETURN. NZ > 0 MEANS THAT THE LAST NZ COMPONENTS WERE SET TO ZERO DUE TO UNDERFLOW. NZ < 0 MEANS UNDERFLOW OCCURRED, BUT THE CONDITION CABS(Z) <= 2*SQRT(FNU+1) WAS VIOLATED AND THE COMPUTATION MUST BE COMPLETED IN ANOTHER ROUTINE WITH N=N-ABS(NZ).

func ZSHCH

func ZSHCH(ZR float64, ZI float64) (float64, float64, float64, float64)

ZSHCH computes the complex hyperbolic functions SINH(x+i*y) and COSH(x+i*y) returns real and imag parts of sinh and cosh

func ZSQRT

func ZSQRT(AR float64, AI float64) (float64, float64)

ZSQRT returns complex square root of a complex number real and imag parts

func ZUCHK

func ZUCHK(YR float64, YI float64, NZ int, ASCLE float64, TOL float64) (float64, float64, int, float64, float64)

ZUCHK takes a scaled quantity Y whose magnitude is greater than

EXP(-ALIM)=ASCLE=1000*D1MACH[1]/TOL

and tests if the magnitude of the real or imaginary part would underflow when y is scaled by TOL to its proper value. Y is accepted if the underflow is at least one precision below the magnitude of the largest component; otherwise the phase angle does not have absolute accuracy and an underflow is assumed.

func ZUNHJ

func ZUNHJ(ZR float64, ZI float64, FNU float64, IPMTR int, TOL float64, PHIR float64, PHII float64,
	ARGR float64, ARGI float64, ZETA1R float64, ZETA1I float64, ZETA2R float64, ZETA2I float64,
	ASUMR float64, ASUMI float64, BSUMR float64, BSUMI float64) (float64, float64, float64, int, float64, float64, float64,
	float64, float64, float64, float64, float64, float64, float64, float64, float64, float64)

ZUNHJ COMPUTES PARAMETERS FOR BESSEL FUNCTIONS C(FNU,Z) = J(FNU,Z), Y(FNU,Z) OR H(I,FNU,Z) I=1,2 FOR LARGE ORDERS FNU BY MEANS OF THE UNIFORM ASYMPTOTIC EXPANSION

C(FNU,Z)=C1*PHI*( ASUM*AIRY(ARG) + C2*BSUM*DAIRY(ARG) )

FOR PROPER CHOICES OF C1, C2, AIRY AND DAIRY WHERE AIRY IS AN AIRY FUNCTION AND DAIRY IS ITS DERIVATIVE.

(2/3)*FNU*ZETA**1.5 = ZETA1-ZETA2,

ZETA1=0.5*FNU*CLOG((1+W)/(1-W)), ZETA2=FNU*W FOR SCALING PURPOSES IN AIRY FUNCTIONS FROM CAIRY OR CBIRY.

MCONJ=SIGN OF AIMAG(Z), BUT IS AMBIGUOUS WHEN Z IS REAL AND MUST BE SPECIFIED. IPMTR=0 RETURNS ALL PARAMETERS. IPMTR= 1 COMPUTES ALL EXCEPT ASUM AND BSUM.

HANDBOOK OF MATHEMATICAL FUNCTIONS BY M. ABRAMOWITZ AND I.A. STEGUN, AMS55, NATIONAL BUREAU OF STANDARDS, 1965, CHAPTER 9. ASYMPTOTICS AND SPECIAL FUNCTIONS BY F.W.J. OLVER, ACADEMIC PRESS, N.Y., 1974, PAGE 420

func ZUNI1

func ZUNI1(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, YI []float64, NZ int, NLAST int, FNUL float64, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, []float64, []float64, int, int, float64, float64, float64, float64)

ZUNI1 COMPUTES I(FNU,Z) BY MEANS OF THE UNIFORM ASYMPTOTIC EXPANSION FOR I(FNU,Z) IN -PI/3 <= ARG Z <= PI/3.

FNUL IS THE SMALLEST ORDER PERMITTED FOR THE ASYMPTOTIC EXPANSION. NLAST=0 MEANS ALL OF THE Y VALUES WERE SET. NLAST.NE.0 IS THE NUMBER LEFT TO BE COMPUTED BY ANOTHER FORMULA FOR ORDERS FNU TO FNU+NLAST-1 BECAUSE FNU+NLAST-1 < FNUL. Y(I)=CZERO FOR I=NLAST+1,N

func ZUNI2

func ZUNI2(ZR float64, ZI float64, FNU float64, KODE int, N int, YR []float64, YI []float64, NZ int, NLAST int, FNUL float64, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, []float64, []float64, int, int, float64, float64, float64, float64)

ZUNI2 COMPUTES I(FNU,Z) IN THE RIGHT HALF PLANE BY MEANS OF UNIFORM ASYMPTOTIC EXPANSION FOR J(FNU,ZN) WHERE ZN IS Z*I OR -Z*I AND ZN IS IN THE RIGHT HALF PLANE ALSO.

FNUL IS THE SMALLEST ORDER PERMITTED FOR THE ASYMPTOTIC EXPANSION. NLAST=0 MEANS ALL OF THE Y VALUES WERE SET. NLAST != 0 IS THE NUMBER LEFT TO BE COMPUTED BY ANOTHER FORMULA FOR ORDERS FNU TO FNU+NLAST-1 BECAUSE FNU+NLAST-1 < FNUL. Y(I)=CZERO FOR I=NLAST+1,N

func ZUNIK

func ZUNIK(ZRR float64, ZRI float64, FNU float64, IKFLG int, IPMTR int, TOL float64, INIT int, PHIR float64,
	PHII float64, ZETA1R float64, ZETA1I float64, ZETA2R float64, ZETA2I float64,
	SUMR float64, SUMI float64, CWRKR []float64, CWRKI []float64) (float64, float64, float64, int, int, float64, int, float64,
	float64, float64, float64, float64, float64, float64, float64, []float64, []float64)

ZUNIK COMPUTES PARAMETERS FOR THE UNIFORM ASYMPTOTIC EXPANSIONS OF THE I AND K FUNCTIONS ON IKFLG= 1 OR 2 RESPECTIVELY BY

W(FNU,ZR) = PHI*EXP(ZETA)*SUM

WHERE ZETA=-ZETA1 + ZETA2 OR

ZETA1 - ZETA2

THE FIRST CALL MUST HAVE INIT=0. SUBSEQUENT CALLS WITH THE SAME ZR AND FNU WILL return THE I OR K FUNCTION ON IKFLG= 1 OR 2 WITH NO CHANGE IN INIT. CWRK IS A COMPLEX WORK ARRAY. IPMTR=0 COMPUTES ALL PARAMETERS. IPMTR=1 COMPUTES PHI, ZETA1,ZETA2.

func ZUNK1

func ZUNK1(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, YI []float64, NZ int, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, float64, float64)

ZUNK1 COMPUTES K(FNU,Z) AND ITS ANALYTIC CONTINUATION FROM THE RIGHT HALF PLANE TO THE LEFT HALF PLANE BY MEANS OF THE UNIFORM ASYMPTOTIC EXPANSION. MR INDICATES THE DIRECTION OF ROTATION FOR ANALYTIC CONTINUATION. NZ=-1 MEANS AN OVERFLOW WILL OCCUR

func ZUNK2

func ZUNK2(ZR float64, ZI float64, FNU float64, KODE int, MR int, N int, YR []float64, YI []float64, NZ int, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, int, []float64, []float64, int, float64, float64, float64)

ZUNK2 COMPUTES K(FNU,Z) AND ITS ANALYTIC CONTINUATION FROM THE RIGHT HALF PLANE TO THE LEFT HALF PLANE BY MEANS OF THE UNIFORM ASYMPTOTIC EXPANSIONS FOR H(KIND,FNU,ZN) AND J(FNU,ZN) WHERE ZN IS IN THE RIGHT HALF PLANE, KIND=(3-MR)/2, MR=+1 OR -1. HERE ZN=ZR*I OR -ZR*I WHERE ZR=Z IF Z IS IN THE RIGHT HALF PLANE OR ZR=-Z IF Z IS IN THE LEFT HALF PLANE. MR INDIC- ATES THE DIRECTION OF ROTATION FOR ANALYTIC CONTINUATION. NZ=-1 MEANS AN OVERFLOW WILL OCCUR

func ZUOIK

func ZUOIK(ZR float64, ZI float64, FNU float64, KODE int, IKFLG int, N int, YR []float64, YI []float64, NUF int, TOL float64, ELIM float64, ALIM float64) (
	float64, float64, float64, int, int, int, []float64, []float64, int, float64, float64, float64)

ZUOIK COMPUTES THE LEADING TERMS OF THE UNIFORM ASYMPTOTIC EXPANSIONS FOR THE I AND K FUNCTIONS AND COMPARES THEM (IN LOGARITHMIC FORM) TO ALIM AND ELIM FOR OVER AND UNDERFLOW WHERE ALIM < ELIM. if THE MAGNITUDE, BASED ON THE LEADING EXPONENTIAL, IS LESS THAN ALIM OR GREATER THAN -ALIM, THEN THE RESULT IS ON SCALE. if NOT, THEN A REFINED TEST USING OTHER MULTIPLIERS (IN LOGARITHMIC FORM) IS MADE BASED ON ELIM. HERE EXP(-ELIM)=SMALLEST MACHINE NUMBER*1.0E+3 AND EXP(-ALIM)= EXP(-ELIM)/TOL

IKFLG=1 MEANS THE I SEQUENCE IS TESTED

=2 MEANS THE K SEQUENCE IS TESTED

NUF = 0 MEANS THE LAST MEMBER OF THE SEQUENCE IS ON SCALE

=-1 MEANS AN OVERFLOW WOULD OCCUR

IKFLG=1 AND NUF > 0 MEANS THE LAST NUF Y VALUES WERE SET TO ZERO

THE FIRST N-NUF VALUES MUST BE SET BY ANOTHER ROUTINE

IKFLG=2 AND NUF == N MEANS ALL Y VALUES WERE SET TO ZERO IKFLG=2 AND 0 < NUF < N NOT CONSIDERED. Y MUST BE SET BY

ANOTHER ROUTINE

func ZWRSK

func ZWRSK(ZRR float64, ZRI float64, FNU float64, KODE int, N int, YR []float64, YI []float64, NZ int, CWR []float64, CWI []float64, TOL float64, ELIM float64, ALIM float64) (float64, float64, float64, int, int, []float64, []float64, int, []float64, []float64, float64, float64, float64)

ZWRSK COMPUTES THE I BESSEL FUNCTION FOR RE(Z) >= 0.0 BY NORMALIZING THE I FUNCTION RATIOS FROM ZRATI BY THE WRONSKIAN