C specifications:

— function: void acmlversion (int *major, int *minor, int *patch);
— function: void acmlinfo (void);

CFFT1D Routine Documentation

SUBROUTINE: CFFT1D (MODE,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CFFT1D.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward transform is performed. Initializations are assumed to have been performed by a prior call to CFFT1D.
  • MODE=1 : a backward (reverse) transform is performed. Initializations are assumed to have been performed by a prior call to CFFT1D.
  • MODE=-2 : (default) initializations and a forward transform are performed.
  • MODE=2 : (default) initializations and a backward transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER N

On input: N is the length of the complex sequence X

— Input/Output: COMPLEX X(N)

On input: X contains the complex sequence of length N to be transformed.
On output: X contains the transformed sequence.

— Input/Output: COMPLEX COMM(5*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL CFFT1D(0,N,X,COMM,INFO)
             CALL CFFT1D(-1,N,X,COMM,INFO)
             CALL CFFT1D(-1,N,Y,COMM,INFO)
             DO 10 I = 1, N
                X(I) = X(I)*CONJG(Y(I))
        10   CONTINUE
             CALL CFFT1D(1,N,X,COMM,INFO)

CFFT1DX Routine Documentation

SUBROUTINE: CFFT1DX (MODE,SCALE,INPL,N,X,INCX,Y,INCY,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CFFT1DX.
On input:

  • MODE=0 : only initializations (specific to the value of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward transform is performed. Initializations are assumed to have been performed by a prior call to CFFT1DX.
  • MODE=1 : a backward (reverse) transform is performed. Initializations are assumed to have been performed by a prior call to CFFT1DX.
  • MODE=-2 : (default) initializations and a forward transform are performed.
  • MODE=2 : (default) initializations and a backward transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: REAL SCALE

On input: SCALE is the scaling factor to apply to the output sequence

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequence; otherwise the output sequence is returned in Y.

— Input: INTEGER N

On input: N is the number of elements to be transformed

— Input/Output: COMPLEX X(1+(N-1)*INCX)

On input: X contains the complex sequence of length N to be transformed, with the ith element stored in X(1+(i-1)*INCX).
On output: if INPL is .TRUE. then X contains the transformed sequence in the same locations as on input; otherwise X remains unchanged.

— Input: INTEGER INCX

On input: INCX is the increment used to store successive elements of a sequence in X.
Constraint: INCX > 0.

— Output: COMPLEX Y(1+(N-1)*INCY)

On output: if INPL is .FALSE. then Y contains the transformed sequence, with the ith element stored in Y(1+(i-1)*INCY); otherwise Y is not referenced.

— Input: INTEGER INCY

On input: INCY is the increment used to store successive elements of a sequence in Y. If INPL is .TRUE. then INCY is not referenced.
Constraint: INCY > 0.

— Input/Output: COMPLEX COMM(5*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward FFTs are performed unscaled and in-place on contiguous
     C     vectors X and Y following initialization. Manipulations on
     C     resultant Fourier coefficients are stored in X which is then
     C     transformed back.
     C
             SCALE = 1.0
             INPL = .TRUE.
             CALL CFFT1DX(0,SCALE,INPL,N,X,1,DUM,1,COMM,INFO)
             CALL CFFT1DX(-1,SCALE,INPL,N,X,1,DUM,1,COMM,INFO)
             CALL CFFT1DX(-1,SCALE,INPL,N,Y,1,DUM,1,COMM,INFO)
             DO 10 I = 1, N
                X(I) = X(I)*CONJG(Y(I))/REAL(N)
        10   CONTINUE
             CALL CFFT1DX(1,SCALE,INPL,N,X,1,DUM,1,COMM,INFO)

CFFT1M Routine Documentation

SUBROUTINE:CFFT1M (MODE,M,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CFFT1M.
On input:

  • MODE=0 : only initializations (specific to the value of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : forward transforms are performed. Initializations are assumed to have been performed by a prior call to CFFT1M.
  • MODE=1 : backward (reverse) transforms are performed. Initializations are assumed to have been performed by a prior call to CFFT1M.
  • MODE=-2 : (default) initializations and forward transforms are performed.
  • MODE=2 : (default) initializations and backward transforms are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is the length of the complex sequences in X

— Input/Output: COMPLEX X(N*M)

On input: X contains the M complex sequences of length N to be transformed. Element i of sequence j is stored in location i+(j-1)*N of X.
On output: X contains the transformed sequences.

— Input/Output: COMPLEX COMM(5*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL CFFT1M(0,1,N,X,COMM,INFO)
             CALL CFFT1M(-1,2,N,X,COMM,INFO)
             DO 10 I = 1, N
                X(I,3) = X(I,1)*CONJG(X(I,2))
                X(I,2) = CMPLX(0.0D0,1.0D0)*X(I,2)
        10   CONTINUE
             CALL CFFT1M(1,2,N,X(1,2),COMM,INFO)

CFFT1MX Routine Documentation

SUBROUTINE:CFFT1MX (MODE,SCALE,INPL,NSEQ,N,X,INCX1,INCX2,Y,INCY1,INCY2,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CFFT1MX.
On input:

  • MODE=0 : only initializations (specific to the value of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward transform is performed. Initializations are assumed to have been performed by a prior call to CFFT1MX.
  • MODE=1 : a backward (reverse) transform is performed. Initializations are assumed to have been performed by a prior call to CFFT1MX.
  • MODE=-2 : (default) initializations and a forward transform are performed.
  • MODE=2 : (default) initializations and a backward transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: REAL SCALE

On input: SCALE is the scaling factor to apply to the output sequences

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequences; otherwise the output sequences are returned in Y.

— Input: INTEGER NSEQ

On input: NSEQ is the number of sequences to be transformed

— Input: INTEGER N

On input: N is the number of elements in each sequence to be transformed

— Input/Output: COMPLEX X(1+(N-1)*INCX1+(NSEQ-1)*INCX2)

On input: X contains the NSEQ complex sequences of length N to be transformed; the ith element of sequence j is stored in X(1+(i-1)*INCX1+(j-1)*INCX2).
On output: if INPL is .TRUE. then X contains the transformed sequences in the same locations as on input; otherwise X remains unchanged.

— Input: INTEGER INCX1

On input: INCX1 is the increment used to store successive elements of a given sequence in X (INCX1=1 for contiguous data).
Constraint: INCX1 > 0.

— Input: INTEGER INCX2

On input: INCX2 is the increment used to store corresponding elements of successive sequences in X (INCX2=N for contiguous data).
Constraint: INCX2 > 0.

— Output: COMPLEX Y(1+(N-1)*INCY1+(NSEQ-1)*INCY2)

On output: if INPL is .FALSE. then Y contains the transformed sequences with the ith element of sequence j stored in Y(1+(i-1)*INCY1+(j-1)*INCY2); otherwise Y is not referenced.

— Input: INTEGER INCY1

On input: INCY1 is the increment used to store successive elements of a given sequence in Y. If INPL is .TRUE. then INCY1 is not referenced.
Constraint: INCY1 > 0.

— Input: INTEGER INCY2

On input: INCY2 is the increment used to store corresponding elements of successive sequences in Y (INCY2=N for contiguous data). If INPL is .TRUE. then INCY2 is not referenced.
Constraint: INCY2 > 0.

— Input/Output: COMPLEX COMM(5*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward FFTs are performed unscaled and in-place on two
     C     contiguous vectors stored in the first two columns of X.
     C     Manipulations are stored in 2nd and 3rd columns of X which are
     C     then transformed back.
     C
             COMPLEX X(N,3)
             SCALE = 1.0
             INPL = .TRUE.
             CALL CFFT1MX(0,SCALE,INPL,2,N,X,1,N,DUM,1,N,COMM,INFO)
             CALL CFFT1MX(-1,SCALE,INPL,2,N,X,1,N,DUM,1,N,COMM,INFO)
             DO 10 I = 1, N
                X(I,3) = X(I,1)*CONJG(X(I,2))/REAL(N)
                X(I,2) = CMPLX(0.0D0,1.0D0)*X(I,2)/REAL(N)
        10   CONTINUE
             CALL CFFT1MX(1,SCALE,INPL,2,N,X(1,2),1,N,DUM,1,N,COMM,INFO)

CFFT2D Routine Documentation

SUBROUTINE:CFFT2D (MODE,M,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the direction of transform to be performed by CFFT2D.
On input:

  • MODE=-1 : a forward 2D transform is performed.
  • MODE=1 : a backward (reverse) 2D transform is performed.

— Input: INTEGER M

On input: M is the number of rows in the 2D array of data to be transformed. If X is declared as a 2D array then M is the first dimension of X.

— Input: INTEGER N

On input: N is the number of columns in the 2D array of data to be transformed. If X is declared as a 2D array then M is the second dimension of X.

— Input/Output: COMPLEX X(M*N)

On input: X contains the M by N complex 2D array to be transformed. Element ij is stored in location i+(j-1)*M of X.
On output: X contains the transformed sequence.

— Input/Output: COMPLEX COMM(M*N+5*(M+N))

COMM is a communication array used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL CFFT2D(-1,M,N,X,COMM,INFO)
             DO 20 J = 1, N
                DO 10 I = 1, MIN(J-1,M)
                   X(I,J) = 0.5D0*(X(I,J) + X(J,I))
                   X(J,I) = CONJG(X(I,J))
        10      CONTINUE
        20   CONTINUE
             CALL CFFT2D(1,M,N,X,COMM,INFO)

CFFT2DX Routine Documentation

SUBROUTINE:CFFT2DX (MODE,SCALE,LTRANS,INPL,M,N,X,INCX1,INCX2,Y,INCY1,INCY2,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CFFT2DX.
On input:

  • MODE=0 : only initializations (specific to the value of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward 2D transform is performed. Initializations are assumed to have been performed by a prior call to CFFT2DX.
  • MODE=1 : a backward (reverse) 2D transform is performed. Initializations are assumed to have been performed by a prior call to CFFT2DX.
  • MODE=-2 : (default) initializations and a forward 2D transform are performed.
  • MODE=2 : (default) initializations and a backward 2D transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the values of N and M) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the values of N and M.

— Input: REAL SCALE

On input: SCALE is the scaling factor to apply to the output sequences

— Input: LOGICAL LTRANS

On input: if LTRANS is .TRUE. then a normal final transposition is performed internally to return transformed data consistent with the values for arguments INPL, INCX1, INCX2, INCY1 and INCY2. If LTRANS is .FALSE. then the final transposition is not performed explicitly; the storage format on output is determined by whether the output data is stored contiguously or not – please see the output specifications for X and Y for details.

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequences; otherwise the output sequences are returned in Y.

— Input: INTEGER M

On input: M is the first dimension of the 2D transform.

— Input: INTEGER N

On input: N is the second dimension of the 2D transform.

— Input/Output: COMPLEX X(1+(M-1)*INCX1+(N-1)*INCX2)

On input: X contains the M by N complex 2D data array to be transformed; the (ij)th element is stored in X(1+(i-1)*INCX1+(j-1)*INCX2).
On output: if INPL is .TRUE. then X contains the transformed data, either in the same locations as on input when LTRANS=.TRUE.; in locations X((i-1)*N+j) when LTRANS=.FALSE., INCX1=1 and INCX2=M; and otherwise in the same locations as on input. If INPL is .FALSE. X remains unchanged.

— Input: INTEGER INCX1

On input: INCX1 is the increment used to store, in X, successive elements in the first dimension (INCX1=1 for contiguous data).
Constraint: INCX1 > 0.

— Input: INTEGER INCX2

On input: INCX2 is the increment used to store, in X, successive elements in the second dimension (INCX2=M for contiguous data).
Constraint: INCX2 > 0;
INCX2 > (M-1)*INCX1 if N > 1.

— Output: COMPLEX Y(1+(M-1)*INCY1+(N-1)*INCY2)

On output: if INPL is .FALSE. then Y contains the transformed data. If LTRANS=.TRUE. then the (ij)th data element is stored in Y(1+(i-1)*INCY1+(j-1)*INCY2); if LTRANS=.FALSE., INCY1=1 and INCY2=N then the (ij)th data element is stored in Y((i-1)*N+j); and otherwise the (ij)th element is stored in Y(1+(i-1)*INCY1+(j-1)*INCY2). If INPL is .TRUE. then Y is not referenced.

— Input: INTEGER INCY1

On input: INCY1 is the increment used to store successive elements in the first dimension in Y (INCY1=1 for contiguous data). If INPL is .TRUE. then INCY1 is not referenced.
Constraint: INCY1 > 0.

— Input: INTEGER INCY2

On input: INCY2 is the increment used to store successive elements in the second dimension in Y (for contiguous data, INCY2=M when LTRANS is .TRUE. or INCY2=N when LTRANS is .FALSE.). If INPL is .TRUE. then INCY2 is not referenced.
Constraints: INCY2 > 0;
INCY2 > (M-1)*INCY1 if N > 1 and LTRANS is .TRUE.;
INCY2 = N if M > 1 and LTRANS is .FALSE..

— Input/Output: COMPLEX COMM(M*N+5*M+5*N+200)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same dimensions M and N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward 2D FFT is performed unscaled, without final transpose
     C     and out-of-place on data stored in array X and output to Y.
     C     Manipulations are stored in vector Y which is then transformed
     C     back, with scaling, into the first M rows of X.
     C
             COMPLEX X(M,N), Y(N,M)
             SCALE = 1.0
             INPL = .FALSE.
             LTRANS = .FALSE.
             CALL CFFT2DX(0,SCALE,LTRANS,INPL,M,N,X,1,M,Y,1,N,COMM,INFO)
             CALL CFFT2DX(-1,SCALE,LTRANS,INPL,M,N,X,1,M,Y,1,N,COMM,INFO)
             DO 20 I = M
                DO 10 J = 1, N
                   Y(J,I) = 0.5*Y(J,I)*EXP(-0.001*REAL(I+J-2))
                   IY = IY + 1
        10      CONTINUE
        20   CONTINUE
             SCALE = 1.0/REAL(M*N)
             CALL CFFT2DX(1,SCALE,LTRANS,INPL,N,M,Y,1,N,X,1,M,COMM,INFO)

CFFT3D Routine Documentation

SUBROUTINE:CFFT3D (MODE,L,M,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the direction of transform to be performed by CFFT3D.
On input:

  • MODE=-1 : forward 3D transform is performed.
  • MODE=1 : backward (reverse) 3D transform is performed.

— Input: INTEGER L

On input: the length of the first dimension of the 3D array of data to be transformed. If X is declared as a 3D array then L is the first dimension of X.

— Input: INTEGER M

On input: the length of the second dimension of the 3D array of data to be transformed. If X is declared as a 3D array then M is the second dimension of X.

— Input: INTEGER N

On input: the length of the third dimension of the 3D array of data to be transformed. If X is declared as a 3D array then N is the third dimension of X.

— Input/Output: COMPLEX X(L*M*N)

On input: X contains the L by M by N complex 3D array to be transformed. Element ijk is stored in location i+(j-1)*L+(k-1)*L*M of X.
On output: X contains the transformed sequence.

— Input/Output: COMPLEX COMM(5*(L+M+N)+4)

COMM is a communication array used as temporary store. Note that the amount of store explicitly required here is less than in some versions prior to this release (version 4.1 and older). Some further workspace will be allocated internally; the amount of allocated memory requested will be MAX(L*N+N, L*M + L + M, N).

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL CFFT3D(-1,L,M,N,X,COMM,INFO)
             DO 30 K = 1, N
               DO 20 J = 1, M
                 DO 10 I = 1, L
                   X(I,J) = X(I,J)*EXP(-0.001D0*REAL(I+J+K))
        10       CONTINUE
        20     CONTINUE
        30   CONTINUE
             CALL CFFT3D(1,L,M,N,X,COMM,INFO)

CFFT3DX Routine Documentation

SUBROUTINE:CFFT3DX (MODE,SCALE,LTRANS,INPL,L,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CFFT3DX.
On input:

  • MODE=0 : only initializations (specific to the values of L, M and N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward 3D transform is performed. Initializations are assumed to have been performed by a prior call to CFFT3DX.
  • MODE=1 : a backward (reverse) 3D transform is performed. Initializations are assumed to have been performed by a prior call to CFFT3DX.
  • MODE=100 : similar to MODE=0; only initializations (specific to the values of L, M and M) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the values of L, M and N.

— Input: REAL SCALE

On input: SCALE is the scaling factor to apply to the output sequences

— Input: LOGICAL LTRANS

On input: if LTRANS is .TRUE. then a normal final transposition is performed internally to return transformed data using the same storage format as the input data. If LTRANS is .FALSE. then the final transposition is not performed and transformed data is stored, in X or Y, in transposed form.

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequences; otherwise the output sequences are returned in Y.

— Input: INTEGER L

On input: L is the first dimension of the 3D transform.

— Input: INTEGER M

On input: M is the second dimension of the 3D transform.

— Input: INTEGER N

On input: N is the third dimension of the 3D transform.

— Input/Output: COMPLEX X(L*M*N)

On input: X contains the L by M by N complex 3D data array to be transformed; the (ijk)th element is stored in X(i+(j-1)*L+(k-1)*L*M).
On output: if INPL is .TRUE. then X contains the transformed data, either in the same locations as on input when LTRANS=.TRUE.; or in locations X(k+(j-1)*N+(i-1)*N*M) when LTRANS=.FALSE. If INPL is .FALSE. X remains unchanged.

— Output: COMPLEX Y(L*M*N)

On output: if INPL is .FALSE. then Y contains the three-dimensional transformed data. If LTRANS=.TRUE. then the (ijk)th data element is stored in Y(i+(j-1)*L+(k-1)*L*M); otherwise, the (ijk)th data element is stored in Y(k+(j-1)*N+(k-1)*N*M). If INPL is .TRUE. then Y is not referenced.

— Input/Output: COMPLEX COMM(*)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence dimensions. The remainder is used as temporary store. The amount of store required depends on the values of the arguments L, M, N and LTRANS. If LTRANS=.TRUE. then for a genuine 3D transform (all of L, M, N greater than 1) the dimension of COMM need only be 5*(L+M+N) + 150; in this case some further workspace will be allocated internally; the amount of allocated memory requested will be MAX(L*N+N, L*M + L + M, N). If LTRANS=.FALSE. then for a genuine 3D transform the workspace requirement is considerably more to allow for the storage of an intermediate 3D transposed array; in this case the dimension of COMM must be at least L*M*N+5*(L+M+N)+150. It is recommended that the appropriate 1D or 2D FFT routine be called when at least one of L, M or N is 1.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward 3D FFT is performed unscaled, without final transpose
     C     and out-of-place on data stored in array X and output to Y.
     C     Manipulations are stored in vector Y which is then transformed
     C     back, with scaling, into the first M rows of X.
     C
             SCALE = 1.0
             INPL = .FALSE.
             LTRANS = .FALSE.
             CALL CFFT3DX(0,SCALE,LTRANS,INPL,L,M,N,X,Y,COMM,INFO)
             CALL CFFT3DX(-1,SCALE,LTRANS,INPL,L,M,N,X,Y,COMM,INFO)
             IY = 1
             DO 20 I = 1, L
                DO 40 J = 1, M
                   DO 10 K = 1, N
                      Y(IY) = Y(IY)*EXP(-0.001*REAL(I+J+K-3))
                      IY = IY + 1
        10      CONTINUE
        20   CONTINUE
             SCALE = 1.0/REAL(L*M*N)
             CALL CFFT3DX(1,SCALE,LTRANS,INPL,N,M,L,Y,X,COMM,INFO)

CFFT3DY Routine Documentation

SUBROUTINE:CFFT3DY (MODE,SCALE,INPL,L,M,N,X,INCX1,INCX2,INCX3,Y,INCY1,INCY2,INCY3,COMM,LCOMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CFFT3DY.
On input:

  • MODE=0 : only initializations (specific to the values of L, M and N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward 3D transform is performed. Initializations are assumed to have been performed by a prior call to CFFT3DY.
  • MODE=1 : a backward (reverse) 3D transform is performed. Initializations are assumed to have been performed by a prior call to CFFT3DY.
  • MODE=-2 : (default) initializations and a forward 3D transform are performed.
  • MODE=2 : (default) initializations and a backward 3D transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the values of L, M and M) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the values of L, M and N.

— Input: REAL SCALE

On input: SCALE is the scaling factor to apply to the output sequences

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequences; otherwise the output sequences are returned in Y.

— Input: INTEGER L

On input: L is the first dimension of the 3D transform.

— Input: INTEGER M

On input: M is the second dimension of the 3D transform.

— Input: INTEGER N

On input: N is the third dimension of the 3D transform.

— Input/Output: COMPLEX X(*)

On input: X contains the L by M by N complex 3D data array to be transformed; the (ijk)th element is stored in X(1+(i-1)*INCX1+(j-1)*INCX2+(k-1)*INCX3).
On output: if INPL is .TRUE. then X contains the transformed data in the same locations as on input. If INPL is .FALSE. X remains unchanged.

— Input: INTEGER INCX1

On input: INCX1 is the step in index of X between successive data elements in the first dimension of the 3D data. Usually INCX1=1 so that succesive elements in the first dimension are stored contiguously.
Constraint: INCX1 > 0.

— Input: INTEGER INCX2

On input: INCX2 is the step in index of X between successive data elements in the second dimension of the 3D data. For completely contiguous data (no gaps in X) INCX2 should be set to L.
Constraint: INCX2 > 0;
INCX2 > (L-1)*INCX1 if max(M,N) > 1.

— Input: INTEGER INCX3

On input: INCX3 is the step in index of X between successive data elements in the third dimension of the 3D data. For completely contiguous data (no gaps in X) INCX3 should be set to L*M.
Constraint: INCX3 > 0;
INCX3 > (L-1)*INCX1+(M-1)*INCX2 if N > 1.

— Output: COMPLEX Y(*)

On output: if INPL is .FALSE. then Y contains the three-dimensional transformed data. If LTRANS=.TRUE. then the the (ijk)th element is stored in Y(1+(i-1)*INCY1+(j-1)*INCY2+(k-1)*INCY3).
If INPL is .TRUE. then Y is not referenced.

— Input: INTEGER INCY1

On input: if INPL is .FALSE. then INCY1 is the step in index of Y between successive data elements in the first dimension of the 3D transformed data. Usually INCY1=1 so that succesive elements in the first dimension are stored contiguously.
If INPL is .TRUE. then INCY1 is not referenced. Constraint: If INPL is .FALSE. then INCY1 > 0.

— Input: INTEGER INCY2

On input: if INPL is .FALSE. then INCY2 is the step in index of Y between successive data elements in the second dimension of the 3D transformed data. For completely contiguous data (no gaps in Y) INCY2 should be set to L.
Constraint: INCY2 > 0 if INPL is .FALSE.;
INCY2 > (L-1)*INCY1, if INPL is .FALSE. and max(M,N) > 1.

— Input: INTEGER INCY3

On input: if INPL is .FALSE. then INCY3 is the step in index of Y between successive data elements in the third dimension of the 3D transformed data. For completely contiguous data (no gaps in Y) INCY3 should be set to L*M.
Constraint: INCY3 > 0 if INPL is .FALSE.;
INCY3 > (L-1)*INCY1+(M-1)*INCY2, if INPL is .FALSE. and N > 1.

— Input/Output: COMPLEX COMM(LCOMM)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence dimensions. The remainder is used as temporary store; if this is not sufficient for the requirements of the routine then temporary storage space will be dynamically allocated internally.

— Input: INTEGER LCOMM

On input: LCOMM is the length of the communication array COMM. The amount of internal dynamic allocation of temporary storage is dependent on the values of the increment arguments for arrays X and Y. The amount is minimized when the increments for the output array are 1, L and L*M respectively, since this represents a contiguous array in which calculations can be performedd in-place.
Constraint: If min(L,M,N) > 1, LCOMM >= 5*(L+M+N) + 150 (otherwise see see CFFT2DX).

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward 3D FFT is performed unscaled and in-place, on the leading
     C     10x10x10 submatrix of a larger 100x100x100 array of data.
     C     The result is transformed back with scaling.
     C
             SCALE = 1.0
             INPL = .TRUE.
             L = 10
             M = 10
             N = 10
             LCOMM = 2000000
             CALL CFFT3DY(0,SCALE,INPL,L,M,N,X,1,100,10000,Y,1,1,1,
          *               COMM,LCOMM,INFO)
             CALL CFFT3DY(-1,SCALE,INPL,L,M,N,X,1,100,10000,Y,1,1,1,
          *               COMM,LCOMM,INFO)
             IY = 1
             DO 20 I = 1, L
                DO 40 J = 1, M
                   DO 10 K = 1, N
                      X(I,J,K) = X(I,J,K)*EXP(-0.001*REAL(I+J+K-3))
        10      CONTINUE
        20   CONTINUE
             SCALE = 1.0/REAL(L*M*N)
             CALL CFFT3DY(1,SCALE,INPL,L,M,N,X,1,100,10000,Y,1,1,1,
          *               COMM,LCOMM,INFO)

CSFFT Routine Documentation

SUBROUTINE:CSFFT (MODE,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CSFFT.
On input:

  • MODE=0 : only initializations (specific to the values of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real transform is performed. Initializations are assumed to have been performed by a prior call to CSFFT.
  • MODE=2 : (default) initializations and a real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: INTEGER N

On input: N is the length of the sequence in X

— Input/Output: REAL X(N)

On input: X contains the Hermitian sequence of length N to be transformed.
On output: X contains the transformed real sequence.

— Input/Output: REAL COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT(0,N,X,COMM,INFO)
             CALL SCFFT(1,N,X,COMM,INFO)
             DO 10 I = N/2+2, N
                X(I) = -X(I)
        10   CONTINUE
             CALL CSFFT(2,N,X,COMM,INFO)

CSFFT1D Routine Documentation

SUBROUTINE:CSFFT1D (MODE,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CSFFT1D.
On input:

  • MODE=0 : only initializations (specific to the values of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a complex-to-real transform is performed. Initializations are assumed to have been performed by a prior call to CSFFT1D.
  • MODE=2 : (default) initializations and a complex-to-real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: INTEGER N

On input: N is length of the sequence in Y

— Input: COMPLEX X(N/2+1)

On input: X contains the complex sequence to be transformed.

— Output: REAL Y(N)

On output: Y contains the transformed real sequence of length N.

— Input/Output: REAL COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT1D(0,N,X,Y,COMM,INFO)
             CALL SCFFT1D(1,N,X,Y,COMM,INFO)
             DO 10 I = 1, N/2+1
                Y(I) = -Y(I)*EXP(-REAL(I-1)/REAL(N))
             10 CONTINUE
             CALL CSFFT1D(2,N,Y,X,COMM,INFO)

CSFFT1M Routine Documentation

SUBROUTINE:CSFFT1M (MODE,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CSFFT1M.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a complex-to-real transform is performed. Initializations are assumed to have been performed by a prior call to CSFFT1M.
  • MODE=2 : (default) initializations and a complex-to-real transform a re performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is length of the sequence in Y

— Input: COMPLEX X((N/2+1)*M)

On input: X contains the Hermitian sequences to be transformed in complex-Hermitian storage.

— Output: REAL Y(N*M)

On output: Y contains the M real sequences of length N transformed. Element i of sequence j is stored in location i + (j – 1) * N of Y.

— Input/Output: REAL COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT1M(0,M,N,X,Y,COMM,INFO)
             CALL SCFFT1M(1,M,N,X,Y,COMM,INFO)
             DO J = 1, M
                DO I = 1, N/2+1
                   Y(I,J) = -Y(I,J)*EXP(-REAL(I-1)/REAL(N))
                END DO
             EMD DO
             CALL CSFFT1M(2,M,N,Y,X,COMM,INFO)

CSFFT2D Routine Documentation

SUBROUTINE:CSFFT2D (MODE,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CSFFT2D.
On input:

  • MODE=0 : only default initializations (specific to M and N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a complex-to-real transform is performed. Initializations are assumed to have been performed by a prior call to CSFFT2D.
  • MODE=2 : (default) initializations and a complex-to-real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER M

On input: M is the number of rows in the 2D array of the real data obtained from the transform. If Y is declared as a 2D array then M is the first dimension of Y.

— Input: INTEGER N

On input: N is the number of columns in the 2D array of the real data obtained from the transform. If Y is declared as a 2D array then N is the second dimension of Y.

— Input: COMPLEX X((M/2+1)*N)

On input: X contains the Hermitian sequences in complex-Hermitian storage to be transformed.

— Output: REAL Y(M*N)

On output: Y contains the M by N real 2D array obtained from the transform. Element ij is stored in location i + (j – 1) * M of Y.

— Input/Output: REAL COMM(4*M+10*N+300)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length M*N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT2D(0,M,N,X,Y,COMM,INFO)
             CALL SCFFT2D(1,M,N,X,Y,COMM,INFO)
             DO J = 1, N
                DO I = 1, M/2+1
                   Y(I,J) = -Y(I,J)/SQRT(REAL(M*N))
                END DO
             EMD DO
             CALL CSFFT2D(2,M,N,Y,X,COMM,INFO)

CSFFT3D Routine Documentation

SUBROUTINE:CSFFT3D (MODE,L,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by CSFFT3D.
On input:

  • MODE=0 : only default initializations (specific to L, M and N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a complex-to-real transform is performed. Initializations are assumed to have been performed by a prior call to CSFFT3D.
  • MODE=2 : (default) initializations and a complex-to-real transform a re performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER L

On input: L is the length of the first dimension of the 3D array of data obtained from the transform. If Y is declared as a 3D array then L is the first dimension of Y.

— Input: INTEGER M

On input: M is the length of the second dimension of the 3D array of data obtained from the transform. If Y is declared as a 3D array then L is the second dimension of Y.

— Input: INTEGER N

On input: N is the length of the third dimension of the 3D array of data obtained from the transform. If Y is declared as a 3D array then L is the third dimension of Y.

— Input: COMPLEX X((L/2+1)*M*N)

On input: X contains the Hermitian sequences in complex-Hermitian storage to be transformed.

— Output: REAL Y(L*M*N)

On output: Y contains the L by M by N real 3D array obtained from the transform. Element ijk is stored in location i + (j – 1) * L + (k – 1) * L * M of Y.

— Input/Output: REAL COMM(4*L+10*M+10*N+500)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length L*M*N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT3D(0,L,M,N,X,Y,COMM,INFO)
             CALL SCFFT3D(1,L,M,N,X,Y,COMM,INFO)
             DO K = 1, N
                DO J = 1, M
                   DO I = 1, L/2+1
                      Y(I,J,K) = -Y(I,J,K)/SQRT(REAL(L*M*N))
                   END DO
                END DO
             EMD DO
             CALL CSFFT3D(2,L,M,N,Y,X,COMM,INFO)

CSFFTM Routine Documentation

SUBROUTINE:CSFFTM (M,N,X,COMM,INFO)

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is the length of the sequences in X

— Input/Output: REAL X(N*M)

On input: X contains the M Hermitian sequences of length N to be transformed. Element i of sequence j is stored in location i+(j-1)*N of X.

On output: X contains the transformed real sequences.

— Input/Output: REAL COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFTM(1,N,X,COMM,INFO)
             CALL SCFFTM(2,N,X,COMM,INFO)
             DO 10 I = 1, N
                X(I,3) = X(I,1)*X(N-I+1,2)
        10   CONTINUE
             CALL CSFFTM(1,N,X(1,3),COMM,INFO)

DRANDBETA / SRANDBETA

Generates a vector of random variates from a beta distribution with probability density function, f(X), where: f(X) = [Gamma(A + B) * X^(A – 1) * (1 – X)^(B – 1)] / [Gamma(A) * Gamma(B)], if 0 <= X <= 1; A,B > 0.0, otherwise f(X) = 0.

(Note that SRANDBETA is the single precision version of DRANDBETA. The argument lists of both routines are identical except that any double precision arguments of DRANDBETA are replaced in SRANDBETA by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDBETA (N,A,B,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION A

On input: first parameter for the distribution.
Constraint: A>0.

— Input: DOUBLE PRECISION B

On input: second parameter for the distribution.
Constraint: B>0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDBETA STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Beta distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION A,B
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) A,B
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Beta distribution
				CALL DRANDBETA(N,A,B,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDBINOMIAL / SRANDBINOMIAL

Generates a vector of random variates from a Binomial distribution with probability, f(X), defined by: f(X) = [M! * P^X * (1 – P)^(M – X)] / [X! * (M – 1)!], I=0,1,…,M

(Note that SRANDBINOMIAL is the single precision version of DRANDBINOMIAL. The argument lists of both routines are identical except that any double precision arguments of DRANDBINOMIAL are replaced in SRANDBINOMIAL by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDBINOMIAL (N,M,P,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER M

On input: number of trials.
Constraint: M>=0.

— Input: DOUBLE PRECISION P

On input: probability of success.
Constraint: 0<=P <1.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDBINOMIAL STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: INTEGER X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Binomial distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER M
				DOUBLE PRECISION P
				INTEGER X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M,P
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Binomial distribution
				CALL DRANDBINOMIAL(N,M,P,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDBINOMIALREFERENCE / SRANDBINOMIALREFERENCE

Initializes a reference vector for use with DRANDGENERALDISCRETE. Reference vector is for a Binomial distribution with probability, f(X), defined by: f(X) = [M! * P^X * (1 – P)^(M – X)] / [X! * (M – 1)!], I=0,1,…,M

(Note that SRANDBINOMIALREFERENCE is the single precision version of DRANDBINOMIALREFERENCE. The argument lists of both routines are identical except that any double precision arguments of DRANDBINOMIALREFERENCE are replaced in SRANDBINOMIALREFERENCE by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDBINOMIALREFERENCE (M,P,REF,LREF,INFO)

— Input: INTEGER M

On input: number of trials.
Constraint: M>=0.

— Input: DOUBLE PRECISION P

On input: probability of success.
Constraint: 0<=P <1.

— Output: DOUBLE PRECISION REF(LREF)

On output: if INFO returns with a value of 0 then REF contains reference information required to generate values from a Binomial distribution using DRANDGENERALDISCRETE.

— Input/Output: INTEGER LREF

On input: either the length of the reference vector REF, or -1.
On output: if LREF=-1 on input, then LREF is set to the recommended length of the reference vector and the routine returns. Otherwise LREF is left unchanged.

— Output: INTEGER INFO

On output: INFO is an error indicator. If INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then LREF has been set to the recommended length for the reference vector REF. If INFO = 0 then the reference vector, REF, has been successfully initialized.

Example:

		C	Generate 100 values from the Binomial distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER M
				DOUBLE PRECISION P
				INTEGER X(N)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M,P
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDBINOMIALREFERENCE(M,P,REF,LREF,INFO)
          
		C	Generate N variates from the Binomial distribution
				CALL DRANDGENERALDISCRETE(N,REF,STATE,X,INFO)
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDBLUMBLUMSHUB / SRANDBLUMBLUMSHUB

Allows direct access to the bit stream generated by the Blum-Blum-Shub generator.

(Note that SRANDBLUMBLUMSHUB is the single precision version of DRANDBLUMBLUMSHUB. The argument lists of both routines are identical except that any double precision arguments of DRANDBLUMBLUMSHUB are replaced in SRANDBLUMBLUMSHUB by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDBLUMBLUMSHUB (N,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required. The total number of bits generated is 24N.
Constraint: N>=0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDBLUMBLUMSHUB STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: INTEGER X(N)

On output: vector holding the bit stream. The least significant 24 bits of each of the X(i) contain the bit stream as generated by the Blum-Blum-Shub generator. The least significant bit of X(1) is the first bit generated, the second least significant bit of X(1) is the second bit generated etc.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.


DRANDCAUCHY / SRANDCAUCHY

Generates a vector of random variates from a Cauchy distribution with probability density function, f(X), where: f(X) = 1 / [Pi * B * (1 + (X – A / B)^2)].

(Note that SRANDCAUCHY is the single precision version of DRANDCAUCHY. The argument lists of both routines are identical except that any double precision arguments of DRANDCAUCHY are replaced in SRANDCAUCHY by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDCAUCHY (N,A,B,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION A

On input: median of the distribution.

— Input: DOUBLE PRECISION B

On input: semi-quartile range of the distribution.
Constraint: B>=0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDCAUCHY STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Cauchy distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION A,B
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) A,B
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Cauchy distribution
				CALL DRANDCAUCHY(N,A,B,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDCHISQUARED / SRANDCHISQUARED

Generates a vector of random variates from a chi-squared distribution with probability density function, f(X), where: f(X) = [X^(v / 2 – 1) * exp(-X / 2)] / [2^(v / 2) * (v / 2 – 1)!], if X > 0, otherwise f(X) = 0. Here v is the degrees of freedom, DF.

(Note that SRANDCHISQUARED is the single precision version of DRANDCHISQUARED. The argument lists of both routines are identical except that any double precision arguments of DRANDCHISQUARED are replaced in SRANDCHISQUARED by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDCHISQUARED (N,DF,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER DF

On input: degrees of freedom of the distribution.
Constraint: DF>0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDCHISQUARED STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Chi-squared distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER DF
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) DF
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Chi-squared distribution
				CALL DRANDCHISQUARED(N,DF,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDDISCRETEUNIFORM / SRANDDISCRETEUNIFORM

Generates a vector of random variates from a Uniform distribution with probability f(X) defined by: f(X) = 1 / [B – A], X=A,A+1,…,B

(Note that SRANDDISCRETEUNIFORM is the single precision version of DRANDDISCRETEUNIFORM. The argument lists of both routines are identical except that any double precision arguments of DRANDDISCRETEUNIFORM are replaced in SRANDDISCRETEUNIFORM by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDDISCRETEUNIFORM (N,A,B,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER A

On input: minimum for the distribution.

— Input: INTEGER B

On input: maximum for the distribution.
Constraint: B>=A.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDDISCRETEUNIFORM STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: INTEGER X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Uniform distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER A,B
				INTEGER X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) A,B
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Uniform distribution
				CALL DRANDDISCRETEUNIFORM(N,A,B,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDEXPONENTIAL / SRANDEXPONENTIAL

Generates a vector of random variates from an exponential distribution with probability density function, f(X), where: f(X) = exp(-X / A) / A, if X>0, otherwise f(X)=0.

(Note that SRANDEXPONENTIAL is the single precision version of DRANDEXPONENTIAL. The argument lists of both routines are identical except that any double precision arguments of DRANDEXPONENTIAL are replaced in SRANDEXPONENTIAL by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDEXPONENTIAL (N,A,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION A

On input: exponential parameter.
Constraint: A>=0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDEXPONENTIAL STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Exponential distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION A
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) A
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Exponential distribution
				CALL DRANDEXPONENTIAL(N,A,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDF / SRANDF

Generates a vector of random variates from an F distribution, also called the Fisher’s variance ratio distribution, with probability density function, f(X), where: f(X) = [((m + v – 2) / 2)! * X^(m / 2 – 1) * m^(m / 2)] / [(m / 2 – 1)! * (v / 2 – 1)! * (1 + m * X / v)^((m + v) / 2) * v^(m / 2)], if X > 0, otherwise f(X) = 0. Here m is the first degrees of freedom, (DF1) and v is the second degrees of freedom, (DF2).

(Note that SRANDF is the single precision version of DRANDF. The argument lists of both routines are identical except that any double precision arguments of DRANDF are replaced in SRANDF by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDF (N,DF1,DF2,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER DF1

On input: first degrees of freedom.
Constraint: DF1>=0.

— Input: INTEGER DF2

On input: second degrees of freedom.
Constraint: DF2>=0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDF STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the F distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER DF1,DF2
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) DF1,DF2
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the F distribution
				CALL DRANDF(N,DF1,DF2,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDGAMMA / SRANDGAMMA

Generates a vector of random variates from a Gamma distribution with probability density function, f(X), where: f(X) = [X^(A – 1) * exp(-X / B)] / [B^A * Gamma(A)], if X >= 0 and A,B > 0.0, otherwise f(X) = 0.

(Note that SRANDGAMMA is the single precision version of DRANDGAMMA. The argument lists of both routines are identical except that any double precision arguments of DRANDGAMMA are replaced in SRANDGAMMA by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDGAMMA (N,A,B,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION A

On input: first parameter of the distribution.
Constraint: A>0.

— Input: DOUBLE PRECISION B

On input: second parameter of the distribution.
Constraint: B>0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDGAMMA STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Gamma distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION A,B
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) A,B
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Gamma distribution
				CALL DRANDGAMMA(N,A,B,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDGAUSSIAN / DRANDGAUSSIAN

Generates a vector of random variates from a Gaussian distribution with probability density function, f(X), where: f(X) = [exp([-(X – m)^2] / [2s^2])] / [s * sqrt(2 * Pi)]. Here m is the mean, (XMU), and s^2 is the variance, (VAR) of the distribution.

(Note that SRANDGAUSSIAN is the single precision version of DRANDGAUSSIAN. The argument lists of both routines are identical except that any double precision arguments of DRANDGAUSSIAN are replaced in SRANDGAUSSIAN by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDGAUSSIAN (N,XMU,VAR,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION XMU

On input: mean of the distribution.

— Input: DOUBLE PRECISION VAR

On input: variance of the distribution.
Constraint: VAR>=0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDGAUSSIAN STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Gaussian distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION XMU,VAR
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) XMU,VAR
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Gaussian distribution
				CALL DRANDGAUSSIAN(N,XMU,VAR,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDGENERALDISCRETE / SRANDGENERALDISCRETE

Takes a reference vector initialized via one of DRANDBINOMIALREFERENCE, DRANDGEOMETRICREFERENCE, DRANDHYPERGEOMETRICREFERENCE, DRANDNEGATIVEBINOMIALREFERENCE, DRANDPOISSONREFERENCE and generates a vector of random variates from it.

(Note that SRANDGENERALDISCRETE is the single precision version of DRANDGENERALDISCRETE. The argument lists of both routines are identical except that any double precision arguments of DRANDGENERALDISCRETE are replaced in SRANDGENERALDISCRETE by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDGENERALDISCRETE (N,REF,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION REF(*)

On input: reference vector generated by one of the following: DRANDBINOMIALREFERENCE, DRANDGEOMETRICREFERENCE, DRANDHYPERGEOMETRICREFERENCE, DRANDNEGATIVEBINOMIALREFERENCE, DRANDPOISSONREFERENCE.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDGENERALDISCRETE STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: INTEGER X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Binomial distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER M
				DOUBLE PRECISION P
				INTEGER X(N)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M,P
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDBINOMIALREFERENCE(M,P,REF,LREF,INFO)
          
		C	Generate N variates from the Binomial distribution
				CALL DRANDGENERALDISCRETE(N,REF,STATE,X,INFO)
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDGEOMETRIC / SRANDGEOMETRIC

Generates a vector of random variates from a Geometric distribution with probability, f(X), defined by: f(X) = P * (1 – P)^X, X=0,1,….

(Note that SRANDGEOMETRIC is the single precision version of DRANDGEOMETRIC. The argument lists of both routines are identical except that any double precision arguments of DRANDGEOMETRIC are replaced in SRANDGEOMETRIC by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDGEOMETRIC (N,P,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION P

On input: distribution parameter.
Constraint: 0<=P <1.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDGEOMETRIC STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: INTEGER X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Geometric distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION P
				INTEGER X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) P
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Geometric distribution
				CALL DRANDGEOMETRIC(N,P,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDGEOMETRICREFERENCE / SRANDGEOMETRICREFERENCE

Initializes a reference vector for use with DRANDGENERALDISCRETE. Reference vector is for a Geometric distribution with probability, f(X), defined by: f(X) = P * (1 – P)^X, X=0,1,….

(Note that SRANDGEOMETRICREFERENCE is the single precision version of DRANDGEOMETRICREFERENCE. The argument lists of both routines are identical except that any double precision arguments of DRANDGEOMETRICREFERENCE are replaced in SRANDGEOMETRICREFERENCE by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDGEOMETRICREFERENCE (P,REF,LREF,INFO)

— Input: DOUBLE PRECISION P

On input: distribution parameter.
Constraint: 0<=P <1.

— Output: DOUBLE PRECISION REF(LREF)

On output: if INFO returns with a value of 0 then REF contains reference information required to generate values from a Geometric distribution using DRANDGENERALDISCRETE.

— Input/Output: INTEGER LREF

On input: either the length of the reference vector REF, or -1.
On output: if LREF=-1 on input, then LREF is set to the recommended length of the reference vector and the routine returns. Otherwise LREF is left unchanged.

— Output: INTEGER INFO

On output: INFO is an error indicator. If INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then LREF has been set to the recommended length for the reference vector REF. If INFO = 0 then the reference vector, REF, has been successfully initialized.

Example:

		C	Generate 100 values from the Geometric distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION P
				INTEGER X(N)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) P
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDGEOMETRICREFERENCE(P,REF,LREF,INFO)
          
		C	Generate N variates from the Geometric distribution
				CALL DRANDGENERALDISCRETE(N,REF,STATE,X,INFO)
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDHYPERGEOMETRIC / SRANDHYPERGEOMETRIC

Generates a vector of random variates from a Hypergeometric distribution with probability, f(X), defined by: f(X) = [s! * m! * (p – s)!] / [X! * (s – X)! * (m – X)! * (p – m – s + X)! * p!], if X = max(0,m+s-p),…,min(l,m), otherwise f(X) = 0. Here p is the size of the population, (NP), s is the size of the sample taken from the population, (NS) and m is the number of labeled, or specified, items in the population, (M).

(Note that SRANDHYPERGEOMETRIC is the single precision version of DRANDHYPERGEOMETRIC. The argument lists of both routines are identical except that any double precision arguments of DRANDHYPERGEOMETRIC are replaced in SRANDHYPERGEOMETRIC by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDHYPERGEOMETRIC (N,NP,NS,M,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER NP

On input: size of population.
Constraint: NP>=0.

— Input: INTEGER NS

On input: size of sample being taken from population.
Constraint: 0<=NS <=NP.

— Input: INTEGER M

On input: number of specified items in the population.
Constraint: 0<=M <=NP.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDHYPERGEOMETRIC STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: INTEGER X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Hypergeometric distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER NP,NS,M
				INTEGER X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) NP,NS,M
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Hypergeometric distribution
				CALL DRANDHYPERGEOMETRIC(N,NP,NS,M,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDHYPERGEOMETRICREFERENCE / SRANDHYPERGEOMETRICREFERENCE

Initializes a reference vector for use with DRANDGENERALDISCRETE. Reference vector is for a Hypergeometric distribution with probability, f(X), defined by: f(X) = [s! * m! * (p – s)!] / [X! * (s – X)! * (m – X)! * (p – m – s + X)! * p!], if X = max(0,m+s-p),…,min(l,m), otherwise f(X) = 0. Here p is the size of the population, (NP), s is the size of the sample taken from the population, (NS) and m is the number of labeled, or specified, items in the population, (M).

(Note that SRANDHYPERGEOMETRICREFERENCE is the single precision version of DRANDHYPERGEOMETRICREFERENCE. The argument lists of both routines are identical except that any double precision arguments of DRANDHYPERGEOMETRICREFERENCE are replaced in SRANDHYPERGEOMETRICREFERENCE by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDHYPERGEOMETRICREFERENCE (NP,NS,M,REF,LREF,INFO)

— Input: INTEGER NP

On input: size of population.
Constraint: NP>=0.

— Input: INTEGER NS

On input: size of sample being taken from population.
Constraint: 0<=NS <=NP.

— Input: INTEGER M

On input: number of specified items in the population.
Constraint: 0<=M <=NP.

— Output: DOUBLE PRECISION REF(LREF)

On output: if INFO returns with a value of 0 then REF contains reference information required to generate values from a Hypergeometric distribution using DRANDGENERALDISCRETE.

— Input/Output: INTEGER LREF

On input: either the length of the reference vector REF, or -1.
On output: if LREF=-1 on input, then LREF is set to the recommended length of the reference vector and the routine returns. Otherwise LREF is left unchanged.

— Output: INTEGER INFO

On output: INFO is an error indicator. If INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then LREF has been set to the recommended length for the reference vector REF. If INFO = 0 then the reference vector, REF, has been successfully initialized.

Example:

		C	Generate 100 values from the Hypergeometric distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER NP, NS,M
				INTEGER X(N)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) NP, NS,M
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDHYPERGEOMETRICREFERENCE(NP, NS,M,REF,LREF,INFO)
          
		C	Generate N variates from the Hypergeometric distribution
				CALL DRANDGENERALDISCRETE(N,REF,STATE,X,INFO)
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDINITIALIZE / SRANDINITIALIZE

Initialize one of the five supplied base generators; NAG basic generator, Wichmann-Hill generator, Mersenne Twister, L’Ecuyer’s combined recursive generator (MRG32k3a) or the Blum-Blum-Shub generator.

(Note that SRANDINITIALIZE is the single precision version of DRANDINITIALIZE. The argument lists of both routines are identical except that any double precision arguments of DRANDINITIALIZE are replaced in SRANDINITIALIZE by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDINITIALIZE (GENID,SUBID,SEED,LSEED,STATE, LSTATE,INFO)

— Input: INTEGER GENID

On input: a numerical code indicating which of the five base generators to initialize.

Constraint: 1<=GENID <=5.

— Input: INTEGER SUBID

On input: if GENID=2, then SUBID indicates which of the 273 Wichmann-Hill generators to use. If GENID=5 then SUBID indicates the number of bits to use (v) from each of iteration of the Blum-Blum-Shub generator. In all other cases SUBID is not referenced.
Constraint: If GENID=2 then 1<=SUBID <=273 .

— Input: INTEGER SEED(LSEED)

On input: if GENID is not equal to 5 , then SEED is a vector of initial values for the base generator. These values must be positive integers. The number of values required depends on the base generator being used. The NAG basic generator requires one initial value, the Wichmann-Hill generator requires four initial values, the L’Ecuyer combined recursive generator requires six initial values and the Mersenne Twister requires 624 initial values. If the number of seeds required by the chosen generator is >LSEED then SEED(1) is used to initialize the NAG basic generator. This is then used to generate all of the remaining seed values required. In general it is best not to set all the elements of SEED to anything too obvious, such as a single repeated value or a simple sequence. Using such a seed array may lead to several similar values being created in a row when the generator is subsequently called. This is particularly true for the Mersenne Twister generator.

In order to initialize the Blum-Blum-Shub generator two large prime values, p and q are required as well as an initial value s. As p, q and s can be of an arbitrary size, these values are factored into a expressed as a polynomial in B, where B=2^(24). For example, p can be factored into a polynomial of order l_p, with p = p_1 + p_2B + p_3B^2 + … + p_(l_p) * B^(l_p-1). The elements of SEED should then be set to the following:

  • SEED(1) = l_p
  • SEED(2) to SEED(l_p+1) = p_1 to p_(l_p)
  • SEED(l_p+2) = l_q
  • SEED(l_p+3) to SEED(l_p+l_q+2) = q_1 to q_(l_q)
  • SEED(l_p+l_q+3) = l_s
  • SEED(l_p+l_q+4) to SEED(l_p+l_q+l_s+3) = s_1 to s_(l_s)

Constraint: If GENID is not equal to 5 then SEED(i)>0 for i=1,2,…. If GENID=5 then SEED must take the values described above.

— Input/Output: INTEGER LSEED

On input: either the length of the seed vector, SEED, or a value <=0 .
On output: if LSEED<=0 on input, then LSEED is set to the number of initial values required by the selected generator, and the routine returns. Otherwise LSEED is left unchanged.

— Output: INTEGER STATE(LSTATE)

On output: the state vector required by all of the supplied distributional and base generators.

— Input/Output: INTEGER LSTATE

On input: either the length of the state vector, STATE, or a value <=0 .
On output: if LSTATE<=0 on input, then LSTATE is set to the minimum length of the state vector STATE for the base generator chosen, and the routine returns. Otherwise LSTATE is left unchanged.
Constraint: LSTATE<=0 or the minimum length for the chosen base generator, given by:

  • GENID=1: LSTATE>=16,
  • GENID=2: LSTATE>=20,
  • GENID=3: LSTATE>=633,
  • GENID=4: LSTATE>=61,
  • GENID=5: LSTATE>=l_p+l_q+l_s+6 , where l_p, l_q and l_s are the order of the polynomials used to express the parameters p,q and s respectively.

— Output: INTEGER INFO

On output: INFO is an error indicator. If INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then either, or both of LSEED and / or LSTATE have been set to the required length for vectors SEED and STATE respectively. Of the two variables LSEED and LSTATE, only those which had an input value <=0 will have been set. The STATE vector will not have been initialized. If INFO = 0 then the state vector, STATE, has been successfully initialized.

Example:

		C	Generate 100 values from the Beta distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION A,B
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) A,B
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Beta distribution
				CALL DRANDBETA(N,A,B,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)


DRANDINITIALIZEBBS / SRANDINITIALIZEBBS

Alternative initialization routine for the Blum-Blum-Shub generator. Unlike the other base generators supplied with the ACML, the Blum-Blum-Shub generator requires two additional parameters, p and q as well as an initial state, s. As p, q and s can be of an arbitrary size. In order to avoid overflow these values are expressed as a polynomial in B, where B=2^(24). For example, p can be factored into a polynomial of order l_p, with p = p_1 + p_2B + p_3B^2 + … + p_(l_p) * B^(l_p-1), similarly q = q_1 + q_2B + q_3B^2 + … + q_(l_q) * B^(l_q-1) and s = s_1 + s_2B + s_3B^2 + \cdots + s_(l_s) * B^(l_s-1).

(Note that SRANDINITIALIZEBBS is the single precision version of DRANDINITIALIZEBBS. The argument lists of both routines are identical except that any double precision arguments of DRANDINITIALIZEBBS are replaced in SRANDINITIALIZEBBS by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDINITIALIZEBBS (NBITS,LP,P,LQ,Q,LS,S,STATE,LSTATE, INFO)

— Input: INTEGER NBITS

On input: the number of bits, v, to use from each iteration of the Blum-Blum-Shub generator. If NBITS<1 then NBITS=1. If NBITS>15 then NBITS=15.

— Input: INTEGER LP

On input: the order of the polynomial used to express p (l_p).
Constraint: 1 <=LP <=25.

— Input: INTEGER P(LP)

On input: the coefficients of the polynomial used to express p. P(i)=p_i, i=1 to l_p.
Constraint: 0 <=P (i) < 2^(24)

— Input: INTEGER LQ

On input: the order of the polynomial used to express q (l_q).
Constraint: 1 <=LQ <=25.

— Input: INTEGER Q(LQ)

On input: the coefficients of the polynomial used to express q. Q(i)=q_i, i=1 to l_q.
Constraint: 0 <=Q (i) < 2^(24)

— Input: INTEGER LS

On input: the order of the polynomial used to express s (l_s).
Constraint: 1 <=LS <=25.

— Input: INTEGER S(LS)

On input: the coefficients of the polynomial used to express s. S(i)=s_i, i=1 to l_s.
Constraint: 0 <=S (i) < 2^(24)

— Output: INTEGER STATE(*)

On output: the initial state for the Blum-Blum-Shub generator with parameters P,Q,S and NBITS.

— Input/Output: INTEGER LSTATE

On input: either the length of the state vector, STATE, or a value <=0 .
On output: if LSTATE<=0 on input, then LSTATE is set to the minimum length of the state vector STATE for the parameters chosen, and the routine returns. Otherwise LSTATE is left unchanged.
Constraint: LSTATE<=0 or LSTATE >=l_p+l_q+l_s+6

— Output: INTEGER INFO

On output: INFO is an error indicator. If INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then LSTATE has been set to the required length for the STATE vector. If INFO = 0 then the state vector, STATE, has been successfully initialized.


DRANDINITIALIZEUSER / SRANDINITIALIZEUSER

Registers a user supplied base generator so that it can be used with the ACML distributional generators.

(Note that SRANDINITIALIZEUSER is the single precision version of DRANDINITIALIZEUSER. The argument lists of both routines are identical except that any double precision arguments of DRANDINITIALIZEUSER are replaced in SRANDINITIALIZEUSER by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDINITIALIZEUSER (UINI,UGEN,GENID,SUBID,SEED,LSEED, STATE,LSTATE,INFO)

— Input: SUBROUTINE UINI

On input: routine that will be used to initialize the user supplied generator, UGEN.

— Input: SUBROUTINE UGEN

On input: user supplied base generator.

— Input: INTEGER GENID

On input: parameter is passed directly to UINI. Its function therefore depends on that routine.

— Input: INTEGER SUBID

On input: parameter is passed directly to UINI. Its function therefore depends on that routine.

— Input: INTEGER SEED(LSEED)

On input: parameter is passed directly to UINI. Its function therefore depends on that routine.

— Input/Output: INTEGER LSEED

On input: length of the vector SEED. This parameter is passed directly to UINI and therefore its required value depends on that routine.
On output: whether LSEED changes will depend on UINI.

— Output: INTEGER STATE(LSTATE)

On output: the state vector required by all of the supplied distributional generators. The value of STATE returned by UINI has some housekeeping elements appended to the end before being returned by DRANDINITIALIZEUSER. See User Supplied Generators for details about the form of STATE.

— Input/Output: INTEGER LSTATE

On input: length of the vector STATE. This parameter is passed directly to UINI and therefore its required value depends on that routine.
On output: whether LSTATE changes will depend on UINI. If LSTATE<=0 then it is assumed that a request for the required length of STATE has been made. The value of LSTATE returned from UINI is therefore adjusted to allow for housekeeping elements to be added to the end of the STATE vector. This results in the value of LSTATE returned by DRANDINITIALIZEUSER being 3 larger than that returned by UINI.

— Output: INTEGER INFO

On output: INFO is an error indicator. DRANDINITIALIZEUSER will return a value of -6 if the value of LSTATE is between 1 and 3. Otherwise INFO is passed directly back from UINI. It is recommended that the value of INFO returned by UINI is kept consistent with the rest of the ACML, that is if INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then either, or both of LSEED and / or LSTATE have been set to the required length for vectors SEED and STATE respectively and the STATE vector has not have been initialized. If INFO = 0 then the state vector, STATE, has been successfully initialized.

Example:

          C     Generate 100 values from the Uniform distribution using
          C     a user supplied base generator
                INTEGER LSTATE,N
                PARAMETER (LSTATE=16,N=100)
                INTEGER I,INFO,NSKIP,SEED(1),STATE(LSTATE)
                INTEGER X(N)
                DOUBLE PRECISION A,B
          
          C     Set the seed
                SEED(1) = 1234
          
          C     Set the distributional parameters
                A = 0.0D0
                B = 1.0D0
          
          C     Initialize the base generator. Here ACMLRNGNB0GND is a user
          C     supplied generator and ACMLRNGNB0INI its initializer
                CALL DRANDINITIALIZEUSER(ACMLRNGNB0INI,ACMLRNGNB0GND,1,0,SEED,
               *    LSEED,STATE,LSTATE,INFO)
          
          C     Generate N variates from the Univariate distribution
                CALL DRANDUNIFORM(N,A,B,STATE,X,LDX,INFO)
          
          C     Print the results
                WRITE(6,*) (X(I),I=1,N)

DRANDLEAPFROG / SRANDLEAPFROG

Amend a generator so that it will generate every Kth value.

(Note that SRANDLEAPFROG is the single precision version of DRANDLEAPFROG. The argument lists of both routines are identical except that any double precision arguments of DRANDLEAPFROG are replaced in SRANDLEAPFROG by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDLEAPFROG (N,K,STATE,INFO)

— Input: INTEGER N

On input: total number of streams being used.
Constraint: N>0.

— Input: INTEGER K

On input: number of the current stream
Constraint: 0<K <=N.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDLEAPFROG STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.


On output: The STATE vector for a generator that has been advanced K-1 places and will return every Nth value.
Constraint: The STATE array must be for either the NAG basic, Wichmann-Hill or L’Ecuyer Combined Recursive base generators.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Generate 3 * 100 values from the Uniform distribution
     C     Multiple streams generated using the Leap Frog method
           INTEGER LSTATE,N
           PARAMETER (LSTATE=16,N=100)
           INTEGER I,INFO
           INTEGER SEED(1),STATE1(LSTATE),STATE2(LSTATE),STATE3(LSTATE)
           INTEGER  X1(N),X2(N),X3(N)
           DOUBLE PRECISION A,B
     
     C     Set the seed
           SEED(1) = 1234
     
     C     Set the distributional parameters
           A = 0.0D0
           B = 1.0D0
     
     C     Initialize the STATE1 vector
           CALL DRANDINITIALIZE(1,1,SEED,1,STATE1,LSTATE,INFO)
     
     C     Copy the STATE1 vector into other state vectors
           DO 20 I = 1,LSTATE
             STATE2(I) = STATE1(I)
             STATE3(I) = STATE1(I)
     20    CONTINUE
     
     C     Update each stream so they generate every 3rd value
           CALL DRANDLEAPFROG(3,1,STATE1,INFO)
           CALL DRANDLEAPFROG(3,2,STATE2,INFO)
           CALL DRANDLEAPFROG(3,3,STATE3,INFO)
     
     C     Generate 3 sets of N variates from the Univariate distribution
           CALL DRANDUNIFORM(N,A,B,STATE1,X1,LDX,INFO)
           CALL DRANDUNIFORM(N,A,B,STATE2,X2,LDX,INFO)
           CALL DRANDUNIFORM(N,A,B,STATE3,X3,LDX,INFO)
     
     C     Print the results
           DO 40 I = 1,N
             WRITE(6,*) X1(I),X2(I),X3(I)
     40    CONTINUE

DRANDLOGISTIC / SRANDLOGISTIC

Generates a vector of random variates from a logistic distribution with probability density function, f(X), where: f(X) = [exp((X – A) / B)] / [B * (1 + exp((X – A) / B))^2].

(Note that SRANDLOGISTIC is the single precision version of DRANDLOGISTIC. The argument lists of both routines are identical except that any double precision arguments of DRANDLOGISTIC are replaced in SRANDLOGISTIC by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDLOGISTIC (N,A,B,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION A

On input: mean of the distribution.

— Input: DOUBLE PRECISION B

On input: spread of the distribution. B = SQRT(3) * SIGMA / PI, where SIGMA is the standard deviation of the distribution.
Constraint: B>0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDLOGISTIC STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Logistic distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION A,B
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) A,B
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Logistic distribution
				CALL DRANDLOGISTIC(N,A,B,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDLOGNORMAL / SRANDLOGNORMAL

Generates a vector of random variates from a lognormal distribution with probability density function, f(X), where: f(X) = [exp([-(log(X) – m)^2] / [2s^2])] / [X * s * sqrt(2 * Pi)], if X > 0, otherwise f(X) = 0. Here m is the mean, (XMU), and s^2 is the variance, (VAR) of the underlying Gaussian distribution.

(Note that SRANDLOGNORMAL is the single precision version of DRANDLOGNORMAL. The argument lists of both routines are identical except that any double precision arguments of DRANDLOGNORMAL are replaced in SRANDLOGNORMAL by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDLOGNORMAL (N,XMU,VAR,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION XMU

On input: mean of the underlying Gaussian distribution.

— Input: DOUBLE PRECISION VAR

On input: variance of the underlying Gaussian distribution.
Constraint: VAR>=0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDLOGNORMAL STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Lognormal distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION XMU,VAR
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) XMU,VAR
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Lognormal distribution
				CALL DRANDLOGNORMAL(N,XMU,VAR,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDMULTINOMIAL / SRANDMULTINOMIAL

Generates a matrix of random variates from a Multinomial distribution with probability, f(X), defined by: f(X) = [M! * P1^(X1) * P2^(X2) * … * PK^(XK)] / [X1! * X2! * … * XK!], where X = (X1,X2,…,XK), P = (P1,P2,…,PK), X1 + X2 + … + XK = 1 and P1 + P2 + … + PK = 1.

(Note that SRANDMULTINOMIAL is the single precision version of DRANDMULTINOMIAL. The argument lists of both routines are identical except that any double precision arguments of DRANDMULTINOMIAL are replaced in SRANDMULTINOMIAL by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDMULTINOMIAL (N,M,P,K,STATE,X,LDX,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER M

On input: number of trials.
Constraint: M>=0.

— Input: DOUBLE PRECISION P(K)

On input: vector of probabilities for each of the K possible outcomes.
Constraint: 0 <= P <= 1, for all elements of P and P1 + P2 + … + PK = 1.

— Input: INTEGER K

On input: number of possible outcomes.
Constraint: K>=2.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDBINOMIAL STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: INTEGER X(LDX,K)

On output: matrix of variates from the specified distribution.

— Input: INTEGER LDX

On input: leading dimension of X in the calling routine.
Constraint: LDX>=N.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C Generate 100 values from the Multinomial distribution
           INTEGER LSTATE,N, MK
           PARAMETER (LSTATE=16,N=100,MK=10)
           INTEGER I,J,INFO,SEED(1),STATE(LSTATE)
           INTEGER LDC,LDX,K,M
           INTEGER  X(N,MK)
           DOUBLE PRECISION P(MK)
     
     C Set array sizes
           LDX = N
     
     C Set the seed
           SEED(1) = 1234
     
     C Read in the distributional parameters
           READ(5,*) K
           READ(5,*) (P(I),I=1,K)
     
     C Initialize the STATE vector
           CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
     
     C Generate N variates from the Multinomial distribution
           CALL DRANDMULTINOMIAL(N,M,P,K,STATE,X,LDX,INFO)
     
     C Print the results
           DO 20 I = 1,N
             WRITE(6,*) (X(I,J),J=1,K)
           20 CONTINUE

DRANDMULTINORMAL / SRANDMULTINORMAL

Generates an array of random variates from a Multivariate Normal distribution with probability density function, f(X), where: f(X) = sqrt([det(C^(-1))] / [(2 * Pi)^M]) * exp(-(X – u)^T * C^(-1) * (X – u)), where u is the vector of means, XMU.

(Note that SRANDMULTINORMAL is the single precision version of DRANDMULTINORMAL. The argument lists of both routines are identical except that any double precision arguments of DRANDMULTINORMAL are replaced in SRANDMULTINORMAL by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDMULTINORMAL (N,M,XMU,C,LDC,STATE,X,LDX,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER M

On input: number of dimensions for the distribution.
Constraint: M>=1.

— Input: DOUBLE PRECISION XMU(M)

On input: vector of means for the distribution.

— Input: DOUBLE PRECISION C(LDC,M)

On input: variance / covariance matrix for the distribution.

— Input: INTEGER LDC

On input: leading dimension of C in the calling routine.
Constraint: LDC>=M.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDMULTINORMAL STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(LDX,M)

On output: matrix of variates from the specified distribution.

— Input: INTEGER LDX

On input: leading dimension of X in the calling routine.
Constraint: LDX>=N.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the
		C	Multivariate Normal distribution
				INTEGER LSTATE,N, MM
				PARAMETER (LSTATE=16,N=100,MM=10)
				INTEGER I,J,INFO,SEED(1),STATE(LSTATE)
				INTEGER LDC,LDX,M
				DOUBLE PRECISION X(N,MM),XMU(MM),C(MM,MM)
		C	Set array sizes
				LDC = MM
				LDX = N
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M
				READ(5,*) (XMU(I),I=1,M)
				DO 20 I = 1,M
			READ(5,*) (C(I,J),J=1,M)
		20	CONTINUE
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the
		C	Multivariate Normal distribution
				CALL DRANDMULTINORMAL(N,M,XMU,C,LDC,STATE,X,LDX,INFO)
          
		C	Print the results
				DO 40 I = 1,N
			WRITE(6,*) (X(I,J),J=1,M)
		40	CONTINUE

DRANDMULTINORMALR / SRANDMULTINORMALR

Generates an array of random variates from a Multivariate Normal distribution using a reference vector initialized by DRANDMULTINORMALREFERENCE.

(Note that SRANDMULTINORMALR is the single precision version of DRANDMULTINORMALR. The argument lists of both routines are identical except that any double precision arguments of DRANDMULTINORMALR are replaced in SRANDMULTINORMALR by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDMULTINORMALR (N,REF,STATE,X,LDX,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION REF(*)

On input: a reference vector generated by DRANDMULTINORMALREFERENCE.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDMULTINORMALR STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(LDX,M)

On output: matrix of variates from the specified distribution.

— Input: INTEGER LDX

On input: leading dimension of X in the calling routine.
Constraint: LDX>=N.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the
		C	Multivariate Normal distribution
				INTEGER LSTATE,N, MM
				PARAMETER (LSTATE=16,N=100,MM=10)
				INTEGER I,J,INFO,SEED(1),STATE(LSTATE)
				INTEGER LDC,LDX,M
				DOUBLE PRECISION X(N,MM),XMU(MM),C(MM,MM)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set array sizes
				LDC = MM
				LDX = N
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M
				READ(5,*) (XMU(I),I=1,M)
				DO 20 I = 1,M
				  READ(5,*) (C(I,J),J=1,M)
          20     CONTINUE
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDMULTINORMALREFERENCE(M,XMU,C,LDC,REF,LREF,INFO)
          
		C	Generate N variates from the
		C	Multivariate Normal distribution
				CALL DRANDMULTINORMALR(N,REF,STATE,X,LDX,INFO)
		C	Print the results
				DO 40 I = 1,N
				  WRITE(6,*) (X(I,J),J=1,M)
          40     CONTINUE

DRANDMULTINORMALREFERENCE / SRANDMULTINORMALREFERENCE

Initializes a reference vector for use with DRANDMULTINORMALR. Reference vector is for a Multivariate Normal distribution with probability density function, f(X), where: f(X) = sqrt([det(C^(-1))] / [(2 * Pi)^M]) * exp(-(X – u)^T * C^(-1) * (X – u)), where u is the vector of means, XMU.

(Note that SRANDMULTINORMALREFERENCE is the single precision version of DRANDMULTINORMALREFERENCE. The argument lists of both routines are identical except that any double precision arguments of DRANDMULTINORMALREFERENCE are replaced in SRANDMULTINORMALREFERENCE by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDMULTINORMALREFERENCE (M,XMU,C,LDC,REF,LREF,INFO)

— Input: INTEGER M

On input: number of dimensions for the distribution.
Constraint: M>=1.

— Input: DOUBLE PRECISION XMU(M)

On input: vector of means for the distribution.

— Input: DOUBLE PRECISION C(LDC,M)

On input: variance / covariance matrix for the distribution.

— Input: INTEGER LDC

On input: leading dimension of C in the calling routine.
Constraint: LDC>=M.

— Output: DOUBLE PRECISION REF(LREF)

On output: if INFO returns with a value of 0 then REF contains reference information required to generate values from a Multivariate Normal distribution using DRANDMULTINORMALR.

— Input/Output: INTEGER LREF

On input: either the length of the reference vector REF, or -1.
On output: if LREF=-1 on input, then LREF is set to the recommended length of the reference vector and the routine returns. Otherwise LREF is left unchanged.

— Output: INTEGER INFO

On output: INFO is an error indicator. If INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then LREF has been set to the recommended length for the reference vector REF. If INFO = 0 then the reference vector, REF, has been successfully initialized.

Example:

		C	Generate 100 values from the
		C	Multivariate Normal distribution
				INTEGER LSTATE,N, MM
				PARAMETER (LSTATE=16,N=100,MM=10)
				INTEGER I,J,INFO,SEED(1),STATE(LSTATE)
				INTEGER LDC,LDX,M
				DOUBLE PRECISION X(N,MM),XMU(MM),C(MM,MM)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set array sizes
				LDC = MM
				LDX = N
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M
				READ(5,*) (XMU(I),I=1,M)
				DO 20 I = 1,M
				  READ(5,*) (C(I,J),J=1,M)
          20     CONTINUE
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDMULTINORMALREFERENCE(M,XMU,C,LDC,REF,LREF,INFO)
          
		C	Generate N variates from the
		C	Multivariate Normal distribution
				CALL DRANDMULTINORMALR(N,REF,STATE,X,LDX,INFO)
		C	Print the results
				DO 40 I = 1,N
				  WRITE(6,*) (X(I,J),J=1,M)
          40     CONTINUE

DRANDMULTISTUDENTST / SRANDMULTISTUDENTST

Generates an array of random variates from a Multivariate Students T distribution with probability density function, f(X), where: f(X) = [Gamma([v + M] / 2)] / [(Pi * v)^(M / 2) * Gamma(v / 2) * det(C)^(1 / 2)] * (1 + [(X – u)^T * C^(-1) * (X – u)] / v)^(-[v + M] / 2), where u us the vector of means, XMU and v is the degrees of freedom, DF.

(Note that SRANDMULTISTUDENTST is the single precision version of DRANDMULTISTUDENTST. The argument lists of both routines are identical except that any double precision arguments of DRANDMULTISTUDENTST are replaced in SRANDMULTISTUDENTST by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDMULTISTUDENTST (N,M,DF,XMU,C,LDC,STATE,X,LDX,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER M

On input: number of dimensions for the distribution.
Constraint: M>=1.

— Input: INTEGER DF

On input: degrees of freedom.
Constraint: DF>2.

— Input: DOUBLE PRECISION XMU(M)

On input: vector of means for the distribution.

— Input: DOUBLE PRECISION C(LDC,M)

On input: matrix defining the variance / covariance for the distribution. The variance / covariance matrix is given by [DF * C]/ [DF – 2].

— Input: INTEGER LDC

On input: leading dimension of C in the calling routine.
Constraint: LDC>=M.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDMULTISTUDENTST STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(LDX,M)

On output: matrix of variates from the specified distribution.

— Input: INTEGER LDX

On input: leading dimension of X in the calling routine.
Constraint: LDX>=N.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the
		C	Multivariate Students T distribution
				INTEGER LSTATE,N, MM
				PARAMETER (LSTATE=16,N=100,MM=10)
				INTEGER I,J,INFO,SEED(1),STATE(LSTATE)
				INTEGER LDC,LDX,M,DF
				DOUBLE PRECISION X(N,MM),XMU(MM),C(MM,MM)
		C	Set array sizes
				LDC = MM
				LDX = N
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M,DF
				READ(5,*) (XMU(I),I=1,M)
				DO 20 I = 1,M
				  READ(5,*) (C(I,J),J=1,M)
          20     CONTINUE
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the
		C	Multivariate Students T distribution
				CALL DRANDMULTISTUDENTST(N,M,DF,XMU,C,LDC,STATE,X,LDX,INFO)
          
		C	Print the results
				DO 40 I = 1,N
				  WRITE(6,*) (X(I,J),J=1,M)
          40     CONTINUE

DRANDMULTISTUDENTSTR / SRANDMULTISTUDENTSTR

Generates an array of random variates from a Multivariate Students T distribution using a reference vector initialized by DRANDMULTISTUDENTSTREFERENCE.

(Note that SRANDMULTISTUDENTSTR is the single precision version of DRANDMULTISTUDENTSTR. The argument lists of both routines are identical except that any double precision arguments of DRANDMULTISTUDENTSTR are replaced in SRANDMULTISTUDENTSTR by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDMULTISTUDENTSTR (N,REF,STATE,X,LDX,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION REF(*)

On input: a reference vector generated by DRANDMULTISTUDENTSTREFERENCE.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDMULTISTUDENTSTR STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(LDX,M)

On output: matrix of variates from the specified distribution.

— Input: INTEGER LDX

On input: leading dimension of X in the calling routine.
Constraint: LDX>=N.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the
		C	Multivariate Students T distribution
				INTEGER LSTATE,N, MM
				PARAMETER (LSTATE=16,N=100,MM=10)
				INTEGER I,J,INFO,SEED(1),STATE(LSTATE)
				INTEGER LDC,LDX,M,DF
				DOUBLE PRECISION X(N,MM),XMU(MM),C(MM,MM)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set array sizes
				LDC = MM
				LDX = N
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M,DF
				READ(5,*) (XMU(I),I=1,M)
				DO 20 I = 1,M
				  READ(5,*) (C(I,J),J=1,M)
          20     CONTINUE
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDMULTISTUDENTSTREFERENCE(M,DF,XMU,C,LDC,REF,LREF,INFO)
          
		C	Generate N variates from the
		C	Multivariate Students T distribution
				CALL DRANDMULTISTUDENTSTR(N,REF,STATE,X,LDX,INFO)
		C	Print the results
				DO 40 I = 1,N
				  WRITE(6,*) (X(I,J),J=1,M)
          40     CONTINUE

DRANDMULTISTUDENTSTREFERENCE / SRANDMULTISTUDENTSTREFERENCE

Initializes a reference vector for use with DRANDMULTISTUDENTSTR. Reference vector is for a Multivariate Students T distribution with probability density function, f(X), where: f(X) = [Gamma([v + M] / 2)] / [(Pi * v)^(M / 2) * Gamma(v / 2) * det(C)^(1 / 2)] * (1 + [(X – u)^T * C^(-1) * (X – u)] / v)^(-[v + M] / 2), where u us the vector of means, XMU and v is the degrees of freedom, DF.

(Note that SRANDMULTISTUDENTSTREFERENCE is the single precision version of DRANDMULTISTUDENTSTREFERENCE. The argument lists of both routines are identical except that any double precision arguments of DRANDMULTISTUDENTSTREFERENCE are replaced in SRANDMULTISTUDENTSTREFERENCE by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDMULTISTUDENTSREFERENCE (M,DF,XMU,C,LDC,REF,LREF,INFO)

— Input: INTEGER M

On input: number of dimensions for the distribution.
Constraint: M>=1.

— Input: INTEGER DF

On input: degrees of freedom.
Constraint: DF>2.

— Input: DOUBLE PRECISION XMU(M)

On input: vector of means for the distribution.

— Input: DOUBLE PRECISION C(LDC,M)

On input: matrix defining the variance / covariance for the distribution. The variance / covariance matrix is given by [DF * C]/ [DF – 2].

— Input: INTEGER LDC

On input: leading dimension of C in the calling routine.
Constraint: LDC>=M.

— Output: DOUBLE PRECISION REF(LREF)

On output: if INFO returns with a value of 0 then REF contains reference information required to generate values from a Multivariate Students T distribution using DRANDMULTISTUDENTSTR.

— Input/Output: INTEGER LREF

On input: either the length of the reference vector REF, or -1.
On output: if LREF=-1 on input, then LREF is set to the recommended length of the reference vector and the routine returns. Otherwise LREF is left unchanged.

— Output: INTEGER INFO

On output: INFO is an error indicator. If INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then LREF has been set to the recommended length for the reference vector REF. If INFO = 0 then the reference vector, REF, has been successfully initialized.

Example:

		C	Generate 100 values from the
		C	Multivariate Students T distribution
				INTEGER LSTATE,N, MM
				PARAMETER (LSTATE=16,N=100,MM=10)
				INTEGER I,J,INFO,SEED(1),STATE(LSTATE)
				INTEGER LDC,LDX,M,DF
				DOUBLE PRECISION X(N,MM),XMU(MM),C(MM,MM)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set array sizes
				LDC = MM
				LDX = N
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M,DF
				READ(5,*) (XMU(I),I=1,M)
				DO 20 I = 1,M
				  READ(5,*) (C(I,J),J=1,M)
          20     CONTINUE
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDMULTISTUDENTSTREFERENCE(M,DF,XMU,C,LDC,REF,LREF,INFO)
          
		C	Generate N variates from the
		C	Multivariate Students T distribution
				CALL DRANDMULTISTUDENTSTR(N,REF,STATE,X,LDX,INFO)
		C	Print the results
				DO 40 I = 1,N
				  WRITE(6,*) (X(I,J),J=1,M)
          40     CONTINUE


DRANDNEGATIVEBINOMIAL / SRANDNEGATIVEBINOMIAL

Generates a vector of random variates from a Negative Binomial distribution with probability f(X) defined by: f(X) = [(M + X – 1)! * P^X * (1 – P)^M] / [X! * (M – 1)!], X=0,1,…

(Note that SRANDNEGATIVEBINOMIAL is the single precision version of DRANDNEGATIVEBINOMIAL. The argument lists of both routines are identical except that any double precision arguments of DRANDNEGATIVEBINOMIAL are replaced in SRANDNEGATIVEBINOMIAL by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDNEGATIVEBINOMIAL (N,M,P,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER M

On input: number of failures.
Constraint: M>=0.

— Input: DOUBLE PRECISION P

On input: probability of success.
Constraint: 0<=P <1.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDNEGATIVEBINOMIAL STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: INTEGER X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Negative Binomial distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER M
				DOUBLE PRECISION P
				INTEGER X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M,P
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Negative Binomial distribution
				CALL DRANDNEGATIVEBINOMIAL(N,M,P,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDNEGATIVEBINOMIALREFERENCE / SRANDNEGATIVEBINOMIALREFERENCE

Initializes a reference vector for use with DRANDGENERALDISCRETE. Reference vector is for a Negative Binomial distribution with probability f(X) defined by: f(X) = [(M + X – 1)! * P^X * (1 – P)^M] / [X! * (M – 1)!], X=0,1,…

(Note that SRANDNEGATIVEBINOMIALREFERENCE is the single precision version of DRANDNEGATIVEBINOMIALREFERENCE. The argument lists of both routines are identical except that any double precision arguments of DRANDNEGATIVEBINOMIALREFERENCE are replaced in SRANDNEGATIVEBINOMIALREFERENCE by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDNEGATIVEBINOMIALREFERENCE (M,P,REF,LREF,INFO)

— Input: INTEGER M

On input: number of failures.
Constraint: M>=0.

— Input: DOUBLE PRECISION P

>

On input: probability of success.
Constraint: 0<=P <1.

— Output: DOUBLE PRECISION REF(LREF)

On output: if INFO returns with a value of 0 then REF contains reference information required to generate values from a Negative Binomial distribution using DRANDGENERALDISCRETE.

— Input/Output: INTEGER LREF

On input: either the length of the reference vector REF, or -1.
On output: if LREF=-1 on input, then LREF is set to the recommended length of the reference vector and the routine returns. Otherwise LREF is left unchanged.

— Output: INTEGER INFO

On output: INFO is an error indicator. If INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then LREF has been set to the recommended length for the reference vector REF. If INFO = 0 then the reference vector, REF, has been successfully initialized.

Example:

		C	Generate 100 values from the Negative Binomial distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER M
				DOUBLE PRECISION P
				INTEGER X(N)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) M,P
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDNEGATIVEBINOMIALREFERENCE(M,P,REF,LREF,INFO)
          
		C	Generate N variates from the Negative Binomial distribution
				CALL DRANDGENERALDISCRETE(N,REF,STATE,X,INFO)
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDPOISSON / SRANDPOISSON

Generates a vector of random variates from a Poisson distribution with probability f(X) defined by: f(X) = [l^X * exp(-l)] / X!, X=0,1,…, where l is the mean of the distribution, LAMBDA.

(Note that SRANDPOISSON is the single precision version of DRANDPOISSON. The argument lists of both routines are identical except that any double precision arguments of DRANDPOISSON are replaced in SRANDPOISSON by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDPOISSON (N,LAMBDA,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER M

On input: number of failures.
Constraint: M>=0.

— Input: DOUBLE PRECISION LAMBDA

On input: mean of the distribution.
Constraint: LAMBDA>=0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDPOISSON STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: INTEGER X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Poisson distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION LAMBDA
				INTEGER X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) LAMBDA
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Poisson distribution
				CALL DRANDPOISSON(N,LAMBDA,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDPOISSONREFERENCE / SRANDPOISSONREFERENCE

Initializes a reference vector for use with DRANDGENERALDISCRETE. Reference vector is for a Poisson distribution with probability f(X) defined by: f(X) = [l^X * exp(-l)] / X!, X=0,1,…, where l is the mean of the distribution, LAMBDA.

(Note that SRANDPOISSONREFERENCE is the single precision version of DRANDPOISSONREFERENCE. The argument lists of both routines are identical except that any double precision arguments of DRANDPOISSONREFERENCE are replaced in SRANDPOISSONREFERENCE by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDPOISSONREFERENCE (LAMBDA,REF,LREF,INFO)

— Input: INTEGER M

On input: number of failures.
Constraint: M>=0.

— Input: DOUBLE PRECISION LAMBDA

On input: mean of the distribution.
Constraint: LAMBDA>=0.

— Output: DOUBLE PRECISION REF(LREF)

On output: if INFO returns with a value of 0 then REF contains reference information required to generate values from a Poisson distribution using DRANDGENERALDISCRETE.

— Input/Output: INTEGER LREF

On input: either the length of the reference vector REF, or -1.
On output: if LREF=-1 on input, then LREF is set to the recommended length of the reference vector and the routine returns. Otherwise LREF is left unchanged.

— Output: INTEGER INFO

On output: INFO is an error indicator. If INFO = -i on exit, the i-th argument had an illegal value. If INFO =1 on exit, then LREF has been set to the recommended length for the reference vector REF. If INFO = 0 then the reference vector, REF, has been successfully initialized.

Example:

		C	Generate 100 values from the Poisson distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION LAMBDA
				INTEGER X(N)
				INTEGER LREF
				DOUBLE PRECISION REF(1000)
          
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) LAMBDA
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
          
		C	Initialize the reference vector
				LREF = 1000
				CALL DRANDPOISSONREFERENCE(LAMBDA,REF,LREF,INFO)
          
		C	Generate N variates from the Poisson distribution
				CALL DRANDGENERALDISCRETE(N,REF,STATE,X,INFO)
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDSKIPAHEAD / SRANDSKIPAHEAD

Advance a generator N places.

(Note that SRANDSKIPAHEAD is the single precision version of DRANDSKIPAHEAD. The argument lists of both routines are identical except that any double precision arguments of DRANDSKIPAHEAD are replaced in SRANDSKIPAHEAD by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDSKIPAHEAD (N,STATE,INFO)

— Input: INTEGER N

On input: number of places to skip ahead.
Constraint: N>=0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDSKIPAHEAD STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.


On output: The STATE vector for a generator that has been advanced N places.
Constraint: The STATE vector must be for either the NAG basic, Wichmann-Hill or L’Ecuyer Combined Recursive base generators.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.


Example:

          C     Generate 3 * 100 values from the Uniform distribution
          C     Multiple streams generated using the Skip Ahead method
                INTEGER LSTATE,N
                PARAMETER (LSTATE=16,N=100)
                INTEGER I,INFO,NSKIP
                INTEGER SEED(1),STATE1(LSTATE),STATE2(LSTATE),STATE3(LSTATE)
                INTEGER X1(N),X2(N),X3(N)
                DOUBLE PRECISION A,B
          
          C     Set the seed
                SEED(1) = 1234
          
          C     Set the distributional parameters
                A = 0.0D0
                B = 1.0D0
          
          C     Initialize the STATE1 vector
                CALL DRANDINITIALIZE(1,1,SEED,1,STATE1,LSTATE,INFO)
          
          C     Copy the STATE1 vector into other state vectors
                DO 20 I = 1,LSTATE
                  STATE2(I) = STATE1(I)
                  STATE3(I) = STATE1(I)
          20    CONTINUE
          
          C     Calculate how many places we want to skip, this
          C     should be >> than the number of variates we
          C     wish to generate from each stream
                NSKIP = N * N
          
          C     Advance each stream, first does not need changing
                CALL DRANDSKIPAHEAD(NSKIP,STATE2,INFO)
                CALL DRANDSKIPAHEAD(2*NSKIP,STATE3,INFO)
          
          C     Generate 3 sets of N variates from the Univariate distribution
                CALL DRANDUNIFORM(N,A,B,STATE1,X1,LDX,INFO)
                CALL DRANDUNIFORM(N,A,B,STATE2,X2,LDX,INFO)
                CALL DRANDUNIFORM(N,A,B,STATE3,X3,LDX,INFO)
          
          C     Print the results
                DO 40 I = 1,N
                  WRITE(6,*) X1(I),X2(I),X3(I)
          40    CONTINUE


DRANDSTUDENTST / SRANDSTUDENTST

Generates a vector of random variates from a Students T distribution with probability density function, f(X), where: f(X) = [((v – 1) / 2)!] / [(v / 2 – 1)! * sqrt(Pi * v) * (1 + X^2 / v)^((v + 1) / 2). Here v is the degrees of freedom, DF.

(Note that SRANDSTUDENTST is the single precision version of DRANDSTUDENTST. The argument lists of both routines are identical except that any double precision arguments of DRANDSTUDENTST are replaced in SRANDSTUDENTST by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDSTUDENTST (N,DF,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: INTEGER DF

On input: degrees of freedom.
Constraint: DF>0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDSTUDENTST STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Students T distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				INTEGER DF
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) DF
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Students T distribution
				CALL DRANDSTUDENTST(N,DF,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDTRIANGULAR / SRANDTRIANGULAR

Generates a vector of random variates from a Triangular distribution with probability density function, f(X), where: f(X) = [2 * (X – XMIN)] / [(XMAX – XMIN) * (XMED – XMIN)], if XMIN < X <= XMED, else f(X) = [2 * (XMAX – X)] / [(XMAX – XMIN) * (XMAX – XMED)], if XMED < X <= XMAX, otherwise f(X) = 0.

(Note that SRANDTRIANGULAR is the single precision version of DRANDTRIANGULAR. The argument lists of both routines are identical except that any double precision arguments of DRANDTRIANGULAR are replaced in SRANDTRIANGULAR by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDTRIANGULAR (N,XMIN,XMED,XMAX,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION XMIN

On input: minimum value for the distribution.

— Input: DOUBLE PRECISION XMED

On input: median value for the distribution.
Constraint: XMIN<=XMED <=XMAX.

— Input: DOUBLE PRECISION XMAX

On input: maximum value for the distribution.
Constraint: XMAX>=XMIN.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDTRIANGULAR STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Triangular distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION XMIN,XMAX,XMED
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) XMIN,XMAX,XMED
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Triangular distribution
				CALL DRANDTRIANGULAR(N,XMIN,XMAX,XMED,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDUNIFORM / SRANDUNIFORM

Generates a vector of random variates from a Uniform distribution with probability density function, f(X), where: f(X) = 1 / (B – A).

(Note that SRANDUNIFORM is the single precision version of DRANDUNIFORM. The argument lists of both routines are identical except that any double precision arguments of DRANDUNIFORM are replaced in SRANDUNIFORM by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDUNIFORM (N,A,B,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION A

On input: minimum value for the distribution.

— Input: DOUBLE PRECISION B

On input: maximum value for the distribution.
Constraint: B>=A.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDUNIFORM STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Uniform distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION A,B
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) A,B
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Uniform distribution
				CALL DRANDUNIFORM(N,A,B,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDVONMISES / SRANDVONMISES

Generates a vector of random variates from a Von Mises distribution with probability density function, f(X), where: f(X) = [exp(k * cos(X))] / [2 * Pi * I0(k)], where X is reduced modulo 2 * Pi, so that it lies between +/- Pi, and k is the concentration parameter VK.

(Note that SRANDVONMISES is the single precision version of DRANDVONMISES. The argument lists of both routines are identical except that any double precision arguments of DRANDVONMISES are replaced in SRANDVONMISES by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDVONMISES (N,VK,,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION VK

On input: concentration parameter.
Constraint: VK>0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDVONMISES STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Von Mises distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION VK
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) VK
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Von Mises distribution
				CALL DRANDVONMISES(N,VK,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DRANDWEIBULL / SRANDWEIBULL

Generates a vector of random variates from a Weibull distribution with probability density function, f(X), where: f(X) = [A * X^(A – 1) * exp(-(X^A) / B)] / B, if X > 0, otherwise f(X)=0.

(Note that SRANDWEIBULL is the single precision version of DRANDWEIBULL. The argument lists of both routines are identical except that any double precision arguments of DRANDWEIBULL are replaced in SRANDWEIBULL by single precision arguments – type REAL in FORTRAN or type float in C).

SUBROUTINE:DRANDWEIBULL (N,A,B,STATE,X,INFO)

— Input: INTEGER N

On input: number of variates required.
Constraint: N>=0.

— Input: DOUBLE PRECISION A

On input: shape parameter for the distribution.
Constraint: A>0.

— Input: DOUBLE PRECISION B

On input: scale parameter for the distribution.
Constraint: B>0.

— Input/Output: INTEGER STATE(*)

The STATE vector holds information on the state of the base generator being used and as such its minimum length varies. Prior to calling DRANDWEIBULL STATE must have been initialized. See Initialization of the Base Generators for information on initialization of the STATE variable.
On input: the current state of the base generator.
On output: the updated state of the base generator.

— Output: DOUBLE PRECISION X(N)

On output: vector of variates from the specified distribution.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

		C	Generate 100 values from the Weibull distribution
				INTEGER LSTATE,N
				PARAMETER (LSTATE=16,N=100)
				INTEGER I,INFO,SEED(1),STATE(LSTATE)
				DOUBLE PRECISION A,B
				DOUBLE PRECISION X(N)
		C	Set the seed
				SEED(1) = 1234
          
		C	Read in the distributional parameters
				READ(5,*) A,B
          
		C	Initialize the STATE vector
				CALL DRANDINITIALIZE(1,1,SEED,1,STATE,LSTATE,INFO)
		C	Generate N variates from the Weibull distribution
				CALL DRANDWEIBULL(N,A,B,STATE,X,INFO)
          
		C	Print the results
				WRITE(6,*) (X(I),I=1,N)

DZFFT Routine Documentation

SUBROUTINE:DZFFT (MODE,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by DZFFT.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real transform is performed. Initializations are assumed to have been performed by a prior call to DZFFT.
  • MODE=2 : (default) initializations and a real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER N

On input: N is the length of the real sequence X

— Input/Output: DOUBLE PRECISION X(N)

On input: X contains the real sequence of length N to be transformed.
On output: X contains the transformed Hermitian sequence.

— Input/Output: DOUBLE PRECISION COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT(0,N,X,COMM,INFO)
             CALL DZFFT(1,N,X,COMM,INFO)
             DO 10 I = N/2+2, N
                X(I) = -X(I)
        10   CONTINUE
             CALL ZDFFT(2,N,X,COMM,INFO)

DZFFT1D Routine Documentation

SUBROUTINE:DZFFT1D (MODE,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by DZFFT1D.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real-to-complex transform is performed. Initializations are assumed to have been performed by a prior call to DZFFT1D.
  • MODE=2 : (default) initializations and a real-to-complex transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER N

On input: N is length of the sequence in X

— Input: DOUBLE PRECISION X(N)

On input: X contains the real sequence of length N to be transformed.

— Output: DOUBLE COMPLEX Y(N/2+1)

On output: Y contains the tansformed complex sequence with roughly half redundant information due to complex conjugate removed.

— Input/Output: DOUBLE PRECISION COMM(4*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT1D(0,N,X,Y,COMM,INFO)
             CALL DZFFT1D(1,N,X,Y,COMM,INFO)
             DO 10 I = 1, N/2+1
                Y(I) = -Y(I)*EXP(-DBLE(I-1)/DBLE(N))
             10 CONTINUE
             CALL ZDFFT1D(2,N,Y,X,COMM,INFO)

DZFFT1M Routine Documentation

SUBROUTINE:DZFFT1M (MODE,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by DZFFT1M.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real-to-complex transform is performed. Initializations are assumed to have been performed by a prior call to DZFFT1M.
  • MODE=2 : (default) initializations and a real-to-complex transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is length of the sequence in X

— Input: DOUBLE PRECISION X(N*M)

On input: X contains the M real sequences of length N to be transformed. Element i of sequence j is stored in location i + (j – 1) * N of X.

— Output: DOUBLE COMPLEX Y((N/2+1)*M)

On output: Y contains the transformed Hermitian sequences in complex-Hermitian storage.

— Input/Output: DOUBLE PRECISION COMM(4*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT1M(0,M,N,X,Y,COMM,INFO)
             CALL DZFFT1M(1,M,N,X,Y,COMM,INFO)
             DO J = 1, M
                DO I = 1, N/2+1
                   Y(I,J) = -Y(I,J)*EXP(-DBLE(I-1)/DBLE(N))
                END DO
             EMD DO
             CALL ZDFFT1M(2,M,N,Y,X,COMM,INFO)

DZFFT2D Routine Documentation

SUBROUTINE:DZFFT2D (MODE,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by DZFFT2D.
On input:

  • MODE=0 : only default initializations (specific to M and N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real-to-complex transform is performed. Initializations are assumed to have been performed by a prior call to DZFFT2D.
  • MODE=2 : (default) initializations and a real-to-complex transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER M

On input: M is the number of rows in the 2D array of data to be transformed. If X is declared as a 2D array then M is the first dimension of X.

— Input: INTEGER N

On input: N is the number of columns in the 2D array of data to be transformed. If X is declared as a 2D array then N is the second dimension of X.

— Input: DOUBLE PRECISION X(M*N)

On input: X contains the M by N real 2D array to be transformed. Element ij is stored in location i + (j – 1) * M of X.

— Output: DOUBLE COMPLEX Y((M/2+1)*N)

On output: Y contains the transformed Hermitian sequences in complex-Hermitian storage.

— Input/Output: DOUBLE PRECISION COMM(4*M+6*N+300)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length M*N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT2D(0,M,N,X,Y,COMM,INFO)
             CALL DZFFT2D(1,M,N,X,Y,COMM,INFO)
             DO J = 1, N
                DO I = 1, M/2+1
                   Y(I,J) = -Y(I,J)/SQRT(DBLE(M*N))
                END DO
             EMD DO
             CALL ZDFFT2D(2,M,N,Y,X,COMM,INFO)

DZFFT3D Routine Documentation

SUBROUTINE:DZFFT3D (MODE,L,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by DZFFT3D.
On input:

  • MODE=0 : only default initializations (specific to L, M and N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real-to-complex transform is performed. Initializations are assumed to have been performed by a prior call to DZFFT3D.
  • MODE=2 : (default) initializations and a real-to-complex transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER L

On input: L is the length of the first dimension of the 3D array of data to be transformed. If X is declared as a 3D array then L is the first dimension of X.

— Input: INTEGER M

On input: M is the length of the second dimension of the 3D array of data to be transformed. If X is declared as a 3D array then L is the second dimension of X.

— Input: INTEGER N

On input: N is the length of the third dimension of the 3D array of data to be transformed. If X is declared as a 3D array then L is the third dimension of X.

— Input: DOUBLE PRECISION X(L*M*N)

On input: X contains the L by M by N real 3D array to be transformed. Element ijk is stored in location i + (j – 1)* L + (k – 1) * L * M of X.

— Output: DOUBLE COMPLEX Y((L/2+1)*M*N)

On output: Y contains the transformed Hermitian sequences in complex-Hermitian storage.

— Input/Output: DOUBLE PRECISION COMM(4*L+6*M+6*N+500)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length L*M*N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT3D(0,L,M,N,X,Y,COMM,INFO)
             CALL DZFFT3D(1,L,M,N,X,Y,COMM,INFO)
             DO K = 1, N
                DO J = 1, M
                   DO I = 1, L/2+1
                      Y(I,J,K) = -Y(I,J,K)/SQRT(DBLE(L*M*N))
                   END DO
                END DO
             EMD DO
             CALL ZDFFT3D(2,L,M,N,Y,X,COMM,INFO)

DZFFTM Routine Documentation

SUBROUTINE:DZFFTM (M,N,X,COMM,INFO)

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is the length of the real sequences in X

— Input/Output: DOUBLE PRECISION X(N*M)

On input: X contains the M real sequences of length N to be transformed. Element i of sequence j is stored in location i+(j-1)*N of X.
On output: X contains the transformed Hermitian sequences.

— Input/Output: DOUBLE PRECISION COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFTM(1,N,X,COMM,INFO)
             CALL DZFFTM(2,N,X,COMM,INFO)
             DO 10 I = 1, N
                X(I,3) = X(I,1)*X(N-I+1,2)
        10   CONTINUE
             CALL ZDFFTM(2,N,X(1,3),COMM,INFO)

FORTRAN specifications:

SUBROUTINE:ACMLVERSION (MAJOR, MINOR, PATCH)

— INTEGER: MAJOR, MINOR, PATCH

SUBROUTINE:ACMLINFO ()


SCFFT Routine Documentation

SUBROUTINE:SCFFT (MODE,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by SCFFT.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real transform is performed. Initializations are assumed to have been performed by a prior call to SCFFT.
  • MODE=2 : (default) initializations and a real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER N

On input: N is the length of the real sequence X

— Input/Output: REAL X(N)

On input: X contains the real sequence of length N to be transformed.
On output: X contains the transformed Hermitian sequence.

— Input/Output: REAL COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT(0,N,X,COMM,INFO)
             CALL SCFFT(1,N,X,COMM,INFO)
             DO 10 I = N/2+2, N
                X(I) = -X(I)
        10   CONTINUE
             CALL CSFFT(2,N,X,COMM,INFO)

SCFFT1D Routine Documentation

SUBROUTINE:SCFFT1D (MODE,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by SCFFT1D.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real-to-complex transform is performed. Initializations are assumed to have been performed by a prior call to SCFFT1D.
  • MODE=2 : (default) initializations and a real-to-complex transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER N

On input: N is length of the sequence in X

— Input: REAL X(N)

On input: X contains the real sequence of length N to be transformed.

— Output: COMPLEX Y(N/2+1)

On output: Y contains the tansformed complex sequence with roughly half redundant information due to complex conjugate removed.

— Input/Output: REAL COMM(4*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT1D(0,N,X,Y,COMM,INFO)
             CALL SCFFT1D(1,N,X,Y,COMM,INFO)
             DO 10 I = 1, N/2+1
                Y(I) = -Y(I)*EXP(-REAL(I-1)/REAL(N))
             10 CONTINUE
             CALL CSFFT1D(2,N,Y,X,COMM,INFO)

SCFFT1M Routine Documentation

SUBROUTINE:SCFFT1M (MODE,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by SCFFT1M.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real-to-complex transform is performed. Initializations are assumed to have been performed by a prior call to SCFFT1M.
  • MODE=2 : (default) initializations and a real-to-complex transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is length of the sequence in X

— Input: REAL X(N*M)

On input: X contains the M real sequences of length N to be transformed. Element i of sequence j is stored in location i + (j – 1) * N of X.

— Output: COMPLEX Y((N/2+1)*M)

On output: Y contains the transformed Hermitian sequences in complex-Hermitian storage.

— Input/Output: REAL COMM(4*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT1M(0,M,N,X,Y,COMM,INFO)
             CALL SCFFT1M(1,M,N,X,Y,COMM,INFO)
             DO J = 1, M
                DO I = 1, N/2+1
                   Y(I,J) = -Y(I,J)*EXP(-REAL(I-1)/REAL(N))
                END DO
             EMD DO
             CALL CSFFT1M(2,M,N,Y,X,COMM,INFO)

SCFFT2D Routine Documentation

SUBROUTINE:SCFFT2D (MODE,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by SCFFT2D.
On input:

  • MODE=0 : only default initializations (specific to M and N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real-to-complex transform is performed. Initializations are assumed to have been performed by a prior call to SCFFT2D.
  • MODE=2 : (default) initializations and a real-to-complex transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER M

On input: M is the number of rows in the 2D array of data to be transformed. If X is declared as a 2D array then M is the first dimension of X.

— Input: INTEGER N

On input: N is the number of columns in the 2D array of data to be transformed. If X is declared as a 2D array then N is the second dimension of X.

— Input: REAL X(M*N)

On input: X contains the M by N real 2D array to be transformed. Element ij is stored in location i + (j – 1) * M of X.

— Output: COMPLEX Y((M/2+1)*N)

On output: Y contains the transformed Hermitian sequences in complex-Hermitian storage.

— Input/Output: REAL COMM(4*M+10*N+300)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length M*N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT2D(0,M,N,X,Y,COMM,INFO)
             CALL SCFFT2D(1,M,N,X,Y,COMM,INFO)
             DO J = 1, N
                DO I = 1, M/2+1
                   Y(I,J) = -Y(I,J)/SQRT(REAL(M*N))
                END DO
             EMD DO
             CALL CSFFT2D(2,M,N,Y,X,COMM,INFO)

SCFFT3D Routine Documentation

SUBROUTINE:SCFFT3D (MODE,L,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by SCFFT3D.
On input:

  • MODE=0 : only default initializations (specific to L, M and N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real-to-complex transform is performed. Initializations are assumed to have been performed by a prior call to SCFFT3D.
  • MODE=2 : (default) initializations and a real-to-complex transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER L

On input: L is the length of the first dimension of the 3D array of data to be transformed. If X is declared as a 3D array then L is the first dimension of X.

— Input: INTEGER M

On input: M is the length of the second dimension of the 3D array of data to be transformed. If X is declared as a 3D array then L is the second dimension of X.

— Input: INTEGER N

On input: N is the length of the third dimension of the 3D array of data to be transformed. If X is declared as a 3D array then L is the third dimension of X.

— Input: REAL X(L*M*N)

On input: X contains the L by M by N real 3D array to be transformed. Element ijk is stored in location i + (j – 1)* L + (k – 1) * L * M of X.

— Output: COMPLEX Y((L/2+1)*M*N)

On output: Y contains the transformed Hermitian sequences in complex-Hermitian storage.

— Input/Output: REAL COMM(4*L+10*M+10*N+500)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length L*M*N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFT3D(0,L,M,N,X,Y,COMM,INFO)
             CALL SCFFT3D(1,L,M,N,X,Y,COMM,INFO)
             DO K = 1, N
                DO J = 1, M
                   DO I = 1, L/2+1
                      Y(I,J,K) = -Y(I,J,K)/SQRT(REAL(L*M*N))
                   END DO
                END DO
             EMD DO
             CALL CSFFT3D(2,L,M,N,Y,X,COMM,INFO)

SCFFTM Routine Documentation

SUBROUTINE:SCFFTM (M,N,X,COMM,INFO)

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is the length of the real sequences in X

— Input/Output: REAL X(N*M)

On input: X contains the M real sequences of length N to be transformed. Element i of sequence j is stored in location i+(j-1)*N of X.
On output: X contains the transformed Hermitian sequences.

— Input/Output: REAL COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL SCFFTM(1,N,X,COMM,INFO)
             CALL SCFFTM(2,N,X,COMM,INFO)
             DO 10 I = 1, N
                X(I,3) = X(I,1)*X(N-I+1,2)
        10   CONTINUE
             CALL CSFFTM(1,N,X(1,3),COMM,INFO)

UGEN

Specification for a user supplied base generator.

SUBROUTINE:UGEN (N,STATE,X,INFO)

— Input: INTEGER N

On input: the number of random numbers to be generated.

— Input/Output: INTEGER STATE(*)

On input: the internal state of your generator.

— Output: DOUBLE PRECISION X(N)

On output: the array of N uniform distributed random numbers, each in the semi-open interval (0.0, 1.0] – i.e. 1.0 is a legitimate return value, but 0.0 is not.

— Output: INTEGER INFO

On output: a flag which you can use to signal an error in the call of UGEN – for example, if UGEN is called without being initialized by UINI.


UINI

Specification for a user supplied initialization routine.

SUBROUTINE:UINI (GENID,SUBID,SEED,LSEED,STATE,LSTATE,INFO)

— Input: INTEGER GENID

On input: the ID associated with the generator. It may be used for anything you like.

— Input: INTEGER SUBID

On input: the sub-ID associated with the generator. It may be used for anything you like.

— Input: INTEGER SEED(LSEED)

On input: an array containing the initial seed for your generator.

— Input/Output: INTEGER LSEED

On input: either the size of the SEED array, or a value < 1.
On output: if LSEED<1 on entry, LSEED must be set to the required size of the SEED array. This allows a caller of UINI to query the required size.

— Output: INTEGER STATE(LSTATE)

On output: if LSTATE<1 on entry, STATE should be unchanged.
Otherwise, STATE is a state vector holding internal details required by your generator. On exit from UINI, the array STATE must hold the following information:

STATE(1) = ESTATE, where ESTATE is your minimum allowed size of array STATE.

STATE(2) = MAGIC, where MAGIC is a magic number of your own choice. This can be used by your routine UGEN as a check that UINI has previously been called.

STATE(3) = GENID

STATE(4) = SUBID

STATE(5)STATE(ESTATE-1) = internal state values required by your generator routine UGEN; for example, the current value of your seed.

STATE(ESTATE) = MAGIC, i.e. the same value as STATE(2).

— Input/Output: INTEGER LSTATE

On input: either the size of the STATE array, or a value < 1.
On output: if LSTATE<1 on entry, LSTATE should be set to the required size of the STATE array, i.e. the value ESTATE as described above. This allows the caller of UINI to query the required size.
Constraint: either LSTATE<1 or LSTATE>=ESTATE .

— Output: INTEGER INFO

On output: an error code, to be used in whatever way you wish; for example to flag an incorrect argument to UINI. If no error is encountered, UINI must set INFO to 0.


ZDFFT Routine Documentation

SUBROUTINE:ZDFFT (MODE,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZDFFT.
On input:

  • MODE=0 : only initializations (specific to the values of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a real transform is performed. Initializations are assumed to have been performed by a prior call to ZDFFT.
  • MODE=2 : (default) initializations and a real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: INTEGER N

On input: N is length of the sequence in X

— Input/Output: DOUBLE PRECISION X(N)

On input: X contains the Hermitian sequence of length N to be transformed.
On output: X contains the transformed real sequence.

— Input/Output: DOUBLE PRECISION COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT(0,N,X,COMM,INFO)
             CALL DZFFT(1,N,X,COMM,INFO)
             DO 10 I = N/2+2, N
                X(I) = -X(I)
        10   CONTINUE
             CALL ZDFFT(2,N,X,COMM,INFO)

ZDFFT1D Routine Documentation

SUBROUTINE:ZDFFT1D (MODE,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZDFFT1D.
On input:

  • MODE=0 : only initializations (specific to the values of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a complex-to-real transform is performed. Initializations are assumed to have been performed by a prior call to ZDFFT1D.
  • MODE=2 : (default) initializations and a complex-to-real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: INTEGER N

On input: N is length of the sequence in Y

— Input: DOUBLE COMPLEX X(N/2+1)

On input: X contains the complex sequence to be transformed.

— Output: DOUBLE PRECISION Y(N)

On output: Y contains the transformed real sequence of length N.

— Input/Output: DOUBLE PRECISION COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT1D(0,N,X,Y,COMM,INFO)
             CALL DZFFT1D(1,N,X,Y,COMM,INFO)
             DO 10 I = 1, N/2+1
                Y(I) = -Y(I)*EXP(-DBLE(I-1)/DBLE(N))
             10 CONTINUE
             CALL ZDFFT1D(2,N,Y,X,COMM,INFO)

ZDFFT1M Routine Documentation

SUBROUTINE:ZDFFT1M (MODE,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZDFFT1M.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a complex-to-real transform is performed. Initializations are assumed to have been performed by a prior call to ZDFFT1M.
  • MODE=2 : (default) initializations and a complex-to-real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is length of the sequence in Y

— Input: DOUBLE COMPLEX X((N/2+1)*M)

On input: X contains the Hermitian sequences to be transformed in complex-Hermitian storage.

— Output: DOUBLE PRECISION Y(N*M)

On output: Y contains the M real sequences of length N transformed. Element i of sequence j is stored in location i + (j – 1) * N of Y.

— Input/Output: DOUBLE PRECISION COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT1M(0,M,N,X,Y,COMM,INFO)
             CALL DZFFT1M(1,M,N,X,Y,COMM,INFO)
             DO J = 1, M
                DO I = 1, N/2+1
                   Y(I,J) = -Y(I,J)*EXP(-DBLE(I-1)/DBLE(N))
                END DO
             EMD DO
             CALL ZDFFT1M(2,M,N,Y,X,COMM,INFO)

ZDFFT2D Routine Documentation

SUBROUTINE:ZDFFT2D (MODE,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZDFFT2D.
On input:

  • MODE=0 : only default initializations (specific to M and N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a complex-to-real transform is performed. Initializations are assumed to have been performed by a prior call to ZDFFT2D.
  • MODE=2 : (default) initializations and a complex-to-real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER M

On input: M is the number of rows in the 2D array of the real data obtained from the transform. If Y is declared as a 2D array then M is the first dimension of Y.

— Input: INTEGER N

On input: N is the number of columns in the 2D array of the real data obtained from the transform. If Y is declared as a 2D array then N is the second dimension of Y.

— Input: DOUBLE COMPLEX X((M/2+1)*N)

On input: X contains the Hermitian sequences in complex-Hermitian storage to be transformed.

— Output: DOUBLE PRECISION Y(M*N)

On output: Y contains the M by N real 2D array obtained from the transform. Element ij is stored in location i + (j – 1) * M of Y.

— Input/Output: DOUBLE PRECISION COMM(4*M+6*N+300)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length M*N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT2D(0,M,N,X,Y,COMM,INFO)
             CALL DZFFT2D(1,M,N,X,Y,COMM,INFO)
             DO J = 1, N
                DO I = 1, M/2+1
                   Y(I,J) = -Y(I,J)/SQRT(DBLE(M*N))
                END DO
             EMD DO
             CALL ZDFFT2D(2,M,N,Y,X,COMM,INFO)

ZDFFT3D Routine Documentation

SUBROUTINE:ZDFFT3D (MODE,L,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZDFFT3D.
On input:

  • MODE=0 : only default initializations (specific to L, M and N) are performed; this is usually followed by calls to the same routine with MODE=1.
  • MODE=1 : a complex-to-real transform is performed. Initializations are assumed to have been performed by a prior call to ZDFFT3D.
  • MODE=2 : (default) initializations and a complex-to-real transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER L

On input: L is the length of the first dimension of the 3D array of data obtained from the transform. If Y is declared as a 3D array then L is the first dimension of Y.

— Input: INTEGER M

On input: M is the length of the second dimension of the 3D array of data obtained from the transform. If Y is declared as a 3D array then L is the second dimension of Y.

— Input: INTEGER N

On input: N is the length of the third dimension of the 3D array of data obtained from the transform. If Y is declared as a 3D array then L is the third dimension of Y.

— Input: DOUBLE COMPLEX X((L/2+1)*M*N)

On input: X contains the Hermitian sequences in complex-Hermitian storage to be transformed.

— Output: DOUBLE PRECISION Y(L*M*N)

On output: Y contains the L by M by N real 3D array obtained from the transform. Element ijk is stored in location i + (j – 1) * L + (k – 1) * L * M of Y.

— Input/Output: DOUBLE PRECISION COMM(4*L+6*M+6*N+500)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length L*M*N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFT3D(0,L,M,N,X,Y,COMM,INFO)
             CALL DZFFT3D(1,L,M,N,X,Y,COMM,INFO)
             DO K = 1, N
                DO J = 1, M
                   DO I = 1, L/2+1
                      Y(I,J,K) = -Y(I,J,K)/SQRT(DBLE(L*M*N))
                   END DO
                END DO
             EMD DO
             CALL ZDFFT3D(2,L,M,N,Y,X,COMM,INFO)

ZDFFTM Routine Documentation

SUBROUTINE:ZDFFTM (M,N,X,COMM,INFO)

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is the length of the sequences in X

— Input/Output: DOUBLE PRECISION X(N*M)

On input: X contains the M Hermitian sequences of length N to be transformed. Element i of sequence j is stored in location i+(j-1)*N of X.
On output: X contains the transformed real sequences.

— Input/Output: DOUBLE PRECISION COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL DZFFTM(1,N,X,COMM,INFO)
             CALL DZFFTM(2,N,X,COMM,INFO)
             DO 10 I = 1, N
                X(I,3) = X(I,1)*X(N-I+1,2)
        10   CONTINUE
             CALL ZDFFTM(1,N,X(1,3),COMM,INFO)

ZFFT1D Routine Documentation

SUBROUTINE:ZFFT1D (MODE,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZFFT1D.
On input:

  • MODE=0 : only default initializations (specific to N) are performed; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT1D.
  • MODE=1 : a backward (reverse) transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT1D.
  • MODE=-2 : initializations and a forward transform are performed.
  • MODE=2 : initializations and a backward transform are performed.
  • MODE=100 : similar to MODE=0; only initializations are performed, but first a plan is generated. This plan is chosen based on the fastest FFT computation for a subset of all possible plans.

— Input: INTEGER N

On input: N is the length of the complex sequence X

— Input/Output: COMPLEX*16 X(N)

On input: X contains the complex sequence of length N to be transformed.
On output: X contains the transformed sequence.

— Input/Output: COMPLEX*16 COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL ZFFT1D(0,N,X,COMM,INFO)
             CALL ZFFT1D(-1,N,X,COMM,INFO)
             CALL ZFFT1D(-1,N,Y,COMM,INFO)
             DO 10 I = 1, N
                X(I) = X(I)*DCONJG(Y(I))
        10   CONTINUE
             CALL ZFFT1D(1,N,X,COMM,INFO)

ZFFT1DX Routine Documentation

SUBROUTINE:ZFFT1DX (MODE,SCALE,INPL,N,X,INCX,Y,INCY,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZFFT1DX.
On input:

  • MODE=0 : only initializations (specific to the value of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT1DX.
  • MODE=1 : a backward (reverse) transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT1DX.
  • MODE=-2 : (default) initializations and a forward transform are performed.
  • MODE=2 : (default) initializations and a backward transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: DOUBLE PRECISION SCALE

On input: SCALE is the scaling factor to apply to the output sequence

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequence; otherwise the output sequence is returned in Y.

— Input: INTEGER N

On input: N is the number of elements to be transformed

— Input/Output: COMPLEX*16 X(1+(N-1)*INCX)

On input: X contains the complex sequence of length N to be transformed, with the ith element stored in X(1+(i-1)*INCX).
On output: if INPL is .TRUE. then X contains the transformed sequence in the same locations as on input; otherwise X remains unchanged.

— Input: INTEGER INCX

On input: INCX is the increment used to store successive elements of a sequence in X.
Constraint: INCX > 0.

— Output: COMPLEX*16 Y(1+(N-1)*INCY)

On output: if INPL is .FALSE. then Y contains the transformed sequence, with the ith element stored in Y(1+(i-1)*INCY); otherwise Y is not referenced.

— Input: INTEGER INCY

On input: INCY is the increment used to store successive elements of a sequence in Y. If INPL is .TRUE. then INCY is not referenced.
Constraint: INCY > 0.

— Input/Output: COMPLEX*16 COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward FFTs are performed unscaled and in-place on contiguous
     C     vectors X and Y following initialization. Manipulations on
     C     resultant Fourier coefficients are stored in X which is then
     C     transformed back.
     C
             SCALE = 1.0D0
             INPL = .TRUE.
             CALL ZFFT1DX(0,SCALE,INPL,N,X,1,DUM,1,COMM,INFO)
             CALL ZFFT1DX(-1,SCALE,INPL,N,X,1,DUM,1,COMM,INFO)
             CALL ZFFT1DX(-1,SCALE,INPL,N,Y,1,DUM,1,COMM,INFO)
             DO 10 I = 1, N
                X(I) = X(I)*DCONJG(Y(I))/DBLE(N)
        10   CONTINUE
             CALL ZFFT1DX(1,SCALE,INPL,N,X,1,DUM,1,COMM,INFO)

ZFFT1M Routine Documentation

SUBROUTINE:ZFFT1M (MODE,M,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZFFT1M.
On input:

  • MODE=0 : only initializations (specific to the value of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : forward transforms are performed. Initializations are assumed to have been performed by a prior call to ZFFT1M.
  • MODE=1 : backward (reverse) transforms are performed. Initializations are assumed to have been performed by a prior call to ZFFT1M.
  • MODE=-2 : (default) initializations and forward transforms are performed.
  • MODE=2 : (default) initializations and backward transforms are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: INTEGER M

On input: M is the number of sequences to be transformed.

— Input: INTEGER N

On input: N is the length of the complex sequences in X

— Input/Output: COMPLEX*16 X(N*M)

On input: X contains the M complex sequences of length N to be transformed. Element i of sequence j is stored in location i+(j-1)*N of X.
On output: X contains the transformed sequences.

— Input/Output: COMPLEX*16 COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL ZFFT1M(0,1,N,X,COMM,INFO)
             CALL ZFFT1M(-1,2,N,X,COMM,INFO)
             DO 10 I = 1, N
                X(I,3) = X(I,1)*DCONJG(X(I,2))
                X(I,2) = DCMPLX(0.0D0,1.0D0)*X(I,2)
        10   CONTINUE
             CALL ZFFT1M(1,2,N,X(1,2),COMM,INFO)

ZFFT1MX Routine Documentation

SUBROUTINE:ZFFT1MX (MODE,SCALE,INPL,NSEQ,N,X,INCX1,INCX2,Y,INCY1,INCY2,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZFFT1MX.
On input:

  • MODE=0 : only initializations (specific to the value of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT1MX.
  • MODE=1 : a backward (reverse) transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT1MX.
  • MODE=-2 : (default) initializations and a forward transform are performed.
  • MODE=2 : (default) initializations and a backward transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the value of N) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the value of N.

— Input: DOUBLE PRECISION SCALE

On input: SCALE is the scaling factor to apply to the output sequences

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequences; otherwise the output sequences are returned in Y.

— Input: INTEGER NSEQ

On input: NSEQ is the number of sequences to be transformed

— Input: INTEGER N

On input: N is the number of elements in each sequence to be transformed

— Input/Output: COMPLEX*16 X(1+(N-1)*INCX1+(NSEQ-1)*INCX2)

On input: X contains the NSEQ complex sequences of length N to be transformed; the ith element of sequence j is stored in X(1+(i-1)*INCX1+(j-1)*INCX2).
On output: if INPL is .TRUE. then X contains the transformed sequences in the same locations as on input; otherwise X remains unchanged.

— Input: INTEGER INCX1

On input: INCX1 is the increment used to store successive elements of a given sequence in X (INCX1=1 for contiguous data).
Constraint: INCX1 > 0.

— Input: INTEGER INCX2

On input: INCX2 is the increment used to store corresponding elements of successive sequences in X (INCX2=N for contiguous data).
Constraint: INCX2 > 0.

— Output: COMPLEX*16 Y(1+(N-1)*INCY1+(NSEQ-1)*INCY2)

On output: if INPL is .FALSE. then Y contains the transformed sequences with the ith element of sequence j stored in Y(1+(i-1)*INCY1+(j-1)*INCY2); otherwise Y is not referenced.

— Input: INTEGER INCY1

On input: INCY1 is the increment used to store successive elements of a given sequence in Y. If INPL is .TRUE. then INCY1 is not referenced.
Constraint: INCY1 > 0.

— Input: INTEGER INCY2

On input: INCY2 is the increment used to store corresponding elements of successive sequences in Y (INCY2=N for contiguous data). If INPL is .TRUE. then INCY2 is not referenced.
Constraint: INCY2 > 0.

— Input/Output: COMPLEX*16 COMM(3*N+100)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence length N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward FFTs are performed unscaled and in-place on two
     C     contiguous vectors stored in the first two columns of X.
     C     Manipulations are stored in 2nd and 3rd columns of X which are
     C     then transformed back.
     C
             COMPLEX *16 X(N,3)
             SCALE = 1.0D0
             INPL = .TRUE.
             CALL ZFFT1MX(0,SCALE,INPL,2,N,X,1,N,DUM,1,N,COMM,INFO)
             CALL ZFFT1MX(-1,SCALE,INPL,2,N,X,1,N,DUM,1,N,COMM,INFO)
             DO 10 I = 1, N
                X(I,3) = X(I,1)*DCONJG(X(I,2))/DBLE(N)
                X(I,2) = DCMPLX(0.0D0,1.0D0)*X(I,2)/DBLE(N)
        10   CONTINUE
             CALL ZFFT1MX(1,SCALE,INPL,2,N,X(1,2),1,N,DUM,1,N,COMM,INFO)

ZFFT2D Routine Documentation

SUBROUTINE:ZFFT2D (MODE,M,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the direction of transform to be performed by ZFFT2D.
On input:

  • MODE=-1 : forward 2D transform is performed.
  • MODE=1 : backward (reverse) 2D transform is performed.

— Input: INTEGER M

On input: M is the number of rows in the 2D array of data to be transformed. If X is declared as a 2D array then M is the first dimension of X.

— Input: INTEGER N

On input: N is the number of columns in the 2D array of data to be transformed. If X is declared as a 2D array then M is the second dimension of X.

— Input/Output: COMPLEX*16 X(M*N)

On input: X contains the M by N complex 2D array to be transformed. Element ij is stored in location i+(j-1)*M of X.
On output: X contains the transformed sequence.

— Input/Output: COMPLEX*16 COMM(M*N+3*(M+N)+100)

COMM is a communication array used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL ZFFT2D(-1,M,N,X,COMM,INFO)
             DO 20 J = 1, N
                DO 10 I = 1, MIN(J-1,M)
                   X(I,J) = 0.5D0*(X(I,J) + X(J,I))
                   X(J,I) = DCONJG(X(I,J))
        10      CONTINUE
        20   CONTINUE
             CALL ZFFT2D(1,M,N,X,COMM,INFO)

ZFFT2DX Routine Documentation

SUBROUTINE:ZFFT2DX (MODE,SCALE,LTRANS,INPL,M,N,X,INCX1,INCX2,Y,INCY1,INCY2,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZFFT2DX.
On input:

  • MODE=0 : only initializations (specific to the value of N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward 2D transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT2DX.
  • MODE=1 : a backward (reverse) 2D transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT2DX.
  • MODE=-2 : (default) initializations and a forward 2D transform are performed.
  • MODE=2 : (default) initializations and a backward 2D transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the values of N and M) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the values of N and M.

— Input: DOUBLE PRECISION SCALE

On input: SCALE is the scaling factor to apply to the output sequences

— Input: LOGICAL LTRANS

On input: if LTRANS is .TRUE. then a normal final transposition is performed internally to return transformed data consistent with the values for arguments INPL, INCX1, INCX2, INCY1 and INCY2. If LTRANS is .FALSE. then the final transposition is not performed explicitly; the storage format on output is determined by whether the output data is stored contiguously or not – please see the output specifications for X and Y for details.

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequences; otherwise the output sequences are returned in Y.

— Input: INTEGER M

On input: M is the first dimension of the 2D transform.

— Input: INTEGER N

On input: N is the second dimension of the 2D transform.

— Input/Output: COMPLEX*16 X(1+(M-1)*INCX1+(N-1)*INCX2)

On input: X contains the M by N complex 2D data array to be transformed; the (ij)th element is stored in X(1+(i-1)*INCX1+(j-1)*INCX2).
On output: if INPL is .TRUE. then X contains the transformed data, either in the same locations as on input when LTRANS=.TRUE.; in locations X((i-1)*N+j) when LTRANS=.FALSE., INCX1=1 and INCX2=M; and otherwise in the same locations as on input. If INPL is .FALSE. X remains unchanged.

— Input: INTEGER INCX1

On input: INCX1 is the increment used to store, in X, successive elements in the first dimension (INCX1=1 for contiguous data).
Constraint: INCX1 > 0.

— Input: INTEGER INCX2

On input: INCX2 is the increment used to store, in X, successive elements in the second dimension (INCX2=M for contiguous data).
Constraint: INCX2 > 0;
INCX2 > (M-1)*INCX1 if N > 1.

— Output: COMPLEX*16 Y(1+(M-1)*INCY1+(N-1)*INCY2)

On output: if INPL is .FALSE. then Y contains the transformed data. If LTRANS=.TRUE. then the (ij)th data element is stored in Y(1+(i-1)*INCY1+(j-1)*INCY2); if LTRANS=.FALSE., INCY1=1 and INCY2=N then the (ij)th data element is stored in Y((i-1)*N+j); and otherwise the (ij)th element is stored in Y(1+(i-1)*INCY1+(j-1)*INCY2). If INPL is .TRUE. then Y is not referenced.

— Input: INTEGER INCY1

On input: INCY1 is the increment used to store successive elements in the first dimension in Y (INCY1=1 for contiguous data). If INPL is .TRUE. then INCY1 is not referenced.
Constraint: INCY1 > 0.

— Input: INTEGER INCY2

On input: INCY2 is the increment used to store successive elements in the second dimension in Y (for contiguous data, INCY2=M when LTRANS is .TRUE. or INCY2=N when LTRANS is .FALSE.). If INPL is .TRUE. then INCY2 is not referenced.
Constraints: INCY2 > 0;
INCY2 > (M-1)*INCY1 if N > 1 and LTRANS is .TRUE.;
INCY2 = N if M > 1 and LTRANS is .FALSE..

— Input/Output: COMPLEX*16 COMM(M*N+3*M+3*N+200)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same dimensions M and N. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward 2D FFT is performed unscaled, without final transpose
     C     and out-of-place on data stored in array X and output to Y.
     C     Manipulations are stored in vector Y which is then transformed
     C     back, with scaling, into the first M rows of X.
     C
             COMPLEX *16 X(M,N), Y(N,M)
             SCALE = 1.0D0
             INPL = .FALSE.
             LTRANS = .FALSE.
             CALL ZFFT2DX(0,SCALE,LTRANS,INPL,M,N,X,1,M,Y,1,N,COMM,INFO)
             CALL ZFFT2DX(-1,SCALE,LTRANS,INPL,M,N,X,1,M,Y,1,N,COMM,INFO)
             DO 20 I = M
                DO 10 J = 1, N
                   Y(J,I) = 0.5D0*Y(J,I)*EXP(0.001D0*(I+J-2))
        10      CONTINUE
        20   CONTINUE
             SCALE = 1.0D0/DBLE(M*N)
             CALL ZFFT2DX(1,SCALE,LTRANS,INPL,N,M,Y,1,N,X,1,M,COMM,INFO)

ZFFT3D Routine Documentation

SUBROUTINE:ZFFT3D (MODE,L,M,N,X,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the direction of transform to be performed by ZFFT3D.
On input:

  • MODE=-1 : forward 3D transform is performed.
  • MODE=1 : backward (reverse) 3D transform is performed.

— Input: INTEGER L

On input: the length of the first dimension of the 3D array of data to be transformed. If X is declared as a 3D array then L is the first dimension of X.

— Input: INTEGER M

On input: the length of the second dimension of the 3D array of data to be transformed. If X is declared as a 3D array then M is the second dimension of X.

— Input: INTEGER N

On input: the length of the third dimension of the 3D array of data to be transformed. If X is declared as a 3D array then N is the third dimension of X.

— Input/Output: COMPLEX*16 X(L*M*N)

On input: X contains the L by M by N complex 3D array to be transformed. Element ijk is stored in location i+(j-1)*L+(k-1)*L*M of X.
On output: X contains the transformed sequence.

— Input/Output: COMPLEX*16 COMM(L*M*N+3*(L+M+N))

COMM is a communication array used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

             CALL ZFFT3D(-1,L,M,N,X,COMM,INFO)
             DO 30 K = 1, N
               DO 20 J = 1, M
                 DO 10 I = 1, L
                   X(I,J) = X(I,J)*EXP(-0.001D0*DBLE(I+J+K))
        10       CONTINUE
        20     CONTINUE
        30   CONTINUE
             CALL ZFFT3D(1,L,M,N,X,COMM,INFO)

ZFFT3DX Routine Documentation

SUBROUTINE:ZFFT3DX (MODE,SCALE,LTRANS,INPL,L,M,N,X,Y,COMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZFFT3DX.
On input:

  • MODE=0 : only initializations (specific to the values of L, M and N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward 3D transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT3DX.
  • MODE=1 : a backward (reverse) 3D transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT3DX.
  • MODE=100 : similar to MODE=0; only initializations (specific to the values of L, M and M) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the values of L, M and N.

— Input: DOUBLE PRECISION SCALE

On input: SCALE is the scaling factor to apply to the output sequences

— Input: LOGICAL LTRANS

On input: if LTRANS is .TRUE. then a normal final transposition is performed internally to return transformed data using the same storage format as the input data. If LTRANS is .FALSE. then the final transposition is not performed and transformed data is stored, in X or Y, in transposed form.

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequences; otherwise the output sequences are returned in Y.

— Input: INTEGER L

On input: L is the first dimension of the 3D transform.

— Input: INTEGER M

On input: M is the second dimension of the 3D transform.

— Input: INTEGER N

On input: N is the third dimension of the 3D transform.

— Input/Output: COMPLEX*16 X(L*M*N)

On input: X contains the L by M by N complex 3D data array to be transformed; the (ijk)th element is stored in X(i+(j-1)*L+(k-1)*L*M).
On output: if INPL is .TRUE. then X contains the transformed data, either in the same locations as on input when LTRANS=.TRUE.; or in locations X(k+(j-1)*N+(i-1)*N*M) when LTRANS=.FALSE. If INPL is .FALSE. X remains unchanged.

— Output: COMPLEX*16 Y(L*M*N)

On output: if INPL is .FALSE. then Y contains the three-dimensional transformed data. If LTRANS=.TRUE. then the (ijk)th data element is stored in Y(i+(j-1)*L+(k-1)*L*M); otherwise, the (ijk)th data element is stored in Y(k+(j-1)*N+(i-1)*N*M). If INPL is .TRUE. then Y is not referenced.

— Input/Output: COMPLEX*16 COMM(L*M*N+3*(L+M+N)+300)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence dimensions. The remainder is used as temporary store.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward 3D FFT is performed unscaled, without final transpose
     C     and out-of-place on data stored in array X and output to Y.
     C     Manipulations are stored in vector Y which is then transformed
     C     back, with scaling, into the first M rows of X.
     C
             COMPLEX *16 X(L*M*N), Y(L*M*N)
             SCALE = 1.0D0
             INPL = .FALSE.
             LTRANS = .FALSE.
             CALL ZFFT3DX(0,SCALE,LTRANS,INPL,L,M,N,X,Y,COMM,INFO)
             CALL ZFFT3DX(-1,SCALE,LTRANS,INPL,L,M,N,X,Y,COMM,INFO)
             IY = 1
             DO 20 I = 1, L
                DO 40 J = 1, M
                   DO 10 K = 1, N
                      Y(IY) = Y(IY)*EXP(-0.001D0*DBLE(I+J+K-3))
                      IY = IY + 1
        10      CONTINUE
        20   CONTINUE
             SCALE = 1.0D0/DBLE(L*M*N)
             CALL ZFFT3DX(1,SCALE,LTRANS,INPL,N,M,L,Y,X,COMM,INFO)

ZFFT3DY Routine Documentation

SUBROUTINE:ZFFT3DY (MODE,SCALE,INPL,L,M,N,X,INCX1,INCX2,INCX3,Y,INCY1,INCY2,INCY3,COMM,LCOMM,INFO)

— Input: INTEGER MODE

The value of MODE on input determines the operation performed by ZFFT3DY.
On input:

  • MODE=0 : only initializations (specific to the values of L, M and N) are performed using a default plan; this is usually followed by calls to the same routine with MODE=-1 or 1.
  • MODE=-1 : a forward 3D transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT3DY.
  • MODE=1 : a backward (reverse) 3D transform is performed. Initializations are assumed to have been performed by a prior call to ZFFT3DY.
  • MODE=-2 : (default) initializations and a forward 3D transform are performed.
  • MODE=2 : (default) initializations and a backward 3D transform are performed.
  • MODE=100 : similar to MODE=0; only initializations (specific to the values of L, M and M) are performed, but these are based on a plan that is first generated by timing a subset of all possible plans and choosing the quickest (i.e. the FFT computation was timed as fastest based on the chosen plan). The plan generation phase may take a significant amount of time depending on the values of L, M and N.

— Input: REAL SCALE

On input: SCALE is the scaling factor to apply to the output sequences

— Input: LOGICAL INPL

On input: if INPL is .TRUE. then X is overwritten by the output sequences; otherwise the output sequences are returned in Y.

— Input: INTEGER L

On input: L is the first dimension of the 3D transform.

— Input: INTEGER M

On input: M is the second dimension of the 3D transform.

— Input: INTEGER N

On input: N is the third dimension of the 3D transform.

— Input/Output: COMPLEX*16 X(*)

On input: X contains the L by M by N complex 3D data array to be transformed; the (ijk)th element is stored in X(1+(i-1)*INCX1+(j-1)*INCX2+(k-1)*INCX3).
On output: if INPL is .TRUE. then X contains the transformed data in the same locations as on input. If INPL is .FALSE. X remains unchanged.

— Input: INTEGER INCX1

On input: INCX1 is the step in index of X between successive data elements in the first dimension of the 3D data. Usually INCX1=1 so that succesive elements in the first dimension are stored contiguously.
Constraint: INCX1 > 0.

— Input: INTEGER INCX2

On input: INCX2 is the step in index of X between successive data elements in the second dimension of the 3D data. For completely contiguous data (no gaps in X) INCX2 should be set to L.
Constraint: INCX2 > 0;
INCX2 > (L-1)*INCX1 if max(M,N) > 1.

— Input: INTEGER INCX3

On input: INCX3 is the step in index of X between successive data elements in the third dimension of the 3D data. For completely contiguous data (no gaps in X) INCX3 should be set to L*M.
Constraint: INCX3 > 0;
INCX3 > (L-1)*INCX1+(M-1)*INCX2 if N > 1.

— Output: COMPLEX*16 Y(*)

On output: if INPL is .FALSE. then Y contains the three-dimensional transformed data. If LTRANS=.TRUE. then the the (ijk)th element is stored in Y(1+(i-1)*INCY1+(j-1)*INCY2+(k-1)*INCY3).
If INPL is .TRUE. then Y is not referenced.

— Input: INTEGER INCY1

On input: if INPL is .FALSE. then INCY1 is the step in index of Y between successive data elements in the first dimension of the 3D transformed data. Usually INCY1=1 so that succesive elements in the first dimension are stored contiguously.
If INPL is .TRUE. then INCY1 is not referenced. Constraint: If INPL is .FALSE. then INCY1 > 0.

— Input: INTEGER INCY2

On input: if INPL is .FALSE. then INCY2 is the step in index of Y between successive data elements in the second dimension of the 3D transformed data. For completely contiguous data (no gaps in Y) INCY2 should be set to L.
Constraint: INCY2 > 0 if INPL is .FALSE.;
INCY2 > (L-1)*INCY1, if INPL is .FALSE. and max(M,N) > 1.

— Input: INTEGER INCY3

On input: if INPL is .FALSE. then INCY3 is the step in index of Y between successive data elements in the third dimension of the 3D transformed data. For completely contiguous data (no gaps in Y) INCY3 should be set to L*M.
Constraint: INCY3 > 0 if INPL is .FALSE.;
INCY3 > (L-1)*INCY1+(M-1)*INCY2, if INPL is .FALSE. and N > 1.

— Input/Output: COMPLEX*16 COMM(LCOMM)

COMM is a communication array. Some portions of the array are used to store initializations for subsequent calls with the same sequence dimensions. The remainder is used as temporary store; if this is not sufficient for the requirements of the routine then temporary storage space will be dynamically allocated internally.

— Input: INTEGER LCOMM

On input: LCOMM is the length of the communication array COMM. The amount of internal dynamic allocation of temporary storage can be reduced significantly by declaring COMM to be of length at least L*M*N + 4*(L+M+N) + 300.
Constraint: LCOMM > 3*(L+M+N) + 150.

— Output: INTEGER INFO

On output: INFO is an error indicator. On successful exit, INFO contains 0. If INFO = -i on exit, the i-th argument had an illegal value.

Example:

     C     Forward 3D FFT is performed unscaled and in-place, on the leading
     C     10x10x10 submatrix of a larger 100x100x100 array of data.
     C     The result is transformed back with scaling.
     C
             SCALE = 1.0D0
             INPL = .TRUE.
             L = 10
             M = 10
             N = 10
             LCOMM = 2000000
             CALL ZFFT3DY(0,SCALE,INPL,L,M,N,X,1,100,10000,Y,1,1,1,
          *               COMM,LCOMM,INFO)
             CALL ZFFT3DY(-1,SCALE,INPL,L,M,N,X,1,100,10000,Y,1,1,1,
          *               COMM,LCOMM,INFO)
             IY = 1
             DO 20 I = 1, L
                DO 40 J = 1, M
                   DO 10 K = 1, N
                      X(I,J,K) = X(I,J,K)*EXP(-1.0D-3*DBLE(I+J+K-3))
        10      CONTINUE
        20   CONTINUE
             SCALE = 1.0/DBLE(L*M*N)
             CALL ZFFT3DY(1,SCALE,INPL,L,M,N,X,1,100,10000,Y,1,1,1,
          *               COMM,LCOMM,INFO)