subroutine dmatvecprod(prod, arg1, arg2, size)
    implicit none
    integer size, i
    real(kind=8) prod(size), arg1(size), arg2(size)

    do i = 1, size
        prod(i) = arg1(i) * arg2(i)
    end do
    return
end
subroutine dotprod(arg1, arg2, size, dot)
    implicit none
    integer size, i
    real(kind=8) dot, arg1(size), arg2(size), ddot
    !      dot=ddot(size,arg1,1,arg2,1)
    dot = 0.d0
    do i = 1, size
        dot = dot + arg1(i) * arg2(i)
    end do
    return
end
subroutine cdotprod(arg1, arg2, size, dot)
    implicit none
    integer size, i
    complex(kind=8) dot, arg1(size), arg2(size), zdotu
    !      dot=zdotu(size,arg1,1,arg2,1)
    dot = (0.d0, 0.d0)
    do i = 1, size
        dot = dot + dconjg(arg1(i)) * arg2(i)
    end do
    return
end


text:
      SUBROUTINE TG01FD( COMPQ, COMPZ, JOBA, L, N, M, P, A, LDA, E, LDE,
     $                   B, LDB, C, LDC, Q, LDQ, Z, LDZ, RANKE, RNKA22,
     $                   TOL, IWORK, DWORK, LDWORK, INFO )
C
C     SLICOT RELEASE 5.7.
C
C     Copyright (c) 2002-2020 NICONET e.V.
C
C     PURPOSE
C
C     To compute for the descriptor system (A-lambda E,B,C)
C     the orthogonal transformation matrices Q and Z such that the
C     transformed system (Q'*A*Z-lambda Q'*E*Z, Q'*B, C*Z) is
C     in a SVD-like coordinate form with
C
C                  ( A11  A12 )             ( Er  0 )
C         Q'*A*Z = (          ) ,  Q'*E*Z = (       ) ,
C                  ( A21  A22 )             (  0  0 )
C
C     where Er is an upper triangular invertible matrix.
C     Optionally, the A22 matrix can be further reduced to the form
C
C                  ( Ar  X )
C            A22 = (       ) ,
C                  (  0  0 )
C
C     with Ar an upper triangular invertible matrix, and X either a full
C     or a zero matrix.
C     The left and/or right orthogonal transformations performed
C     to reduce E and A22 can be optionally accumulated.
C
C     ARGUMENTS
C
C     Mode Parameters
C
C     COMPQ   CHARACTER*1
C             = 'N':  do not compute Q;
C             = 'I':  Q is initialized to the unit matrix, and the
C                     orthogonal matrix Q is returned;
C             = 'U':  Q must contain an orthogonal matrix Q1 on entry,
C                     and the product Q1*Q is returned.
C
C     COMPZ   CHARACTER*1
C             = 'N':  do not compute Z;
C             = 'I':  Z is initialized to the unit matrix, and the
C                     orthogonal matrix Z is returned;
C             = 'U':  Z must contain an orthogonal matrix Z1 on entry,
C                     and the product Z1*Z is returned.
C
C     JOBA    CHARACTER*1
C             = 'N':  do not reduce A22;
C             = 'R':  reduce A22 to a SVD-like upper triangular form.
C             = 'T':  reduce A22 to an upper trapezoidal form.
C
C     Input/Output Parameters
C
C     L       (input) INTEGER
C             The number of rows of matrices A, B, and E.  L >= 0.
C
C     N       (input) INTEGER
C             The number of columns of matrices A, E, and C.  N >= 0.
C
C     M       (input) INTEGER
C             The number of columns of matrix B.  M >= 0.
C
C     P       (input) INTEGER
C             The number of rows of matrix C.  P >= 0.
C
C     A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
C             On entry, the leading L-by-N part of this array must
C             contain the state dynamics matrix A.
C             On exit, the leading L-by-N part of this array contains
C             the transformed matrix Q'*A*Z. If JOBA = 'T', this matrix
C             is in the form
C
C                           ( A11  *   *  )
C                  Q'*A*Z = (  *   Ar  X  ) ,
C                           (  *   0   0  )
C
C             where A11 is a RANKE-by-RANKE matrix and Ar is a
C             RNKA22-by-RNKA22 invertible upper triangular matrix.
C             If JOBA = 'R' then A has the above form with X = 0.
C
C     LDA     INTEGER
C             The leading dimension of array A.  LDA >= MAX(1,L).
C
C     E       (input/output) DOUBLE PRECISION array, dimension (LDE,N)
C             On entry, the leading L-by-N part of this array must
C             contain the descriptor matrix E.
C             On exit, the leading L-by-N part of this array contains
C             the transformed matrix Q'*E*Z.
C
C                      ( Er  0 )
C             Q'*E*Z = (       ) ,
C                      (  0  0 )
C
C             where Er is a RANKE-by-RANKE upper triangular invertible
C             matrix.
C
C     LDE     INTEGER
C             The leading dimension of array E.  LDE >= MAX(1,L).
C
C     B       (input/output) DOUBLE PRECISION array, dimension (LDB,M)
C             On entry, the leading L-by-M part of this array must
C             contain the input/state matrix B.
C             On exit, the leading L-by-M part of this array contains
C             the transformed matrix Q'*B.
C
C     LDB     INTEGER
C             The leading dimension of array B.
C             LDB >= MAX(1,L) if M > 0 or LDB >= 1 if M = 0.
C
C     C       (input/output) DOUBLE PRECISION array, dimension (LDC,N)
C             On entry, the leading P-by-N part of this array must
C             contain the state/output matrix C.
C             On exit, the leading P-by-N part of this array contains
C             the transformed matrix C*Z.
C
C     LDC     INTEGER
C             The leading dimension of array C.  LDC >= MAX(1,P).
C
C     Q       (input/output) DOUBLE PRECISION array, dimension (LDQ,L)
C             If COMPQ = 'N':  Q is not referenced.
C             If COMPQ = 'I':  on entry, Q need not be set;
C                              on exit, the leading L-by-L part of this
C                              array contains the orthogonal matrix Q,
C                              where Q' is the product of Householder
C                              transformations which are applied to A,
C                              E, and B on the left.
C             If COMPQ = 'U':  on entry, the leading L-by-L part of this
C                              array must contain an orthogonal matrix
C                              Q1;
C                              on exit, the leading L-by-L part of this
C                              array contains the orthogonal matrix
C                              Q1*Q.
C
C     LDQ     INTEGER
C             The leading dimension of array Q.
C             LDQ >= 1,        if COMPQ = 'N';
C             LDQ >= MAX(1,L), if COMPQ = 'U' or 'I'.
C
C     Z       (input/output) DOUBLE PRECISION array, dimension (LDZ,N)
C             If COMPZ = 'N':  Z is not referenced.
C             If COMPZ = 'I':  on entry, Z need not be set;
C                              on exit, the leading N-by-N part of this
C                              array contains the orthogonal matrix Z,
C                              which is the product of Householder
C                              transformations applied to A, E, and C
C                              on the right.
C             If COMPZ = 'U':  on entry, the leading N-by-N part of this
C                              array must contain an orthogonal matrix
C                              Z1;
C                              on exit, the leading N-by-N part of this
C                              array contains the orthogonal matrix
C                              Z1*Z.
C
C     LDZ     INTEGER
C             The leading dimension of array Z.
C             LDZ >= 1,        if COMPZ = 'N';
C             LDZ >= MAX(1,N), if COMPZ = 'U' or 'I'.
C
C     RANKE   (output) INTEGER
C             The estimated rank of matrix E, and thus also the order
C             of the invertible upper triangular submatrix Er.
C
C     RNKA22  (output) INTEGER
C             If JOBA = 'R' or 'T', then RNKA22 is the estimated rank of
C             matrix A22, and thus also the order of the invertible
C             upper triangular submatrix Ar.
C             If JOBA = 'N', then RNKA22 is not referenced.
C
C     Tolerances
C
C     TOL     DOUBLE PRECISION
C             The tolerance to be used in determining the rank of E
C             and of A22. If the user sets TOL > 0, then the given
C             value of TOL is used as a lower bound for the
C             reciprocal condition numbers of leading submatrices
C             of R or R22 in the QR decompositions E * P = Q * R of E
C             or A22 * P22 = Q22 * R22 of A22.
C             A submatrix whose estimated condition number is less than
C             1/TOL is considered to be of full rank.  If the user sets
C             TOL <= 0, then an implicitly computed, default tolerance,
C             defined by  TOLDEF = L*N*EPS,  is used instead, where
C             EPS is the machine precision (see LAPACK Library routine
C             DLAMCH). TOL < 1.
C
C     Workspace
C
C     IWORK   INTEGER array, dimension (N)
C
C     DWORK   DOUBLE PRECISION array, dimension (LDWORK)
C             On exit, if INFO = 0, DWORK(1) returns the optimal value
C             of LDWORK.
C
C     LDWORK  INTEGER
C             The length of the array DWORK.
C             LDWORK >= MAX( 1, N+P, MIN(L,N)+MAX(3*N-1,M,L) ).
C             For optimal performance, LDWORK should be larger.
C
C             If LDWORK = -1, then a workspace query is assumed;
C             the routine only calculates the optimal size of the
C             DWORK array, returns this value as the first entry of
C             the DWORK array, and no error message related to LDWORK
C             is issued by XERBLA.
C
C     Error Indicator
C
C     INFO    INTEGER
C             = 0:  successful exit;
C             < 0:  if INFO = -i, the i-th argument had an illegal
C                   value.
C
C     METHOD
C
C     The routine computes a truncated QR factorization with column
C     pivoting of E, in the form
C
C                       ( E11 E12 )
C           E * P = Q * (         )
C                       (  0  E22 )
C
C     and finds the largest RANKE-by-RANKE leading submatrix E11 whose
C     estimated condition number is less than 1/TOL. RANKE defines thus
C     the rank of matrix E. Further E22, being negligible, is set to
C     zero, and an orthogonal matrix Y is determined such that
C
C           ( E11 E12 ) = ( Er  0 ) * Y .
C
C     The overal transformation matrix Z results as Z = P * Y' and the
C     resulting transformed matrices Q'*A*Z and Q'*E*Z have the form
C
C                          ( Er  0 )                      ( A11  A12 )
C         E <- Q'* E * Z = (       ) ,  A <- Q' * A * Z = (          ) ,
C                          (  0  0 )                      ( A21  A22 )
C
C     where Er is an upper triangular invertible matrix.
C     If JOBA = 'R' the same reduction is performed on A22 to obtain it
C     in the form
C
C                  ( Ar  0 )
C            A22 = (       ) ,
C                  (  0  0 )
C
C     with Ar an upper triangular invertible matrix.
C     If JOBA = 'T' then A22 is row compressed using the QR
C     factorization with column pivoting to the form
C
C                  ( Ar  X )
C            A22 = (       )
C                  (  0  0 )
C
C     with Ar an upper triangular invertible matrix.
C
C     The transformations are also applied to the rest of system
C     matrices
C
C          B <- Q' * B, C <- C * Z.
C
C     NUMERICAL ASPECTS
C
C     The algorithm is numerically backward stable and requires
C     0( L*L*N )  floating point operations.
C
C     CONTRIBUTOR
C
C     A. Varga, German Aerospace Center, DLR Oberpfaffenhofen.
C     March 1999. Based on the RASP routine RPDSSV.
C
C     REVISIONS
C
C     V. Sima, Research Institute for Informatics, Bucharest, July 1999,
C     May 2003, Jan. 2009, Apr. 2011, Sep. 2016, Nov. 2016, Feb. 2017.
C
C     KEYWORDS
C
C     Descriptor system, matrix algebra, matrix operations,
C     orthogonal transformation.
C
C     ******************************************************************
C
C     .. Parameters ..
      DOUBLE PRECISION   ONE, ZERO
      PARAMETER          ( ONE = 1.0D0, ZERO = 0.0D0 )
C     .. Scalar Arguments ..
      CHARACTER          COMPQ, COMPZ, JOBA
      INTEGER            INFO, L, LDA, LDB, LDC, LDE, LDQ, LDWORK,
     $                   LDZ, M, N, P, RANKE, RNKA22
      DOUBLE PRECISION   TOL
C     .. Array Arguments ..
      INTEGER            IWORK( * )
      DOUBLE PRECISION   A( LDA, * ), B( LDB, * ), C( LDC, * ),
     $                   DWORK( * ),  E( LDE, * ), Q( LDQ, * ),
     $                   Z( LDZ, * )
C     .. Local Scalars ..
      LOGICAL            ILQ, ILZ, LQUERY, REDA, REDTR, WITHB, WITHC
      INTEGER            I, ICOMPQ, ICOMPZ, IR1, IRE1, J, K, KW, LA22,
     $                   LH, LN, LWR, NA22, WRKOPT
      DOUBLE PRECISION   SVLMAX, TOLDEF
C     .. Local Arrays ..
      DOUBLE PRECISION   SVAL(3)
C     .. External Functions ..
      LOGICAL            LSAME
      DOUBLE PRECISION   DLAMCH, DLANGE
      EXTERNAL           DLAMCH, DLANGE, LSAME
C     .. External Subroutines ..
      EXTERNAL           DLASET, DORMQR, DORMRZ, DSWAP, DTZRZF, MB03OY,
     $                   XERBLA
C     .. Intrinsic Functions ..
      INTRINSIC          DBLE, INT, MAX, MIN
C
C     .. Executable Statements ..
C
C     Decode COMPQ.
C
      IF( LSAME( COMPQ, 'N' ) ) THEN
         ILQ = .FALSE.
         ICOMPQ = 1
      ELSE IF( LSAME( COMPQ, 'U' ) ) THEN
         ILQ = .TRUE.
         ICOMPQ = 2
      ELSE IF( LSAME( COMPQ, 'I' ) ) THEN
         ILQ = .TRUE.
         ICOMPQ = 3
      ELSE
         ICOMPQ = 0
      END IF
C
C     Decode COMPZ.
C
      IF( LSAME( COMPZ, 'N' ) ) THEN
         ILZ = .FALSE.
         ICOMPZ = 1
      ELSE IF( LSAME( COMPZ, 'U' ) ) THEN
         ILZ = .TRUE.
         ICOMPZ = 2
      ELSE IF( LSAME( COMPZ, 'I' ) ) THEN
         ILZ = .TRUE.
         ICOMPZ = 3
      ELSE
         ICOMPZ = 0
      END IF
      REDA  = LSAME( JOBA, 'R' )
      REDTR = LSAME( JOBA, 'T' )
      WITHB = M.GT.0
      WITHC = P.GT.0
      LQUERY = ( LDWORK.EQ.-1 )
C
C     Test the input parameters.
C
      LN = MIN( L, N )
      INFO = 0
      WRKOPT = MAX( 1, N+P, LN + MAX( 3*N-1, M, L ) )
      IF( ICOMPQ.LE.0 ) THEN
         INFO = -1
      ELSE IF( ICOMPZ.LE.0 ) THEN
         INFO = -2
      ELSE IF( .NOT.LSAME( JOBA, 'N' ) .AND. .NOT.REDA .AND.
     $         .NOT.REDTR ) THEN
         INFO = -3
      ELSE IF( L.LT.0 ) THEN
         INFO = -4
      ELSE IF( N.LT.0 ) THEN
         INFO = -5
      ELSE IF( M.LT.0 ) THEN
         INFO = -6
      ELSE IF( P.LT.0 ) THEN
         INFO = -7
      ELSE IF( LDA.LT.MAX( 1, L ) ) THEN
         INFO = -9
      ELSE IF( LDE.LT.MAX( 1, L ) ) THEN
         INFO = -11
      ELSE IF( LDB.LT.1 .OR. ( WITHB .AND. LDB.LT.L ) ) THEN
         INFO = -13
      ELSE IF( LDC.LT.MAX( 1, P ) ) THEN
         INFO = -15
      ELSE IF( ( ILQ .AND. LDQ.LT.L ) .OR. LDQ.LT.1 ) THEN
         INFO = -17
      ELSE IF( ( ILZ .AND. LDZ.LT.N ) .OR. LDZ.LT.1 ) THEN
         INFO = -19
      ELSE IF( TOL.GE.ONE ) THEN
         INFO = -22
      ELSE
         IF( LQUERY ) THEN
            CALL DORMQR( 'Left', 'Transpose', L, N, LN, E, LDE, DWORK,
     $                   A, LDA, DWORK, -1, INFO )
            WRKOPT = MAX( WRKOPT, LN + INT( DWORK(1) ) )
            IF( WITHB ) THEN
               CALL DORMQR( 'Left', 'Transpose', L, M, LN, E, LDE,
     $                      DWORK, B, LDB, DWORK, -1, INFO )
               WRKOPT = MAX( WRKOPT, LN + INT( DWORK(1) ) )
            END IF
            IF( ILQ ) THEN
               CALL DORMQR( 'Right', 'No Transpose', L, L, LN, E, LDE,
     $                      DWORK, Q, LDQ, DWORK, -1, INFO )
               WRKOPT = MAX( WRKOPT, LN + INT( DWORK(1) ) )
            END IF
            CALL DTZRZF( LN, N, E, LDE, DWORK, DWORK, -1, INFO )
            WRKOPT = MAX( WRKOPT, LN + INT( DWORK(1) ) )
            CALL DORMRZ( 'Right', 'Transpose', L, N, LN, N, E, LDE,
     $                   DWORK, A, LDA, DWORK, -1, INFO )
            WRKOPT = MAX( WRKOPT, N + INT( DWORK(1) ) )
            IF( WITHC ) THEN
               CALL DORMRZ( 'Right', 'Transpose', P, N, LN, N, E, LDE,
     $                      DWORK, C, LDC, DWORK, -1, INFO )
               WRKOPT = MAX( WRKOPT, N + INT( DWORK(1) ) )
            END IF
            IF( ILZ ) THEN
               CALL DORMRZ( 'Right', 'Transpose', N, N, LN, N, E, LDE,
     $                      DWORK, Z, LDZ, DWORK, -1, INFO )
               WRKOPT = MAX( WRKOPT, N + INT( DWORK(1) ) )
            END IF
         ELSE IF( LDWORK.LT.WRKOPT ) THEN
            INFO = -25
         END IF
      END IF
      IF( INFO.NE.0 ) THEN
         CALL XERBLA( 'TG01FD', -INFO )
         RETURN
      ELSE IF( LQUERY ) THEN
         DWORK(1) = WRKOPT
         RETURN
      END IF
C
C     Initialize Q and Z if necessary.
C
      IF( ICOMPQ.EQ.3 )
     $   CALL DLASET( 'Full', L, L, ZERO, ONE, Q, LDQ )
      IF( ICOMPZ.EQ.3 )
     $   CALL DLASET( 'Full', N, N, ZERO, ONE, Z, LDZ )
C
C     Quick return if possible.
C
      IF( L.EQ.0 .OR. N.EQ.0 ) THEN
         DWORK(1) = ONE
         RANKE = 0
         IF( REDA .OR. REDTR ) RNKA22 = 0
         RETURN
      END IF
C
      TOLDEF = TOL
      IF( TOLDEF.LE.ZERO ) THEN
C
C        Use the default tolerance for rank determination.
C
         TOLDEF = DBLE( L*N )*DLAMCH( 'EPSILON' )
      END IF
C
      SVLMAX = ZERO
C
C     Compute the rank-revealing QR decomposition of E,
C
C                        ( E11 E12 )
C           E * P = Qr * (         ) ,
C                        (  0  E22 )
C
C     and determine the rank of E using incremental condition
C     estimation.
C     Workspace: MIN(L,N) + 3*N - 1.
C
      LWR = LDWORK - LN
      KW  = LN + 1
C
      CALL MB03OY( L, N, E, LDE, TOLDEF, SVLMAX, RANKE, SVAL, IWORK,
     $             DWORK, DWORK(KW), INFO )
C
C     Apply transformation on the rest of matrices.
C
      IF( RANKE.GT.0 ) THEN
C
C        A <-- Qr' * A.
C        Workspace: need   MIN(L,N) + N;
C                   prefer MIN(L,N) + N*NB.
C
         CALL DORMQR( 'Left', 'Transpose', L, N, RANKE, E, LDE, DWORK,
     $                A, LDA, DWORK(KW), LWR, INFO )
         WRKOPT = MAX( WRKOPT, LN + INT( DWORK(KW) ) )
C
C        B <-- Qr' * B.
C        Workspace: need   MIN(L,N) + M;
C                   prefer MIN(L,N) + M*NB.
C
         IF( WITHB ) THEN
            CALL DORMQR( 'Left', 'Transpose', L, M, RANKE, E, LDE,
     $                   DWORK, B, LDB, DWORK(KW), LWR, INFO )
            WRKOPT = MAX( WRKOPT, LN + INT( DWORK(KW) ) )
         END IF
C
C        Q <-- Q * Qr.
C        Workspace: need   MIN(L,N) + L;
C                   prefer MIN(L,N) + L*NB.
C
         IF( ILQ ) THEN
            CALL DORMQR( 'Right', 'No Transpose', L, L, RANKE, E, LDE,
     $                   DWORK, Q, LDQ, DWORK(KW), LWR, INFO )
            WRKOPT = MAX( WRKOPT, LN + INT( DWORK(KW) ) )
         END IF
C
C        Set lower triangle of E to zero.
C
         IF( L.GE.2 )
     $      CALL DLASET( 'Lower', L-1, RANKE, ZERO, ZERO, E(2,1), LDE )
C
C        Compute A*P, C*P and Z*P by forward permuting the columns of
C        A, C and Z based on information in IWORK.
C
         DO 10 J = 1, N
            IWORK(J) = -IWORK(J)
   10    CONTINUE
         DO 30 I = 1, N
            IF( IWORK(I).LT.0 ) THEN
               J = I
               IWORK(J) = -IWORK(J)
   20          CONTINUE
               K = IWORK(J)
               IF( IWORK(K).LT.0 ) THEN
                  CALL DSWAP( L, A(1,J), 1, A(1,K), 1 )
                  IF( WITHC )
     $               CALL DSWAP( P, C(1,J), 1, C(1,K), 1 )
                  IF( ILZ )
     $               CALL DSWAP( N, Z(1,J), 1, Z(1,K), 1 )
                  IWORK(K) = -IWORK(K)
                  J = K
                  GO TO 20
               END IF
            END IF
   30    CONTINUE
C
C        Determine an orthogonal matrix Y such that
C
C           ( E11 E12 ) = ( Er  0 ) * Y .
C
C        Compute E <-- E*Y', A <-- A*Y', C <-- C*Y', Z <-- Z*Y'.
C
         IF( RANKE.LT.N ) THEN
C
C           Workspace: need   2*N;
C                      prefer N + N*NB.
C
            KW = RANKE + 1
            CALL DTZRZF( RANKE, N, E, LDE, DWORK, DWORK(KW),
     $                   LDWORK-KW+1, INFO )
            WRKOPT = MAX( WRKOPT, INT( DWORK(KW) ) + KW - 1 )
C
C           Workspace: need   N + MAX(L,P,N);
C                      prefer N + MAX(L,P,N)*NB.
C
            LH = N - RANKE
            CALL DORMRZ( 'Right', 'Transpose', L, N, RANKE, LH, E, LDE,
     $                   DWORK, A, LDA, DWORK(KW), LDWORK-KW+1, INFO )
            WRKOPT = MAX( WRKOPT, INT( DWORK(KW) ) + KW - 1 )
            IF( WITHC ) THEN
               CALL DORMRZ( 'Right', 'Transpose', P, N, RANKE, LH, E,
     $                      LDE, DWORK, C, LDC, DWORK(KW), LDWORK-KW+1,
     $                      INFO )
               WRKOPT = MAX( WRKOPT, INT( DWORK(KW) ) + KW - 1 )
            END IF
            IF( ILZ ) THEN
               CALL DORMRZ( 'Right', 'Transpose', N, N, RANKE, LH, E,
     $                      LDE, DWORK, Z, LDZ, DWORK(KW), LDWORK-KW+1,
     $                      INFO )
               WRKOPT = MAX( WRKOPT, INT( DWORK(KW) ) + KW - 1 )
            END IF
C
C           Set E12 and E22 to zero.
C
            CALL DLASET( 'Full', L, LH, ZERO, ZERO, E(1,KW), LDE )
         END IF
      ELSE
         CALL DLASET( 'Full', L, N, ZERO, ZERO, E, LDE )
      END IF
C
C     Reduce A22 if necessary.
C
      IF( REDA .OR. REDTR ) THEN
         LA22 = L - RANKE
         NA22 = N - RANKE
         IF( MIN( LA22, NA22 ).EQ.0 ) THEN
            RNKA22 = 0
         ELSE
C
C           Compute the rank-revealing QR decomposition of A22,
C
C                              ( R11 R12 )
C              A22 * P2 = Q2 * (         ) ,
C                              (  0  R22 )
C
C           and determine the rank of A22 using incremental
C           condition estimation.
C           Workspace: MIN(L,N) + 3*N - 1.
C
C           Set the estimate of maximum singular value of A to detect
C           a negligible A matrix.
C
            SVLMAX = DLANGE( 'Frobenius', L, N, A, LDA, DWORK )
            IR1 = RANKE + 1
            KW  = MIN( LA22, NA22 ) + 1
            CALL MB03OY( LA22, NA22, A(IR1,IR1), LDA, TOLDEF,
     $                   SVLMAX, RNKA22, SVAL, IWORK, DWORK,
     $                   DWORK(KW), INFO )
C
C           Apply transformation on the rest of matrices.
C
            IF( RNKA22.GT.0 ) THEN
C
C              A <-- diag(I, Q2') * A
C              Workspace: need   MIN(L,N) + N;
C                         prefer MIN(L,N) + N*NB.
C
               CALL DORMQR( 'Left', 'Transpose', LA22, RANKE, RNKA22,
     $                      A(IR1,IR1), LDA, DWORK, A(IR1,1), LDA,
     $                      DWORK(KW), LWR, INFO )
C
C              B <-- diag(I, Q2') * B
C              Workspace: need   MIN(L,N) + M;
C                         prefer MIN(L,N) + M*NB.
C
               IF ( WITHB )
     $            CALL DORMQR( 'Left', 'Transpose', LA22, M, RNKA22,
     $                         A(IR1,IR1), LDA, DWORK, B(IR1,1), LDB,
     $                         DWORK(KW), LWR, INFO )
C
C              Q <-- Q * diag(I, Q2)
C              Workspace: need   MIN(L,N) + L;
C                         prefer MIN(L,N) + L*NB.
C
               IF( ILQ )
     $            CALL DORMQR( 'Right', 'No transpose', L, LA22, RNKA22,
     $                         A(IR1,IR1), LDA, DWORK, Q(1,IR1), LDQ,
     $                         DWORK(KW), LWR, INFO )
C
C              Set lower triangle of A22 to zero.
C
               IF( LA22.GE.2 )
     $            CALL DLASET( 'Lower', LA22-1, RNKA22, ZERO, ZERO,
     $                         A(IR1+1,IR1), LDA )
C
C              Compute A*diag(I,P2), C*diag(I,P2) and Z*diag(I,P2)
C              by forward permuting the columns of A, C and Z based
C              on information in IWORK.
C
               DO 40 J = 1, NA22
                  IWORK(J) = -IWORK(J)
   40          CONTINUE
               DO 60 I = 1, NA22
                  IF( IWORK(I).LT.0 ) THEN
                     J = I
                     IWORK(J) = -IWORK(J)
   50                CONTINUE
                     K = IWORK(J)
                     IF( IWORK(K).LT.0 ) THEN
                        CALL DSWAP( RANKE, A(1,RANKE+J), 1,
     $                              A(1,RANKE+K), 1 )
                        IF( WITHC )
     $                     CALL DSWAP( P, C(1,RANKE+J), 1,
     $                                 C(1,RANKE+K), 1 )
                        IF( ILZ )
     $                     CALL DSWAP( N, Z(1,RANKE+J), 1,
     $                                 Z(1,RANKE+K), 1 )
                        IWORK(K) = -IWORK(K)
                        J = K
                        GO TO 50
                     END IF
                  END IF
   60          CONTINUE
C
               IF( REDA .AND. RNKA22.LT.NA22 ) THEN
C
C                 Determine an orthogonal matrix Y2 such that
C
C                 ( R11 R12 ) = ( Ar  0 ) * Y2 .
C
C                 Compute A <-- A*diag(I, Y2'), C <-- C*diag(I, Y2'),
C                         Z <-- Z*diag(I, Y2').
C                 Workspace: need   2*N.
C                            prefer N + N*NB.
C
                  KW = RANKE + 1
                  CALL DTZRZF( RNKA22, NA22, A(IR1,IR1), LDA, DWORK,
     $                         DWORK(KW), LDWORK-KW+1, INFO )
                  WRKOPT = MAX( WRKOPT, INT( DWORK(KW) ) + KW - 1 )
C
C                 Workspace: need   N + MAX(P,N);
C                            prefer N + MAX(P,N)*NB.
C
                  LH = NA22 - RNKA22
                  CALL DORMRZ( 'Right', 'Transpose', RANKE, NA22,
     $                         RNKA22, LH, A(IR1,IR1), LDA, DWORK,
     $                         A(1,IR1), LDA, DWORK(KW), LDWORK-KW+1,
     $                         INFO )
                  IF( WITHC ) THEN
                     CALL DORMRZ( 'Right', 'Transpose', P, NA22, RNKA22,
     $                            LH, A(IR1,IR1), LDA, DWORK, C(1,IR1),
     $                            LDC, DWORK(KW), LDWORK-KW+1, INFO )
                     WRKOPT = MAX( WRKOPT, INT( DWORK(KW) ) + KW - 1 )
                  END IF
                  IF( ILZ ) THEN
                     CALL DORMRZ( 'Right', 'Transpose', N, NA22, RNKA22,
     $                            LH, A(IR1,IR1), LDA, DWORK, Z(1,IR1),
     $                            LDZ, DWORK(KW), LDWORK-KW+1, INFO )
                     WRKOPT = MAX( WRKOPT, INT( DWORK(KW) ) + KW - 1 )
                  END IF
                  IRE1 = RANKE + RNKA22 + 1
C
C                 Set R12 and R22 to zero.
C
                  CALL DLASET( 'Full', LA22, LH, ZERO, ZERO,
     $                         A(IR1,IRE1), LDA )
               END IF
            ELSE
               CALL DLASET( 'Full', LA22, NA22, ZERO, ZERO,
     $                      A(IR1,IR1), LDA)
            END IF
         END IF
      END IF
C
      DWORK(1) = WRKOPT
C
      RETURN
C *** Last line of TG01FD ***
      END

ext: f

max_stars_repo_name: bnavigator/SLICOT-Reference

max_issues_repo_name: RJHKnight/slicotr

max_forks_repo_name: RJHKnight/slicotr


text:
!!  Copyright (C)  Stichting Deltares, 2012-2016.
!!
!!  This program is free software: you can redistribute it and/or modify
!!  it under the terms of the GNU General Public License version 3,
!!  as published by the Free Software Foundation.
!!
!!  This program is distributed in the hope that it will be useful,
!!  but WITHOUT ANY WARRANTY; without even the implied warranty of
!!  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
!!  GNU General Public License for more details.
!!
!!  You should have received a copy of the GNU General Public License
!!  along with this program. If not, see <http://www.gnu.org/licenses/>.
!!
!!  contact: delft3d.support@deltares.nl
!!  Stichting Deltares
!!  P.O. Box 177
!!  2600 MH Delft, The Netherlands
!!
!!  All indications and logos of, and references to registered trademarks
!!  of Stichting Deltares remain the property of Stichting Deltares. All
!!  rights reserved.

!    Date:       22 Oct 1992
!    Time:       21:28
!    Program:    MATINV.FOR
!    Version:    6.00.00
!    Programmer: Nicolaas M de Rooij
!
      subroutine matinv(a, n, b, m, d, w, ising)
!           4-69
!           cholesky's diagonal method -- for symmetric matrices only
!           matrix inversion and / or simultaneous equation solution
!             m =  0    invert only and store inverse in a
!             m = -1    solve simultaneous equations.  solution in b
!             m =  1    combination of m = 0 and m = -1.
!             m = -2    solve simultaneous equations with already
!                       calculated upper triangular and diagonal
!                       inverse stored in a, and solution in b
!           dimension of the arrays should be equal to the value of
!        maxmd in common labeled /intgr/ of other routines in deck.
      implicit real*8 (a-h,o-z)
      dimension a(60,60), b(60), d(60), w(60)
      ising = 0
      if(m.eq.-2)go to 35
!        find upper triangular matrix and diagonal matrix
      do 30 i=1,n
         if (a(i,i).eq.0.) go to 170
         d(i) = a(i,i)
         a(i,i) = 1.
         if (i.eq.n) go to 30
         k = i + 1
         do 20 j=k,n
            do 10 l=j,n
               a(j,l) = a(j,l) - a(i,j)*a(i,l)/d(i)
   10       continue
            a(i,j) = a(i,j)/d(i)
   20    continue
   30 continue
      if (m.eq.0) go to 80
!        solve simultaneous equations
   35 w(1) = b(1)
      do 50 i=2,n
         s = 0.
         l2 = i - 1
         do 40 l=1,l2
            s = s + a(l,i)*w(l)
   40    continue
         w(i) = b(i) - s
   50 continue
      b(n) = w(n)/d(n)
      do 70 ii=2,n
         i = n + 1 - ii
         s = 0.
         l1 = i + 1
         do 60 l=l1,n
            s = s + a(i,l)*b(l)
   60    continue
         b(i) = w(i)/d(i) - s
   70 continue
   80 if (m.lt.0) return
!        invert upper triangular matrix  to get  inverse
      do 110 ii=2,n
         i = n + 1 - ii
         jj = ii - 1
         j3 = i + 1
         do 100 j2=1,jj
            j = n + 1 - j2
            s = 0.
            do 90 k=j3,j
               s = s + a(i,k)*a(k,j)
   90       continue
            a(i,j) = -s
  100    continue
  110 continue
!        multiply inverse * reciprocal diagonal * inverse transposed
      do 140 i=1,n
         do 130 j=i,n
            s = 0.
            do 120 k=j,n
               s = s + a(i,k)*a(j,k)/d(k)
  120       continue
            a(i,j) = s
  130    continue
  140 continue
!        re-symmetrize
      do 160 j=2,n
         j2 = j - 1
         do 150 i=1,j2
            a(j,i) = a(i,j)
  150    continue
  160 continue
      return
  170 ising = i
      return
      end

max_stars_repo_name: liujiamingustc/phd

max_issues_repo_name: liujiamingustc/phd

max_forks_repo_name: liujiamingustc/phd


     subroutine fort_set_homog_bcs(lo,hi,domlo,domhi, &
                                   phi,phi_l1,phi_l2,phi_l3,phi_h1,phi_h2,phi_h3,dx);

     use amrex_fort_module, only : rt => amrex_real
     implicit none

     integer         ,intent(in   ) :: lo(3),hi(3),domlo(3),domhi(3)
     integer         ,intent(in   ) :: phi_l1,phi_l2,phi_l3,phi_h1,phi_h2,phi_h3
     real(rt),intent(  out) :: phi(phi_l1:phi_h1,phi_l2:phi_h2,phi_l3:phi_h3)
     real(rt),intent(in   ) :: dx(3)

     phi = 0.d0

     end subroutine fort_set_homog_bcs


text:
      SUBROUTINE DAALL(IC,L,CCC,NO,NV)
      COMMON /main1/ eps
      END

max_issues_repo_name: sujankh/rose-matlab

max_forks_repo_name: sujankh/rose-matlab


text:
+UVM_NO_RELNOTES
+incdir+${PROJ_HOME}
+incdir+${PROJ_HOME}/sv

-dpi
-DIRUN
-cpost "${PROJ_HOME}/sv/svaunit_vpi_api.cpp" -end
-uvm
-sv
-assert
${PROJ_HOME}/sv/svaunit_pkg.sv
${PROJ_HOME}/sv/svaunit_vpi_interface.sv
-f ${file_name}
-access rwc
-linedebug
-assert
-coverage all
-covoverwrite
-timescale 1ns/1ps

-top ${top_name}


text:
module pf_my_level
  use encap
  use pf_mod_imex_sweeper

  !>  extend the generic level type by defining transfer operators
  type, extends(pf_user_level_t) :: my_level_t
   integer :: nlevels
   contains
     procedure :: restrict
     procedure :: interpolate
  end type my_level_t

  interface

     subroutine HypreRestrict(y_f, y_c, f_level, c_level) bind(c, name="HypreRestrict")
           use iso_c_binding
           type(c_ptr), value :: y_c, y_f
           integer, value :: f_level, c_level
     end subroutine HypreRestrict

     subroutine HypreProlong(y_f, y_c, f_level, c_level) bind(c, name="HypreProlong")
           use iso_c_binding
           type(c_ptr), value :: y_c, y_f
           integer, value :: f_level, c_level
     end subroutine HypreProlong

  end interface

contains

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!>  These are the transfer functions that must be  provided for the level
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  !>  Interpolate from coarse level to fine
  subroutine interpolate(this, f_lev, c_lev, f_vec, c_vec, t, flags)
    class(my_level_t), intent(inout) :: this
    real(pfdp),        intent(in   ) :: t
    class(pf_level_t), intent(inout) :: f_lev, c_lev  !  fine and coarse levels
    class(pf_encap_t), intent(inout) :: f_vec, c_vec  !  fine and coarse vectors
    integer, intent(in), optional :: flags
    integer :: f_level_index, c_level_index
    integer :: restrict_flag

    class(hypre_vector_encap), pointer :: y_f, y_c

    !>  Cast the abstract encap as my data type
    y_f => cast_as_hypre_vector(f_vec)
    y_c => cast_as_hypre_vector(c_vec)

    f_level_index = f_lev%index
    c_level_index = c_lev%index

    !> Here we use the identity map
    call HypreProlong(y_f%c_hypre_vector_ptr, &
                      y_c%c_hypre_vector_ptr, &
                      f_level_index, c_level_index)
  end subroutine interpolate

  !>  Restrict from fine level to coarse
  subroutine restrict(this, f_lev, c_lev, f_vec, c_vec, t, flags)
    class(my_level_t), intent(inout) :: this
    class(pf_level_t), intent(inout) :: f_lev, c_lev  !  fine and coarse levels
    class(pf_encap_t),  intent(inout) :: f_vec, c_vec  !  fine and coarse vectors
    real(pfdp), intent(in) :: t             !<  time of solution
    integer, intent(in), optional :: flags
    integer :: f_level_index, c_level_index
    integer :: restrict_flag

    class(hypre_vector_encap), pointer :: y_f, y_c

    !>  Cast the abstract encap as my data type
    y_f => cast_as_hypre_vector(f_vec)
    y_c => cast_as_hypre_vector(c_vec)

    !f_level_index = abs(f_lev%index - this%nlevels)
    !c_level_index = abs(c_lev%index - this%nlevels)

    f_level_index = f_lev%index
    c_level_index = c_lev%index

    !> Here we use the identity map
    call HypreRestrict(y_f%c_hypre_vector_ptr, &
                       y_c%c_hypre_vector_ptr, &
                       f_level_index, c_level_index)
  end subroutine restrict

  function cast_as_my_level_t(level_polymorph) result(my_level_t_obj)
    class(pf_user_level_t), intent(in), target :: level_polymorph
    type(my_level_t), pointer :: my_level_t_obj

    select type(level_polymorph)
    type is (my_level_t)
       my_level_t_obj => level_polymorph
    end select
  end function cast_as_my_level_t

end module pf_my_level

max_stars_repo_name: schreiberx/LibPFASST

max_issues_repo_name: schreiberx/LibPFASST

max_forks_repo_name: schreiberx/LibPFASST


text:
module Mod_supm_TemporalDerivatives
   use Mod_supm_BaseElmope
   use Mod_ThreeFieldElement
   use Mod_LogarithmicElement
   implicit none
   private
   public SetPointersTemporalDerivatives
   integer(ip), allocatable :: kfl_IsSet

contains

   subroutine SetPointersTemporalDerivatives(itask)
      integer(ip) :: itask
      !procedure() :: NULL()

      select case (itask)
         case(0)
            allocate(kfl_IsSet)
            call a%Memor%allocObj(0,'kfl_IsSet','InitProcedurePointer',1)
            kfl_IsSet = -1

         case(1)
            if (kfl_IsSet == -1) then
               kfl_IsSet = 1
                if (a%LogFormulation==0) then
                  ProcPointer%TemporalDerivatives => sup_temporalpart
                else if (a%LogFormulation==1) then
                  ProcPointer%TemporalDerivatives => sup_temporalpart_LCR
                end if
            endif

         case(100)
            deallocate(kfl_IsSet)
            call a%Memor%deallocObj(0,'kfl_IsSet','InitProcedurePointer',1)

       end select
   end subroutine

   subroutine sup_temporalpart
      implicit none
      call sup_TimeIntegrationToElext(e,auxtens,integrator,auxVE,a%dtinv,gpsig,elextS)
   end subroutine

   subroutine sup_temporalpart_LCR
      implicit none
      call sup_LCR2_TimeIntegrationToElext(e,auxtens,nsteps,integrator,a%dtinv,lambda,lambda0,gpsig,elextS)
   end subroutine

end module

max_stars_repo_name: ciaid-colombia/InsFEM

max_issues_repo_name: ciaid-colombia/InsFEM

max_forks_repo_name: ciaid-colombia/InsFEM


text:
Describe Users/il here.

Nice work. jw

20080401 05:53:11 nbsp dumbass Users/rickrolled
 Thanks. I appreciate open communication. Users/il

20080508 09:28:12 nbsp Hey, I dont know about publicity this year, but in the past everyone knew about the coop party during WEF..?  I think there were even posters up around advertising the bands?  Maybe, if its a problem, the reference to the coop could be removed and the pictures / rest of text remain? (Every year theres a big party accompanying WEF..?) Users/PhilipNeustrom

20080607 18:03:30 nbsp lol, I thought I was the one that accidentally erased the events board cause I tried to add an event while the board was gone and saw what happened. Users/SunjeetBaadkar

20080610 15:37:14 nbsp Its possible; I do remember someone asking me about Naked Lunch.  Say hi next time youre in. Users/SunjeetBaadkar

20090208 18:30:52 nbsp Hey, thanks for updating Board games!  Do you have any more details (time, do people actually show up, etc.) about Go at Delta?  Ive been looking for people to play, but havent had much luck. Users/JonathanLawton

20090407 11:43:40 nbsp So, I never made it to Go night at Delta last quarter.  Are you still meeting there at 7PM on Wednesdays this quarter? Users/JonathanLawton

20091115 20:23:21 nbsp I know its not in Davis, but I have bought lime leaves at SF Market in Sac. Users/JonathanLawton

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


text:
c
c  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
c  .                                                             .
c  .                  copyright (c) 1998 by UCAR                 .
c  .                                                             .
c  .       University Corporation for Atmospheric Research       .
c  .                                                             .
c  .                      all rights reserved                    .
c  .                                                             .
c  .                                                             .
c  .                         SPHEREPACK3.0                       .
c  .                                                             .
c  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
c
c
c
c ... file shagc.f
c
c     this file contains code and documentation for subroutines
c     shagc and shagci
c
c ... files which must be loaded with shagc.f
c
c     sphcom.f, hrfft.f, gaqd.f
c
c
c     subroutine shagc(nlat,nlon,isym,nt,g,idg,jdg,a,b,mdab,ndab,
c    +                 wshagc,lshagc,work,lwork,ierror)
c
c     subroutine shagc performs the spherical harmonic analysis
c     on the array g and stores the result in the arrays a and b.
c     the analysis is performed on a gaussian grid in colatitude
c     and an equally spaced grid in longitude.  the associated
c     legendre functions are recomputed rather than stored as they
c     are in subroutine shags.  the analysis is described below
c     at output parameters a,b.
c
c     input parameters
c
c     nlat   the number of points in the gaussian colatitude grid on the
c            full sphere. these lie in the interval (0,pi) and are compu
c            in radians in theta(1),...,theta(nlat) by subroutine gaqd.
c            if nlat is odd the equator will be included as the grid poi
c            theta((nlat+1)/2).  if nlat is even the equator will be
c            excluded as a grid point and will lie half way between
c            theta(nlat/2) and theta(nlat/2+1). nlat must be at least 3.
c            note: on the half sphere, the number of grid points in the
c            colatitudinal direction is nlat/2 if nlat is even or
c            (nlat+1)/2 if nlat is odd.
c
c     nlon   the number of distinct londitude points.  nlon determines
c            the grid increment in longitude as 2*pi/nlon. for example
c            nlon = 72 for a five degree grid. nlon must be greater
c            than or equal to 4. the efficiency of the computation is
c            improved when nlon is a product of small prime numbers.
c
c     isym   = 0  no symmetries exist about the equator. the analysis
c                 is performed on the entire sphere.  i.e. on the
c                 array g(i,j) for i=1,...,nlat and j=1,...,nlon.
c                 (see description of g below)
c
c            = 1  g is antisymmetric about the equator. the analysis
c                 is performed on the northern hemisphere only.  i.e.
c                 if nlat is odd the analysis is performed on the
c                 array g(i,j) for i=1,...,(nlat+1)/2 and j=1,...,nlon.
c                 if nlat is even the analysis is performed on the
c                 array g(i,j) for i=1,...,nlat/2 and j=1,...,nlon.
c
c
c            = 2  g is symmetric about the equator. the analysis is
c                 performed on the northern hemisphere only.  i.e.
c                 if nlat is odd the analysis is performed on the
c                 array g(i,j) for i=1,...,(nlat+1)/2 and j=1,...,nlon.
c                 if nlat is even the analysis is performed on the
c                 array g(i,j) for i=1,...,nlat/2 and j=1,...,nlon.
c
c     nt     the number of analyses.  in the program that calls shagc,
c            the arrays g,a and b can be three dimensional in which
c            case multiple analyses will be performed.  the third
c            index is the analysis index which assumes the values
c            k=1,...,nt.  for a single analysis set nt=1. the
c            discription of the remaining parameters is simplified
c            by assuming that nt=1 or that the arrays g,a and b
c            have only two dimensions.
c
c     g      a two or three dimensional array (see input parameter
c            nt) that contains the discrete function to be analyzed.
c            g(i,j) contains the value of the function at the gaussian
c            point theta(i) and longitude point phi(j) = (j-1)*2*pi/nlon
c            the index ranges are defined above at the input parameter
c            isym.
c
c     idg    the first dimension of the array g as it appears in the
c            program that calls shagc. if isym equals zero then idg
c            must be at least nlat.  if isym is nonzero then idg must
c            be at least nlat/2 if nlat is even or at least (nlat+1)/2
c            if nlat is odd.
c
c     jdg    the second dimension of the array g as it appears in the
c            program that calls shagc. jdg must be at least nlon.
c
c     mdab   the first dimension of the arrays a and b as it appears
c            in the program that calls shagc. mdab must be at least
c            min0((nlon+2)/2,nlat) if nlon is even or at least
c            min0((nlon+1)/2,nlat) if nlon is odd
c
c     ndab   the second dimension of the arrays a and b as it appears
c            in the program that calls shaec. ndab must be at least nlat
c
c     wshagc an array which must be initialized by subroutine shagci.
c            once initialized, wshagc can be used repeatedly by shagc.
c            as long as nlat and nlon remain unchanged.  wshagc must
c            not be altered between calls of shagc.
c
c     lshagc the dimension of the array wshagc as it appears in the
c            program that calls shagc. define
c
c               l1 = min0(nlat,(nlon+2)/2) if nlon is even or
c               l1 = min0(nlat,(nlon+1)/2) if nlon is odd
c
c            and
c
c               l2 = nlat/2        if nlat is even or
c               l2 = (nlat+1)/2    if nlat is odd
c
c            then lshagc must be at least
c
c                  nlat*(2*l2+3*l1-2)+3*l1*(1-l1)/2+nlon+15
c
c
c     work   a work array that does not have to be saved.
c
c     lwork  the dimension of the array work as it appears in the
c            program that calls shagc. define
c
c               l1 = min0(nlat,(nlon+2)/2) if nlon is even or
c               l1 = min0(nlat,(nlon+1)/2) if nlon is odd
c
c            and
c
c               l2 = nlat/2        if nlat is even or
c               l2 = (nlat+1)/2    if nlat is odd
c
c            if isym is zero then lwork must be at least
c
c                      nlat*(nlon*nt+max0(3*l2,nlon))
c
c            if isym is not zero then lwork must be at least
c
c                      l2*(nlon*nt+max0(3*nlat,nlon))
c
c     **************************************************************
c
c     output parameters
c
c     a,b    both a,b are two or three dimensional arrays (see input
c            parameter nt) that contain the spherical harmonic
c            coefficients in the representation of g(i,j) given in the
c            discription of subroutine shagc. for isym=0, a(m,n) and
c            b(m,n) are given by the equations listed below. symmetric
c            versions are used when isym is greater than zero.
c
c     definitions
c
c     1. the normalized associated legendre functions
c
c     pbar(m,n,theta) = sqrt((2*n+1)*factorial(n-m)/(2*factorial(n+m)))
c                       *sin(theta)**m/(2**n*factorial(n)) times the
c                       (n+m)th derivative of (x**2-1)**n with respect
c                       to x=cos(theta).
c
c     2. the fourier transform of g(i,j).
c
c     c(m,i)          = 2/nlon times the sum from j=1 to j=nlon of
c                       g(i,j)*cos((m-1)*(j-1)*2*pi/nlon)
c                       (the first and last terms in this sum
c                       are divided by 2)
c
c     s(m,i)          = 2/nlon times the sum from j=2 to j=nlon of
c                       g(i,j)*sin((m-1)*(j-1)*2*pi/nlon)
c
c
c     3. the gaussian points and weights on the sphere
c        (computed by subroutine gaqd).
c
c        theta(1),...,theta(nlat) (gaussian pts in radians)
c        wts(1),...,wts(nlat) (corresponding gaussian weights)
c
c     4. the maximum (plus one) longitudinal wave number
c
c            mmax = min0(nlat,(nlon+2)/2) if nlon is even or
c            mmax = min0(nlat,(nlon+1)/2) if nlon is odd.
c
c
c     then for m=0,...,mmax-1 and n=m,...,nlat-1 the arrays a,b
c     are given by
c
c     a(m+1,n+1)     =  the sum from i=1 to i=nlat of
c                       c(m+1,i)*wts(i)*pbar(m,n,theta(i))
c
c     b(m+1,n+1)      = the sum from i=1 to nlat of
c                       s(m+1,i)*wts(i)*pbar(m,n,theta(i))
c
c     ierror = 0  no errors
c            = 1  error in the specification of nlat
c            = 2  error in the specification of nlon
c            = 3  error in the specification of isym
c            = 4  error in the specification of nt
c            = 5  error in the specification of idg
c            = 6  error in the specification of jdg
c            = 7  error in the specification of mdab
c            = 8  error in the specification of ndab
c            = 9  error in the specification of lshagc
c            = 10 error in the specification of lwork
c
c
c ****************************************************************
c
c     subroutine shagci(nlat,nlon,wshagc,lshagc,dwork,ldwork,ierror)
c
c     subroutine shagci initializes the array wshagc which can then
c     be used repeatedly by subroutines shagc. it precomputes
c     and stores in wshagc quantities such as gaussian weights,
c     legendre polynomial coefficients, and fft trigonometric tables.
c
c     input parameters
c
c     nlat   the number of points in the gaussian colatitude grid on the
c            full sphere. these lie in the interval (0,pi) and are compu
c            in radians in theta(1),...,theta(nlat) by subroutine gaqd.
c            if nlat is odd the equator will be included as the grid poi
c            theta((nlat+1)/2).  if nlat is even the equator will be
c            excluded as a grid point and will lie half way between
c            theta(nlat/2) and theta(nlat/2+1). nlat must be at least 3.
c            note: on the half sphere, the number of grid points in the
c            colatitudinal direction is nlat/2 if nlat is even or
c            (nlat+1)/2 if nlat is odd.
c
c     nlon   the number of distinct londitude points.  nlon determines
c            the grid increment in longitude as 2*pi/nlon. for example
c            nlon = 72 for a five degree grid. nlon must be greater
c            than or equal to 4. the efficiency of the computation is
c            improved when nlon is a product of small prime numbers.
c
c     wshagc an array which must be initialized by subroutine shagci.
c            once initialized, wshagc can be used repeatedly by shagc
c            as long as nlat and nlon remain unchanged.  wshagc must
c            not be altered between calls of shagc.
c
c     lshagc the dimension of the array wshagc as it appears in the
c            program that calls shagc. define
c
c               l1 = min0(nlat,(nlon+2)/2) if nlon is even or
c               l1 = min0(nlat,(nlon+1)/2) if nlon is odd
c
c            and
c
c               l2 = nlat/2        if nlat is even or
c               l2 = (nlat+1)/2    if nlat is odd
c
c            then lshagc must be at least
c
c                  nlat*(2*l2+3*l1-2)+3*l1*(1-l1)/2+nlon+15
c
c     dwork   a double precision work array that does not have to be saved.
c
c     ldwork  the dimension of the array dwork as it appears in the
c            program that calls shagci. ldwork must be at least
c
c                nlat*(nlat+4)
c
c     output parameter
c
c     wshagc an array which must be initialized before calling shagc or
c            once initialized, wshagc can be used repeatedly by shagc or
c            as long as nlat and nlon remain unchanged.  wshagc must not
c            altered between calls of shagc.
c
c     ierror = 0  no errors
c            = 1  error in the specification of nlat
c            = 2  error in the specification of nlon
c            = 3  error in the specification of lshagc
c            = 4  error in the specification of ldwork
c            = 5  failure in gaqd to compute gaussian points
c                 (due to failure in eigenvalue routine)
c
c
c ****************************************************************
      subroutine shagc(nlat,nlon,isym,nt,g,idg,jdg,a,b,mdab,ndab,
     1                    wshagc,lshagc,work,lwork,ierror)
c     subroutine shagc performs the spherical harmonic analysis on
c     a gaussian grid on the array(s) in g and returns the coefficients
c     in array(s) a,b. the necessary legendre polynomials are computed
c     as needed in this version.
c
      dimension g(idg,jdg,1),a(mdab,ndab,1),b(mdab,ndab,1),
     1          wshagc(lshagc),work(lwork)
c     check input parameters
      ierror = 1
      if (nlat.lt.3) return
      ierror = 2
      if (nlon.lt.4) return
      ierror = 3
      if (isym.lt.0 .or.isym.gt.2) return
      ierror = 4
      if (nt.lt.1) return
c     set upper limit on m for spherical harmonic basis
      l = min0((nlon+2)/2,nlat)
c     set gaussian point nearest equator pointer
      late = (nlat+mod(nlat,2))/2
c     set number of grid points for analysis/synthesis
      lat = nlat
      if (isym.ne.0) lat = late
      ierror = 5
      if (idg.lt.lat) return
      ierror = 6
      if (jdg.lt.nlon) return
      ierror = 7
      if(mdab .lt. l) return
      ierror = 8
      if(ndab .lt. nlat) return
      l1 = l
      l2 = late
      ierror = 9
c     check permanent work space length
      if (lshagc .lt. nlat*(2*l2+3*l1-2)+3*l1*(1-l1)/2+nlon+15)return
      ierror = 10
c     check temporary work space length
      if (isym.eq.0) then
      if(lwork.lt.nlat*(nlon*nt+max0(3*l2,nlon))) return
      else
c     isym.ne.0
      if(lwork.lt.l2*(nlon*nt+max0(3*nlat,nlon))) return
      end if
      ierror = 0
c     starting address for gaussian wts in shigc and fft values
      iwts = 1
      ifft = nlat+2*nlat*late+3*(l*(l-1)/2+(nlat-l)*(l-1))+1
c     set pointers for internal storage of g and legendre polys
      ipmn = lat*nlon*nt+1
      call shagc1(nlat,nlon,l,lat,isym,g,idg,jdg,nt,a,b,mdab,ndab,
     1wshagc,wshagc(iwts),wshagc(ifft),late,work(ipmn),work)
      return
      end
      subroutine shagc1(nlat,nlon,l,lat,mode,gs,idg,jdg,nt,a,b,mdab,
     1                  ndab,w,wts,wfft,late,pmn,g)
      dimension gs(idg,jdg,nt),a(mdab,ndab,nt),
     1          b(mdab,ndab,nt),g(lat,nlon,nt)
      dimension w(1),wts(nlat),wfft(1),pmn(nlat,late,3)
c     set gs array internally in shagc1
      do 100 k=1,nt
      do 100 j=1,nlon
      do 100 i=1,lat
      g(i,j,k) = gs(i,j,k)
  100 continue
c     do fourier transform
      do 101 k=1,nt
      call hrfftf(lat,nlon,g(1,1,k),lat,wfft,pmn)
  101 continue
c     scale result
      sfn = 2.0/float(nlon)
      do 102 k=1,nt
      do 102 j=1,nlon
      do 102 i=1,lat
      g(i,j,k) = sfn*g(i,j,k)
  102 continue
c     compute using gaussian quadrature
c     a(n,m) = s (ga(theta,m)*pnm(theta)*sin(theta)*dtheta)
c     b(n,m) = s (gb(theta,m)*pnm(theta)*sin(theta)*dtheta)
c     here ga,gb are the cos(phi),sin(phi) coefficients of
c     the fourier expansion of g(theta,phi) in phi.  as a result
c     of the above fourier transform they are stored in array
c     g as follows:
c     for each theta(i) and k= l-1
c     ga(0),ga(1),gb(1),ga(2),gb(2),...,ga(k-1),gb(k-1),ga(k)
c     correspond to (in the case nlon=l+l-2)
c     g(i,1),g(i,2),g(i,3),g(i,4),g(i,5),...,g(i,2l-4),g(i,2l-3),g(i,2l-
c     initialize coefficients to zero
      do 103 k=1,nt
      do 103 np1=1,nlat
      do 103 mp1=1,l
      a(mp1,np1,k) = 0.0
      b(mp1,np1,k) = 0.0
  103 continue
c     set m+1 limit on b(m+1) calculation
      lm1 = l
      if (nlon .eq. l+l-2) lm1 = l-1
      if (mode.eq.0) then
c     for full sphere (mode=0) and even/odd reduction:
c     overwrite g(i) with (g(i)+g(nlat-i+1))*wts(i)
c     overwrite g(nlat-i+1) with (g(i)-g(nlat-i+1))*wts(i)
      nl2 = nlat/2
      do 104 k=1,nt
      do 104 j=1,nlon
      do 105 i=1,nl2
      is = nlat-i+1
      t1 = g(i,j,k)
      t2 = g(is,j,k)
      g(i,j,k) = wts(i)*(t1+t2)
      g(is,j,k) = wts(i)*(t1-t2)
  105 continue
c     adjust equator if necessary(nlat odd)
      if (mod(nlat,2).ne.0) g(late,j,k) = wts(late)*g(late,j,k)
  104 continue
c     set m = 0 coefficients first
      m = 0
      call legin(mode,l,nlat,m,w,pmn,km)
      do 106 k=1,nt
      do 106 i=1,late
      is = nlat-i+1
      do 107 np1=1,nlat,2
c     n even
      a(1,np1,k) = a(1,np1,k)+g(i,1,k)*pmn(np1,i,km)
  107 continue
      do 108 np1=2,nlat,2
c     n odd
      a(1,np1,k) = a(1,np1,k)+g(is,1,k)*pmn(np1,i,km)
  108 continue
  106 continue
c     compute coefficients for which b(m,n) is available
      do 109 mp1=2,lm1
      m = mp1-1
      mp2 = m+2
c     compute pmn for all i and n=m,...,l-1
      call legin(mode,l,nlat,m,w,pmn,km)
      do 110 k=1,nt
      do 111 i=1,late
      is = nlat-i+1
c     n-m even
      do 112 np1=mp1,nlat,2
      a(mp1,np1,k) = a(mp1,np1,k)+g(i,2*m,k)*pmn(np1,i,km)
      b(mp1,np1,k) = b(mp1,np1,k)+g(i,2*m+1,k)*pmn(np1,i,km)
  112 continue
c     n-m odd
      do 113 np1=mp2,nlat,2
      a(mp1,np1,k) = a(mp1,np1,k)+g(is,2*m,k)*pmn(np1,i,km)
      b(mp1,np1,k) = b(mp1,np1,k)+g(is,2*m+1,k)*pmn(np1,i,km)
  113 continue
  111 continue
  110 continue
  109 continue
      if (nlon .eq. l+l-2) then
c     compute a(l,np1) coefficients only
      m = l-1
      call legin(mode,l,nlat,m,w,pmn,km)
      do 114 k=1,nt
      do 114 i=1,late
      is = nlat-i+1
c     n-m even
      do 124 np1=l,nlat,2
      a(l,np1,k) = a(l,np1,k)+0.5*g(i,nlon,k)*pmn(np1,i,km)
  124 continue
      lp1 = l+1
c     n-m odd
      do 125 np1=lp1,nlat,2
      a(l,np1,k) = a(l,np1,k)+0.5*g(is,nlon,k)*pmn(np1,i,km)
  125 continue
  114 continue
      end if
      else
c     half sphere
c     overwrite g(i) with wts(i)*(g(i)+g(i)) for i=1,...,nlate/2
      nl2 = nlat/2
      do 116  k=1,nt
      do 116 j=1,nlon
      do 115 i=1,nl2
      g(i,j,k) = wts(i)*(g(i,j,k)+g(i,j,k))
  115 continue
c     adjust equator separately if a grid point
      if (nl2.lt.late) g(late,j,k) = wts(late)*g(late,j,k)
  116 continue
c     set m = 0 coefficients first
      m = 0
      call legin(mode,l,nlat,m,w,pmn,km)
      ms = 1
      if (mode.eq.1) ms = 2
      do 117 k=1,nt
      do 117 i=1,late
      do 117 np1=ms,nlat,2
      a(1,np1,k) = a(1,np1,k)+g(i,1,k)*pmn(np1,i,km)
  117 continue
c     compute coefficients for which b(m,n) is available
      do 118 mp1=2,lm1
      m = mp1-1
      ms = mp1
      if (mode.eq.1) ms = mp1+1
c     compute pmn for all i and n=m,...,nlat-1
      call legin(mode,l,nlat,m,w,pmn,km)
      do 119 k=1,nt
      do 119 i=1,late
      do 119 np1=ms,nlat,2
      a(mp1,np1,k) = a(mp1,np1,k)+g(i,2*m,k)*pmn(np1,i,km)
      b(mp1,np1,k) = b(mp1,np1,k)+g(i,2*m+1,k)*pmn(np1,i,km)
  119 continue
  118 continue
      if (nlon.eq.l+l-2) then
c     compute coefficient a(l,np1) only
      m = l-1
      call legin(mode,l,nlat,m,w,pmn,km)
      ns = l
      if (mode.eq.1) ns = l+1
      do 120 k=1,nt
      do 120 i=1,late
      do 120 np1=ns,nlat,2
      a(l,np1,k) = a(l,np1,k)+0.5*g(i,nlon,k)*pmn(np1,i,km)
  120 continue
      end if
      end if
      return
      end
      subroutine shagci(nlat,nlon,wshagc,lshagc,dwork,ldwork,ierror)
c     this subroutine must be called before calling shagc with
c     fixed nlat,nlon. it precomputes quantites such as the gaussian
c     points and weights, m=0,m=1 legendre polynomials, recursion
c     recursion coefficients.
      dimension wshagc(lshagc)
      double precision dwork(ldwork)
      ierror = 1
      if (nlat.lt.3) return
      ierror = 2
      if (nlon.lt.4) return
c     set triangular truncation limit for spherical harmonic basis
      l = min0((nlon+2)/2,nlat)
c     set equator or nearest point (if excluded) pointer
      late = (nlat+mod(nlat,2))/2
      l1 = l
      l2 = late
      ierror = 3
c     check permanent work space length
      if (lshagc .lt. nlat*(2*l2+3*l1-2)+3*l1*(1-l1)/2+nlon+15)return
      ierror = 4
      if (ldwork.lt.nlat*(nlat+4))return
      ierror = 0
c     set pointers
      i1 = 1
      i2 = i1+nlat
      i3 = i2+nlat*late
      i4 = i3+nlat*late
      i5 = i4+l*(l-1)/2 +(nlat-l)*(l-1)
      i6 = i5+l*(l-1)/2 +(nlat-l)*(l-1)
      i7 = i6+l*(l-1)/2 +(nlat-l)*(l-1)
c     set indices in temp work for double precision gaussian wts and pts
      idth = 1
      idwts = idth+nlat
      iw = idwts+nlat
      call shagci1(nlat,nlon,l,late,wshagc(i1),wshagc(i2),wshagc(i3),
     1wshagc(i4),wshagc(i5),wshagc(i6),wshagc(i7),dwork(idth),
     2dwork(idwts),dwork(iw),ierror)
      if (ierror.ne.0) ierror = 5
      return
      end
      subroutine shagci1(nlat,nlon,l,late,wts,p0n,p1n,abel,bbel,cbel,
     1                  wfft,dtheta,dwts,work,ier)
      dimension wts(nlat),p0n(nlat,late),p1n(nlat,late),abel(1),bbel(1),
     1 cbel(1),wfft(1)
      double precision pb,dtheta(nlat),dwts(nlat),work(*)
c     compute the nlat  gaussian points and weights, the
c     m=0,1 legendre polys for gaussian points and all n,
c     and the legendre recursion coefficients
c     define index function used in storing
c     arrays for recursion coefficients (functions of (m,n))
c     the index function indx(m,n) is defined so that
c     the pairs (m,n) map to [1,2,...,indx(l-1,l-1)] with no
c     "holes" as m varies from 2 to n and n varies from 2 to l-1.
c     (m=0,1 are set from p0n,p1n for all n)
c     define for 2.le.n.le.l-1
      indx(m,n) = (n-1)*(n-2)/2+m-1
c     define index function for l.le.n.le.nlat
      imndx(m,n) = l*(l-1)/2+(n-l-1)*(l-1)+m-1
c     preset quantites for fourier transform
      call hrffti(nlon,wfft)
c     compute double precision gaussian points and weights
c     lw = 4*nlat*(nlat+1)+2
      lw = nlat*(nlat+2)
      call gaqd(nlat,dtheta,dwts,work,lw,ier)
      if (ier.ne.0) return
c     store gaussian weights single precision to save computation
c     in inner loops in analysis
      do 100 i=1,nlat
      wts(i) = dwts(i)
  100 continue
c     initialize p0n,p1n using double precision dnlfk,dnlft
      do 101 np1=1,nlat
      do 101 i=1,late
      p0n(np1,i) = 0.0
      p1n(np1,i) = 0.0
  101 continue
c     compute m=n=0 legendre polynomials for all theta(i)
      np1 = 1
      n = 0
      m = 0
      call dnlfk(m,n,work)
      do 103 i=1,late
      call dnlft(m,n,dtheta(i),work,pb)
      p0n(1,i) = pb
  103 continue
c     compute p0n,p1n for all theta(i) when n.gt.0
      do 104 np1=2,nlat
      n = np1-1
      m = 0
      call dnlfk(m,n,work)
      do 105 i=1,late
      call dnlft(m,n,dtheta(i),work,pb)
      p0n(np1,i) = pb
  105 continue
c     compute m=1 legendre polynomials for all n and theta(i)
      m = 1
      call dnlfk(m,n,work)
      do 106 i=1,late
      call dnlft(m,n,dtheta(i),work,pb)
      p1n(np1,i) = pb
  106 continue
  104 continue
c     compute and store swarztrauber recursion coefficients
c     for 2.le.m.le.n and 2.le.n.le.nlat in abel,bbel,cbel
      do 107 n=2,nlat
      mlim = min0(n,l)
      do 107 m=2,mlim
      imn = indx(m,n)
      if (n.ge.l) imn = imndx(m,n)
      abel(imn)=sqrt(float((2*n+1)*(m+n-2)*(m+n-3))/
     1               float(((2*n-3)*(m+n-1)*(m+n))))
      bbel(imn)=sqrt(float((2*n+1)*(n-m-1)*(n-m))/
     1               float(((2*n-3)*(m+n-1)*(m+n))))
      cbel(imn)=sqrt(float((n-m+1)*(n-m+2))/
     1               float(((n+m-1)*(n+m))))
  107 continue
      return
      end

max_stars_repo_name: xylar/cdat

max_issues_repo_name: xylar/cdat

max_forks_repo_name: CDAT/uvcdat


text:
Subroutine lulist(mlist)
  Common /lujets/n, k(9000, 5), p(9000, 5), v(9000, 5)
  Save /lujets/
  Common /ludat1/mstu(200), paru(200), mstj(200), parj(200)
  Save /ludat1/
  Common /ludat2/kchg(500, 3), pmas(500, 4), parf(2000), vckm(4, 4)
  Save /ludat2/
  Common /ludat3/mdcy(500, 3), mdme(2000, 2), brat(2000), kfdp(2000, 5)
  Save /ludat3/
  Character chap*16, chac*16, chan*16, chad(5)*16, chmo(12)*3, chdl(7)*4
  Dimension ps(6)
  Data chmo/'Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sep', 'Oct', 'Nov', 'Dec'/, chdl/'(())', ' ', '()', '!!', '<>', '==', '(==)'/
  If (mlist>=1 .And. mlist<=3) Then
    If (mlist==1) Write (mstu(11), 1100)
    If (mlist==2) Write (mstu(11), 1200)
    If (mlist==3) Write (mstu(11), 1300)
    lmx = 12
    If (mlist>=2) lmx = 16
    istr = 0
    imax = n
    If (mstu(2)>0) imax = mstu(2)
    Do i = max(1, mstu(1)), max(imax, n+max(0,mstu(3)))
      If ((i>imax .And. i<=n) .Or. k(i,1)<0) Goto 120
      Call luname(k(i,2), chap)
      len = 0
      Do lem = 1, 16
        If (chap(lem:lem)/=' ') len = lem
      End Do
      mdl = (k(i,1)+19)/10
      ldl = 0
      If (mdl==2 .Or. mdl>=8) Then
        chac = chap
        If (len>lmx) chac(lmx:lmx) = '?'
      Else
        ldl = 1
        If (mdl==1 .Or. mdl==7) ldl = 2
        If (len==0) Then
          chac = chdl(mdl)(1:2*ldl) // ' '
        Else
          chac = chdl(mdl)(1:ldl) // chap(1:min(len,lmx-2*ldl)) // chdl(mdl)(ldl+1:2*ldl) // ' '
          If (len+2*ldl>lmx) chac(lmx:lmx) = '?'
        End If
      End If
      If (k(i,1)==1 .Or. k(i,1)==2 .Or. k(i,1)==11 .Or. k(i,1)==12) Then
        kc = lucomp(k(i,2))
        kcc = 0
        If (kc/=0) kcc = kchg(kc, 2)
        If (kcc/=0 .And. istr==0) Then
          istr = 1
          If (len+2*ldl+3<=lmx) chac(lmx-1:lmx-1) = 'A'
        Else If (kcc/=0 .And. (k(i,1)==2 .Or. k(i,1)==12)) Then
          If (len+2*ldl+3<=lmx) chac(lmx-1:lmx-1) = 'I'
        Else If (kcc/=0) Then
          istr = 0
          If (len+2*ldl+3<=lmx) chac(lmx-1:lmx-1) = 'V'
        End If
      End If
      If (mlist==1 .And. abs(p(i,4))<9999.) Then
        Write (mstu(11), 1400) i, chac(1:12), (k(i,j1), j1=1, 3), (p(i,j2), j2=1, 5)
      Else If (mlist==1 .And. abs(p(i,4))<99999.) Then
        Write (mstu(11), 1500) i, chac(1:12), (k(i,j1), j1=1, 3), (p(i,j2), j2=1, 5)
      Else If (mlist==1) Then
        Write (mstu(11), 1600) i, chac(1:12), (k(i,j1), j1=1, 3), (p(i,j2), j2=1, 5)
      Else If (mstu(5)==10000 .And. (k(i,1)==3 .Or. k(i,1)==13 .Or. k(i,1)==14)) Then
        Write (mstu(11), 1700) i, chac, (k(i,j1), j1=1, 3), k(i, 4)/100000000, mod(k(i,4)/10000, 10000), mod(k(i,4), 10000), k(i, 5)/100000000, mod(k(i,5)/10000, 10000), mod(k(i,5), 10000), (p(i,j2), j2=1, 5)
      Else
        Write (mstu(11), 1800) i, chac, (k(i,j1), j1=1, 5), (p(i,j2), j2=1, 5)
      End If
      If (mlist==3) Write (mstu(11), 1900)(v(i,j), j=1, 5)
      If (mstu(70)>=1) Then
        isep = 0
        Do j = 1, min(10, mstu(70))
          If (i==mstu(70+j)) isep = 1
        End Do
        If (isep==1 .And. mlist==1) Write (mstu(11), 2000)
        If (isep==1 .And. mlist>=2) Write (mstu(11), 2100)
      End If
    120 End Do
    Do j = 1, 6
      ps(j) = plu(0, j)
    End Do
    If (mlist==1 .And. abs(ps(4))<9999.) Then
      Write (mstu(11), 2200) ps(6), (ps(j), j=1, 5)
    Else If (mlist==1 .And. abs(ps(4))<99999.) Then
      Write (mstu(11), 2300) ps(6), (ps(j), j=1, 5)
    Else If (mlist==1) Then
      Write (mstu(11), 2400) ps(6), (ps(j), j=1, 5)
    Else
      Write (mstu(11), 2500) ps(6), (ps(j), j=1, 5)
    End If
  Else If (mlist==11) Then
    Write (mstu(11), 2600)
    Do kf = 1, 40
      Call luname(kf, chap)
      Call luname(-kf, chan)
      If (chap/=' ' .And. chan==' ') Write (mstu(11), 2700) kf, chap
      If (chan/=' ') Write (mstu(11), 2700) kf, chap, -kf, chan
    End Do
    Do kfls = 1, 3, 2
      Do kfla = 1, 8
        Do kflb = 1, kfla - (3-kfls)/2
          kf = 1000*kfla + 100*kflb + kfls
          Call luname(kf, chap)
          Call luname(-kf, chan)
          Write (mstu(11), 2700) kf, chap, -kf, chan
        End Do
      End Do
    End Do
    Do kmul = 0, 5
      kfls = 3
      If (kmul==0 .Or. kmul==3) kfls = 1
      If (kmul==5) kfls = 5
      kflr = 0
      If (kmul==2 .Or. kmul==3) kflr = 1
      If (kmul==4) kflr = 2
      Do kflb = 1, 8
        Do kflc = 1, kflb - 1
          kf = 10000*kflr + 100*kflb + 10*kflc + kfls
          Call luname(kf, chap)
          Call luname(-kf, chan)
          Write (mstu(11), 2700) kf, chap, -kf, chan
        End Do
        kf = 10000*kflr + 110*kflb + kfls
        Call luname(kf, chap)
        Write (mstu(11), 2700) kf, chap
      End Do
    End Do
    kf = 130
    Call luname(kf, chap)
    Write (mstu(11), 2700) kf, chap
    kf = 310
    Call luname(kf, chap)
    Write (mstu(11), 2700) kf, chap
    Do kflsp = 1, 3
      kfls = 2 + 2*(kflsp/3)
      Do kfla = 1, 8
        Do kflb = 1, kfla
          Do kflc = 1, kflb
            If (kflsp==1 .And. (kfla==kflb .Or. kflb==kflc)) Goto 180
            If (kflsp==2 .And. kfla==kflc) Goto 180
            If (kflsp==1) kf = 1000*kfla + 100*kflc + 10*kflb + kfls
            If (kflsp>=2) kf = 1000*kfla + 100*kflb + 10*kflc + kfls
            Call luname(kf, chap)
            Call luname(-kf, chan)
            Write (mstu(11), 2700) kf, chap, -kf, chan
          180 End Do
        End Do
      End Do
    End Do
  Else If (mlist==12) Then
    Write (mstu(11), 2800)
    mstj24 = mstj(24)
    mstj(24) = 0
    kfmax = 20883
    If (mstu(2)/=0) kfmax = mstu(2)
    Do kf = max(1, mstu(1)), kfmax
      kc = lucomp(kf)
      If (kc==0) Goto 220
      If (mstu(14)==0 .And. kf>100 .And. kc<=100) Goto 220
      If (mstu(14)>0 .And. kf>100 .And. max(mod(kf/1000,10),mod(kf/100,10))>mstu(14)) Goto 220
      Call luname(kf, chap)
      If (kf<=100 .And. chap==' ' .And. mdcy(kc,2)==0) Goto 220
      Call luname(-kf, chan)
      pm = ulmass(kf)
      Write (mstu(11), 2900) kf, kc, chap, chan, kchg(kc, 1), kchg(kc, 2), kchg(kc, 3), pm, pmas(kc, 2), pmas(kc, 3), pmas(kc, 4), mdcy(kc, 1)
      If (kf>100 .And. kc<=100) Goto 220
      Do idc = mdcy(kc, 2), mdcy(kc, 2) + mdcy(kc, 3) - 1
        Do j = 1, 5
          Call luname(kfdp(idc,j), chad(j))
        End Do
        Write (mstu(11), 3000) idc, mdme(idc, 1), mdme(idc, 2), brat(idc), (chad(j), j=1, 5)
      End Do
    220 End Do
    mstj(24) = mstj24
  Else If (mlist==13) Then
    Write (mstu(11), 3100)
    Do i = 1, 200
      Write (mstu(11), 3200) i, mstu(i), paru(i), mstj(i), parj(i), parf(i)
    End Do
  End If
  Return
  1100 Format (///28X, 'Event listing (summary)'//4X, 'I  particle/jet KS', 5X, 'KF orig    p_x      p_y      p_z       E        m'/)
  1200 Format (///28X, 'Event listing (standard)'//4X, 'I  particle/jet', '  K(I,1)   K(I,2) K(I,3)     K(I,4)      K(I,5)       P(I,1)', '       P(I,2)       P(I,3)       P(I,4)       P(I,5)'/)
  1300 Format (///28X, 'Event listing (with vertices)'//4X, 'I  particle/j', 'et  K(I,1)   K(I,2) K(I,3)     K(I,4)      K(I,5)       P(I,1)', '       P(I,2)       P(I,3)       P(I,4)       P(I,5)'/73X, 'V(I,1)       V(I,2)       V(I,3)       V(I,4)       V(I,5)'/)
  1400 Format (1X, I4, 2X, A12, 1X, I2, 1X, I6, 1X, I4, 5F9.3)
  1500 Format (1X, I4, 2X, A12, 1X, I2, 1X, I6, 1X, I4, 5F9.2)
  1600 Format (1X, I4, 2X, A12, 1X, I2, 1X, I6, 1X, I4, 5F9.1)
  1700 Format (1X, I4, 2X, A16, 1X, I3, 1X, I8, 2X, I4, 2(3X,I1,2I4), 5F13.5)
  1800 Format (1X, I4, 2X, A16, 1X, I3, 1X, I8, 2X, I4, 2(3X,I9), 5F13.5)
  1900 Format (66X, 5(1X,F12.3))
  2000 Format (1X, 78('='))
  2100 Format (1X, 130('='))
  2200 Format (19X, 'sum:', F6.2, 5X, 5F9.3)
  2300 Format (19X, 'sum:', F6.2, 5X, 5F9.2)
  2400 Format (19X, 'sum:', F6.2, 5X, 5F9.1)
  2500 Format (19X, 'sum charge:', F6.2, 3X, 'sum momentum and inv. mass:', 5F13.5)
  2600 Format (///20X, 'List of KF codes in program'/)
  2700 Format (4X, I6, 4X, A16, 6X, I6, 4X, A16)
  2800 Format (///30X, 'Particle/parton data table'//5X, 'KF', 5X, 'KC', 4X, 'particle', 8X, 'antiparticle', 6X, 'chg  col  anti', 8X, 'mass', 7X, 'width', 7X, 'w-cut', 5X, 'lifetime', 1X, 'decay'/11X, 'IDC', 1X, 'on/off', 1X, 'ME', 3X, 'Br.rat.', 4X, 'decay products')
  2900 Format (/1X, I6, 3X, I4, 4X, A16, A16, 3I5, 1X, F12.5, 2(1X,F11.5), 2X, F12.5, 3X, I2)
  3000 Format (10X, I4, 2X, I3, 2X, I3, 2X, F8.5, 4X, 5A16)
  3100 Format (///20X, 'Parameter value table'//4X, 'I', 3X, 'MSTU(I)', 8X, 'PARU(I)', 3X, 'MSTJ(I)', 8X, 'PARJ(I)', 8X, 'PARF(I)')
  3200 Format (1X, I4, 1X, I9, 1X, F14.5, 1X, I9, 1X, F14.5, 1X, F14.5)
End Subroutine lulist

max_stars_repo_name: xiaohaijin/AMPT

max_issues_repo_name: xiaohaijin/AMPT

max_forks_repo_name: xiaohaijin/AMPT


text:
This is a simple list of the people of Davis, pretty heavily biased to those who use the wiki (since this list is compiled by people who use the wiki). You can also see some User Statistics. When people get together they form Social Structures.

You are part of Davis, so add yourself: make a homepage for yourself then add it to this list! Dont be shy! To make a homepage for yourself, click on your username in the upperright corner of the screen OR enter your username into the search box. Then, choose to create a new page with this title (or maybe if theres already a page for you, then youre awfully lucky!) Then just fill in your page with whatever you want. This is, in general, how you create a new page in the wiki without adding a link to an existing page. Keep format, and add after youve made yourself a nice page, keeps Wanted Pages down!

Then, once youve got a homepage, come back here and add yourself, but keep it alphabetized so Users/AmitVainsencher doesnt cry!

Or add other people. Some people are so notable locally that they have entries even though they havent edited the wiki. Please please please add people if you think theyre relevant to Davis! The more the merrier! The only reason most of the names on there are of people who have wiki accounts is that no one has taken the time to flesh it out more. As with all things on the Wiki, theres a selfselection bias. People with accounts are the ones who add stuff. One cautionary note: if the names youre adding arent those of public figures, its probably better to check with the people first to make sure they dont mind having their name on the page.

There are also some Wiki Animals that use the wiki (or their owners do...)
 >
<
 A
  Users/KyleAaron
  Users/AbcDef
  Users/TomAbeles
  Users/ArlenAbraham
  Users/MichelleAccurso
  Users/JustinAi
  Debbie Aldridge
  Users/MistressAlthaea
  Users/JanelleAlvstadMattson
  Users/PaulAmnuaypayoat
  Users/AllyAnaheim
  Users/AdamAnderson
  Users/EricaAnderson
  Users/ErikAnderson
  Users/TojoAndrianarivo
  Users/BrianAng
  Users/Angel.York
  Users/ThomasAngell
  Users/PeterAnselmo
  Users/AnthonySantos
  Users/AraAroyan
  Users/EthanAshley
  Users/SyedAslam
  Users/AlinaAsrian
  Frank AtaBaah
  Users/PxlAted
  Users/JeffAustin
  Users/AlanAuyeung
  Users/ChelseaAyers

 B
  Users/VishalBhula
  Users/BrianWickstrom
  Users/Jenb
  Users/RonB
  Users/Babhari
  Users/LeslieBach
  Users/BlairBaldwin
  Users/JamesBanks
  Users/AndrewBanta
  wiki:wikispot:users/BaltoCo BaltoCo
  Users/KateBaker
  Users/CarlosBarahona
  Users/BrianBarbier
  Users/CaseyBarc
  Users/DavidBarnum David Barnum
  Users/RobertBaron
  Users/FredBarrett
  Users/BrandonBarrette
  Users/NicholasBarry Nicholas Barry (formerly users/NaBarry)
  Users/ElizabethBarthel
  Users/HilaryBartle
  Users/ToddBartlow
  Users/PeteBasofin
  Users/TonyBataska
  Herbert Bauer
  Users/PayamBavafa
  Users/jbeda
  Users/AndrewBell
  Users/BryanBell
  Jay Belsky
  Users/AndrewBenson
  Users/JillBenciWoodward
  Users/BradBenedict
  Users/DavidBenjamin
  Users/JoeyBennett
  Users/PhilBerghuis
  Users/MarielleBerman
  Eve West Bessier
  Users/AndrewBianchi
  Users/JuliAnnBlanco
  Users/ScottBiggart
  Users/WilliamBinnsSmith
  Users/BradBitler
  Users/KristenBirdsall
  Rebecca Black is a local etiquette coach.
  Users/EmilyBlake
  Users/JosephBleckman
  Users/KenBloom
  Users/BarryBolden
  Users/ClaireBolding
  Users/MarieBoisvert
  Users/BoogabaahWeesnaah
  Users/borax
  Lourdes E. Borrayo
  Users/CaseyBotkin
  Users/jbowlus
  Users/RussBowlus
  Users/JeffBoyce
  Users/BrendanBoyle
  Users/ChrisMoBradley
  Barbara Brady
  Users/CurtisBray
  Users/TeddyBriggs
  Users/ArianaBrill
  Users/BritMSnow
  Users/ErnBro
  Users/BrookeB
  Users/DanielBrown
  Users/MichaelBrown
  Users/RyanBrown
  Users/NicoleLBrowner
  Users/CraigBrozinsky
  Users/DeannaBrusco
  Users/KonradBrynda
  Users/GerritBuddingh
  Users/VinceBuffalo
  Users/ChristineBui
  Users/AdamBurgett
  Users/MikeyBurke
  Users/DavidBurns
  Users/CliftonBurton
  Users/NicoleBush
  Users/BreeButler
  Users/RachelButler
  Users/RichardByron Richard Zoogy Byron

 C
  Users/NumiaCairaguas
  Users/VictorCairel
  Users/JasminCamacho
  Users/BrentCampos
  Users/MinervaCelesteCardona
  Users/GennaCarnes
  Users/AmeliaCarlson
  users/Calvin
  Users/ascapoccia
  Users/MikeCarr
  Users/NoraCary
  Users/JasonCash
  Users/RyanCastellucci
  Users/RaynatoCastro
  Users/BenjaminCatabas
  Users/JoeCattarin
  Users/LizethCazares
  James G. Cecil
  Users/PaulCernea
  Users/IrinaChakraborty
  Users/ChristinaChang
  Users/EvelinaChang
  Users/AndrewChen
  Users/BrianChen
  Users/JoeyChen
  Users/WiniChen
  Users/YawenChen
  Users/AricaChhay
  Users/AliceChoe
  Users/AndrewCholerton
  Users/DorotheaChow
  Users/Chris47
  Users/PopeChristof
  Users/JenniferChu
  Users/IanCipperly
  Users/ChristopherCivil
  Users/CamilleClark
  Users/RogerClark
  Users/AinsleyClarkGaveston
  Users/TerryCliff
  Users/EmmaCoats
  Users/KaceyCoburn
  Users/LauraCochrane
  users/JamesColligan JamesColligan
  Users/condemned2bfree
  Users/DeeConway DeeConway
  Users/BlakeNCooper Blake N. Cooper
  Users/KellyCorcoran
  Users/AnaluciaCordova
  Users/NicoleCorrea
  George Covell
  Users/JackkiCox
  Users/KellyCox
  Users/BrianneCrabtree
  Users/KristineCraig
  Users/LukeCrawford
  Users/LawrenceCuevas
  Users/TorreyaCummings
  Users/PeterCung
  Users/LynneCunningham
  Users/BradCuppy
  Users/CurlyGirl26
  Users/Aaron.Curtin
 D

  Users/BrianDarnell
  Users/EyadDarras
  Users/ghdaskalos
  Users/LisaDavis
  Users/SteveDavison
  Users/ShyDavisite
  Users/ShannonDeArmond
  Users/KatieDelwiche
  Users/RoseAnneDeCristoforo
  Users/PaulDeLeon
  Ashley Detchman
  Users/ErikadelosReyes
  Users/NinadelRosario
  Users/CristinaDeptula
  Users/MandeepDhaliwal
  Users/JakeDiego
  Users/ChrisDietrich
  Users/AlinaDini
  Users/IsaacDivide
  Users/DinaDobkin
  Users/JenDockter
  Users/NathanDodd
  Users/DerekDowney
  Users/PatrickDragon
  Users/BrigitteDriller
  Users/JonDsouza
  Users/JohnDudek
  Users/MichaelDuffey
  Users/MikeDuffy
  Users/JasonDunne
  Users/ZitaDemaree

 E

  Users/CraigEasdale
  Users/LieshaEberst
  Users/AbdolhosseinEdalatiSarayani
  Users/JaneEdberg
  Users/edubs
  Evan Users/JabberWokky Edwards
  Users/JulieEickhof
  Users/JonathanEisen
  Users/RickEle
  Users/AlexEley
  Users/StaceyEllis
  Users/Elly
  Users/KristElohim
  Users/JulienBiewerElstob Julien Elstob
  Users/LaurelE Laurel Engbring
  Users/AikoEnoki
  Users/AllisonEriksen
  Users/AmyErven
  Users/VinceEstremera
  Users/MarionEveridge
  Users/GonzaloEyzaguirre

 F

  Users/FirstStreetRealEstate
  Users/EmyFargey
  Users/JoshFeldman
  Users/FeralFeline
  Users/TylerFelix
  Users/DavidFeliz
  Users/ShiFeng
  Users/JoshFernandez
  Users/ThomasFifield
  Users/JoeFinkel
  Users/KennethFirestein
  Users/DanFisher
  Users/MattFisher
  Users/LucyFitzGibbon Lucy Fitz Gibbon
  Users/BryceFitzsimons
  Users/KyleFlick
  Users/MorganFlores
  Users/SylviaFlores
  Users/AdamFlowers
  Users/PatrickFoley
  Users/DavidFong
  Jack Forbes
  Users/BibbleFoshizzle
  Users/FreelanceFoodCritic
  wiki:wikipedia:Michael_Franti  See Significant Davis Bands
  Users/ColinFreeman
  Users/GrahamFreeman
  Tenaya Freitas
  Users/MelissaFreixas
  Users/KrisFricke
  Users/TobinFricke
  Users/SimonFung
  Users/RussellFurey
  Users/JuliaRoseFurman
  Users/GregFurnish
  Users/PeterFuszard
  Users/WengFF

 G

  Cliff Gamble
  Users/ManiGandham
  Users/TomGarberson
  Users/AlexGarbutt
  Users/MeganGarcia
  Users/ShaunGeer
  Dave GertmenianWong
  Users/AdamGhaziTehrani
  Users/JenniferGiang
  Users/MichaelGiardina
  Users/JeremyGimbel
  Users/MeredithGlaser
  Users/DanGlendening
  Users/JoshuaGo
  Users/AndreyGoder
  Users/SaraGoldband
  Users/IanGolden
  Users/DanielGonzales
  Users/MichaelGoodell
  Users/LaurenGoshen
  Users/TheodoreGould
  Users/LenaGranik
  Users/SamanthaGrant
  Users/TravisGrathwell
  Emil F. Grauel
  Users/AlexGray
  Users/BobbyGray
  Users/JevanGray
  Users/DominikGreen
  Users/AnneGreenley
  Users/perrygregg Perry Gregg
  Gary R. Grelecki
  Ray Grelecki
  Users/TylerGriffith
  Users/MatthiasGropp
  Users/DavidGrundler
  Users/GracielaGuardado
  Users/GiladGurantz
  Users/JaniceGutshall
  Cape Guy

 H

  Users/BriannaHaag
  Users/JamesHaile
  Users/FredHall
  Users/MargieHalloran
  Users/RodolfoHamZhu
  Users/IsaacHamlenGomez
  Users/StephanieHammon
  Users/JasonHammons
  Users/AdamHancock
  Users/OlinHannum
  Bruce Hansen, Users/BruceHansen Users/BruceHansen
  Users/WesHardaker
  Users/JasonHarden
  Users/ChristineHarms
  Users/PaulHarms
  Users/ArielaHaro
  Users/AndreHarris
  Kristin Harris
  Users/MaxHart
  Users/LauraHartman
  Users/JackHaskel
  Users/BryceHathaway
  Users/BobbyHayashi
  Users/BrooksHayes
  Ernie Head
  Users/mheathershaw
  Users/KristyHeidenberger
  Mark Heisey
  Users/JoshuaHeller
  Users/LouHelt
  Users/EdHenn
  Users/SallyHensel
  Users/MichaelHernandez
  Users/MattHh
  Jack Hicks  UCD English professor and author
  Users/ScottHildebrand
  Users/EmilyHill
  Users/LeightonHinkley
  Users/BroodingHippie
  Users/AlexHirsch
  Users/JamesHirsh
  Users/GregHirson
  Users/WendellHise
  Users/StephenHjellum
  Users/DonHo om nom nom nom
  Users/EileenHo
  Users/StephenHo
  Users/AleeceHoffmann
  Users/DannieHolden
  Users/BrettHolycross
  Users/MilesHookey
  Users/NiareeHopelian
  Users/JoshHorne
  Users/ElisaHough
  Users/WesHough
  Von Howard
  Liz Hummer Users/LizHummer
  Users/JohnHumperdinkle
  Users/MarieHuynh
  Users/MinhHuynh
  Users/L.A.Hernandez
  Users/MichelleHV

 I

  Users/IanBurau
  Users/VanessaIrving
  Users/IasriIsrin
  Users/IDoNotExist
  Users/AlishaIsrael
  Users/MikeIvanov
  Users/PaulIvanov

 J
  Hair Saloon Sharon Jackson
  Users/JustinKumar
  Users/NamitJain
  Users/Jeff
  Users/TimJarratt
  Users/NathanielJarrett
  Users/RohiniJasavala
  Users/GavinJensen
  Users/jeanpierre
  Users/JamesJobe
  Users/WillJobe
  Users/SeanJohannessen
  Users/AubreyJohnson
  Users/ColeJohnson
  Users/GeoffJohnson
  Users/JessicaJohnson
  Users/KatiJohnson
  Users/ZachariahJohnson
  Users/MattJojola
  Users/AndyJones
  Users/AnnaJones
  Users/AnsateJones
  Users/ChristaJones
  Users/DaveJones
  Users/DagonJones
  Users/ElliotJorgensen
  Users/GautamJoshi
  Users/MaxJ
  Users/MattJurach
  Users/NikhilJoshi

 K

  Users/kelvinkay
  Users/RobinKahn
  Users/ToddKaiser
  Users/JeffreyKalvass
  Users/MaxKam
  Elise Kane
  Users/DylanKann
  Users/JenKao
  Sahas Katta
  Users/MiriamKaufman
  Users/WillKaufman
  Users/RyanKazoo
  Users/JohnKearney
  Users/BrandonKeene
  Users/KrystenKellum
  Users/DiorKelley
  Users/DonovanKelley
  Users/BillKendrick
  Users/TeresaKenny
  Users/AlisonKeenan
  Users/AlisonKent
  Users/BrandonHKey
  Users/ImranKhan
  Users/OmidKhandan
  Users/JustinKhoo
  Users/YongKil
  Users/hankim
  Users/PhillipKim
  Users/ShawnBKing
  Users/AlexKistner
  Users/NatalieKitchiner
  Users/MarkKlebanov
  Users/EricKlein
  Users/ErikKlemetti
  Users/HeatherKlinger
  Users/MeganKlugh
  Users/FKnochenhauer
  Users/EthanKnoop
  Users/SandyKnowles
  Users/MichelleKoehler
  Users/GrahamKolbeins
  Users/JeremyKoo
  Users/JeyKottalam
  Users/GabeKoulikov
  Users/LindseyKrabbenhoft
  Users/JamesKrause
  Joe Krovoza  the Mayor of Davis
  Users/JenniferKucich
  Users/KimKunaniec
  Users/WinnieKuo
  Users/InnaKurikova
  Users/AlexKwan
  Users/CherylKysh

 L

  Users/BrentLaabs
  Users/JeffLagasca
  Users/JonathanLai
  Users/AndrewLam
  Users/JasonLamar
  Users/EricLambden
  Users/TrishaLambert
  Users/TomLangston
  Users/GlennLangton
  Users/PaulLanzi
  Users/AdamLau
  Users/JasmineLau
  Users/KitLau
  Users/AbbyLawson
  Users/JonathanLawton
  Users/AnneLe
  Users/CarolLe
  Users/AlexLeach
  Users/MichaelLeahy
  Users/JonathonLeathers
  Users/DanielleLee
  Users/HollyLee
  Users/LanaLee
  Users/MatthewLee
  Users/SharonLee
  Users/StevenLee
  Users/DanielLebowitz
  Users/JaniceLee
  Users/BurkhardLehner
  Users/thelonepiper555 Erik Leiken
  Users/AlexisLemus Alexis Lemus
  Users/AndrewLeonard
  Users/HalonaLeung
  Users/DavidLevin
  Users/GeorgeLewis
  Users/SpeedbagLewis
  Users/WilliamLewis
  Users/MaryLieth
  Users/ChristineLim
  Users/SteffanyLim
  Users/RayLin
  Users/RichLindvall
  Users/HilarieLloyd
  Users/ThomasLloyd
  Users/MamaLoebs
  Users/Oleg.Lokhvitsky Oleg Lokhvitsky
  Users/JakeLosh
  Users/MattLow
  Users/AlexLubinski
  Users/MaxLucas
  Users/JessicaLuedtke
  Robert Lugo
  Users/JulieLujano
  Users/BrianLum
  Users/HowardLuong
  Users/VickiLy
  Users/PatrickLynch
  Users/DeniseLiggettMcMerrick

 M
  Users/MelanieMadden
  Pete Madrone
  Users/TonyMagagna
  Users/JamesMah
  Users/AngelinaMalfitano
  Users/JosephMalki
  Users/AlexMandel
  Users/ErinManley
  Users/JakeMann
  Users/MattisManzel
  Users/MattHh
  Users/AlyssaMarkham
  Users/BenMarkus
  Users/SaraMarostica
  Users/LeonardMarque
  Users/RichardMarquez
  Users/ChristyMarsden
  David Martin
  Users/KristenMartin
  Tony Martin
  Users/BrianMartinez
  Users/DanMasiel
  Users/AshleyMatson
  Users/KristenMatsumura
  Users/TimMatthiessen
  Users/mbindrajana
  Users/CarlMcCabe
  Users/GregMcClanahan
  Users/BrianMcFadden
  Users/ShaunaMcKellar
  Users/AllisonMears
  Users/m610 Mike Meier
  Users/MikeMele
  Users/ShawnaMellon
  Users/JessicaMerrill
  Users/CraigMerry
  Users/LizMerry
  Users/MaxMikalonis
  Users/MaxStevenson
  Users/chrisfmichaels
  Users/Michellaneous
  Users/RyanMikulovsky
  Users/AllisonMiller
  Users/BenjaminMiller
  Users/JesseMiller
  Users/ScottMiller
  Users/PaulMin
  Users/Ming
  Users/AdamMoerschell
  Users/KarlMogel
  Users/JohnMolina
  Users/MikeMolina
  Paul Moller
  Users/GruntledMonkey
  Users/MariahMontanez
  Users/HannahMoore
  Users/ZahraMorad
  Users/StasMoreinis
  Users/AlysonMorgan
  Users/ScottMorgan
  Users/ZacMorris
  Users/CaitlinMorrow
  Users/TaylorMoseley
  Users/PaulMulase
  Users/DennisMuldoon
  Users/MasonMurray
  Users/LenoreMyers

 N

  Users/MattNagel
  Users/BrianNeal
  Users/TomNelson
  Users/PhilipNeustrom
  Users/KenNeville
  Users/ShannonNewbold
  Users/JeffNewmiller
  users/Ben
  Peter Ng
  Users/DatNguyen
  Users/TysonNichols
  Users/HarrisonNoah
  Users/ErikNoftle
  Users/MikeyNolan
  Users/BlancheNonken
  Users/JeffreyNonken
  Users/ZacharyNorwood
  Users/AlexNorris
  Users/NotTires
  Users/TimothyNutter

 O

  Users/YooHyunOak
  Users/MarenObrien
  Users/NicoleOcampo
  Users/ReiOkamoto
  Users/HollywoOd
  Users/TimOldenkamp
  Users/CristinaOlaverri
  Users/StarrOleary
  Users/KirstenOnell
  Users/LoriOrf
  Users/AshleyOrsaba
  Users/KennyOrtiz
  Users/CarlosOverstreet

 P
 Users/Phillco
  Users/CecilioPadilla
  Users/SirPallas
  Users/MichellePalmer
  Users/AnthonyParisi
  Users/KirkParrott
  Users/SusanParson
  Users/AndyPastalaniec
  Users/PenguiN42
  Users/GregPera
  Users/EdPereira
  Users/JustinPerkins
  Users/TracyPerkins
  Users/marilynlperry
  Users/JosephdelPesco
  Users/PeteForsyth
  Users/RachelPeters
  Users/AliPezeshkpour
  Daniel Pfister
  Users/CoreyPham
  Users/VivianPham
  Users/TerencePhan
  Users/LizPhillips
  Users/RopponPicha
  Users/GarrettPierce
  Users/ronpina
  Users/AnthonyPineda
  Users/JustinPolanich
  Users/RobertPolyn
  Users/AlexPomeranz
  USERS/JoePomidor
  Users/DavidPoole Dave/DavidPoole
  Users/KemblePope
  Users/DavidPrenatt
  Users/RyanPrendiville
  Users/NoahPretentious
  Users/ChristopherPrice
  Users/EricPrice
  Users/TravisPritchard
  Users/CovertProfessor
  Users/ElenkaProulx
  Users/SandeepPuri
  Users/parisperrault
  Users/PxlAted

 Q

  Users/AlexQuan
  Users/DanielQuincy
  Users/LucasQuass
  Users/GQuincy

 R

  Users/JaimeRaba
  Users/GadielRachelson
  Users/TheRadish
  Users/AllanRae
  Users/DanielRaff
  Users/KarthikRam
  Users/DoctorRama
  Users/IanRamsay
  Users/TusharRawat Tushar
  Users/ShiraRawlinson
  Walter G. Read
  Users/DaveyReed
  Users/SeanReedy
  Users/HannahReff
  Users/MarcRehfuss
  Users/DavidReid
  Users/EfremRensi
  Users/VanessaRenwick
  Users/AynReyes
  Users/BarryRice
  Users/KristinRicherson
  Keith Richman
  Users/KalenRidenour
  Users/RomanRivilis
  Users/RiotInDavis
  Users/JoRo
  Users/StephanieRobinson
  Users/JacobRoche
  Users/BradRodriguez
  Bernard Rogers
  Users/PatrickRogers Patrick Rogers
  Users/ErikRonneberg
  Users/MattRoper
  Users/StevenRoseJr

  Users/AaronRosenberg
  Users/DavidRosenberg
  Users/BenjaminRosenstein
  Users/FernandaGuimaraesRosa
  Users/ErinRoss
  Users/JeremyRoss
  Users/StephanieRoss
  Users/RobRoy
  Users/AndresRuiz
  Users/JoeRunnels
  Users/CaptainRush
  Users/SharonRuth
  Users/JenniferRutherford
  Users/EricRuud
  Users/JessicaRuygt

 S

  Users/Salie
  Users/EdwinS
  Users/ChrisSaito
  Users/AmiraSamaha
  Users/SamSampson
  Users/NoraSandstedt
  Users/ArielSanJose
  Users/GarySaylin
  Users/NickSchmalenberger
  Users/JamesSchwab
  Users/CaseyScott
  Users/AldenSeabolt
  Users/EllenSeeley
  Users/JeniferSelke
  Users/SamSellers
  Users/JulieSetele
  Users/TheShah
  Users/AlliraSharma
  Users/JeffShaw
  Users/MariSherbert
  Users/JoeSherman
  Users/AmberShipherd
  Users/IanSierchio
  Users/MikeSiminitus
  Users/SabrinaSimonton
  Users/PatrickSing
  Users/AviSingh
  Users/JesseSingh
  Users/MarcusSitz
  Users/EmilySiu
  Users/RachelSkytt
  Users/JonathanSlakey
  Users/KaliSlusser
  Users/JordanSmart
  Users/CharlesSmith
  Users/ChristopherSmith
  Users/GregorySmith
  Users/NateSmith
  Users/NathanSmith Nathan Smith
  Users/SaraSmith
  Users/TerraSmithers
  Users/strawberry
  Users/SMCareBare
  Users/SimiSohota
  Users/PanSola
  Users/SummerSong
  Users/DaleSorbello
  Users/BettySousa
  Users/RichardSpencer
  Users/CindySperry
  Users/DougSprankling
  Users/OraSraboy
  Users/KatStambuk
  Users/CaitlinthebestStanley
  Users/AndreaStewart
  Users/BonnieStewart
  Users/TeenaStockert
  Users/TaylorStreet
  Users/DanaSullivan
  Users/KarinaSummers
  Users/JennSuzuki
  Users/JoannaSwan
  Users/SydneySwartz
  Users/BevSykes
  Users/ElliotSzkup

 T

  Lucas Taber
  Users/susant
  Users/ToniTabora
  Users/TaiTakahara
  Users/ChrisTakemura
  Users/NatashaTavakoli
  Users/ConnieTaxiera
  Users/BenTaylor
  Users/JeremyTaylor
  Users/RJTaylor
  Users/TravisTaylor
  Users/JohnTerning
  Users/SorryAboutThat
  Users/BruceThomas
  Users/JonThomas
  Users/MattThornburg
  Users/IsaacDivide IsaacThornton
  Users/ThUn
  Users/thejoedude
  Users/NickTiano
  Users/GarrettTiemann
  Users/KaiTing
  Users/ChristinaTom
  Users/MoTorres
  Users/StacieTownsend
  Users/MeghanTracy Meghan Tracy
  Users/MaianTran
  Users/MarshaTran
  Users/MinhTran
  Users/RishiTrivedi
  Users/AlvinTsao
  Users/EricTsay
  Users/AnnieTsou
  Users/EmilyTung
  NAM TURK
  Users/TusharRawat Tushar
  Users/DawnTibbits Dawn Tibbits
  Users/T.Zukumori T.Zukumori
  Users/TonyKe

 U

  Users/RobertUllrey
  Users/JesseUnger

 V

  Users/MichaelVanderArend
  Users/AmitVainsencher
  Users/AnthonyValenzuela
  Users/RevChad ChadVanSchoelandt
  Users/JenndelaVega
  Users/JerryVigil
  Users/KrystalVillanueva
  Users/maxvidrine
  Users/ClarkVineyard
  Users/PatrickVisperas
  Users/RainbowVogt
  Users/GraceVoong
  Maveric Vu

 W

  Users/EdgarWai
  Users/SamanthaWagner
  Users/LaineWalkerAvina
  Users/willwallace
  Users/TigWallis
  Users/ColinWalsh
  Users/DougWalter
  Users/NingWan
  Users/RebeccaWan
  Users/HelenWang
  Users/TimWaring
  Users/MichelleWarren
  wiki:westsac:users/MasumiWatanabe Masumi Watanabe
  Users/KennethWaters
  Greg Webb
  Users/EricWecker
  Users/KyleWecker
  Users/Weesito
  Users/BoogabaahWeesnaah
  Users/JulieWeiss
  Lee Weissmann
  Users/JillWeinstein
  Users/HilaryWeisert
  Users/ARWENNHOLD
  Users/DebWestergaard
  Users/RobWestergaard
  Users/ShanaWhite
  Users/carl.whithaus
  Users/DarinWick
  Users/JudsonWickham
  Users/PaulWilburn
  Users/JamesWilk
  Users/KatieWillheim
  Users/PrevinWitana
  Users/MichaelWoffenden
  Users/EricWilliamWong
  Users/StevenWong
  Users/LaddWoodland
  Users/CathyWoolery
  Users/MichaelAWoods
  Users/StevenWorker
  Natalie Wormeli
  Users/DanielWorthington
  Bob Wren
  Users/MattWright
  Users/RoyWright
  Users/ColinWu
  Users/EricWu
  Users/FrankWu
  Users/VannyWu
  Users/JeffWurz

 X

  Users/DanXie
  Users/XavierHenry
  Users/XeDty

 Y

  Users/KenjiYamada
  Users/HannahYang
  users/EliYani
  Users/AdamYergovich
  Gary Lee Yoder
  Users/TamaraYurovsky
  Users/KenjiYonemoto
  Users/PamelaYoo
  Users/Kirstenelise
  Users/HaiYu

 Z

  Users/VictorZapanta
  Users/AmyZimmerman
  Users/SharonZimmerman
  Users/TianxiaZhou


text:
C$Procedure      PODIED ( Pod, insert elements, double precision )

      SUBROUTINE PODIED ( ELEMS, N, LOC, POD )
      IMPLICIT NONE

C$ Abstract
C
C     Insert elements at a specified location within the active group
C     of a pod.
C
C$ Disclaimer
C
C     THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE
C     CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S.
C     GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE
C     ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE
C     PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS"
C     TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY
C     WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A
C     PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC
C     SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE
C     SOFTWARE AND RELATED MATERIALS, HOWEVER USED.
C
C     IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA
C     BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT
C     LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND,
C     INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS,
C     REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE
C     REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY.
C
C     RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF
C     THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY
C     CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE
C     ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE.
C
C$ Required_Reading
C
C     PODS
C
C$ Keywords
C
C     ARRAY
C     CELLS
C     PODS
C
C$ Declarations

      INTEGER               LBCELL
      PARAMETER           ( LBCELL = -5 )

      DOUBLE PRECISION      ELEMS   (          * )
      INTEGER               N
      INTEGER               LOC
      DOUBLE PRECISION      POD     ( LBCELL : * )

C$ Brief_I/O
C
C     Variable  I/O  Description
C     --------  ---  --------------------------------------------------
C     ELEMS      I   New elements.
C     N          I   Number of new elements.
C     LOC        I   Location at which elements are to be inserted.
C     POD       I,O  Pod.
C
C$ Detailed_Input
C
C     ELEMS      contains elements to be inserted into the active
C                group of POD.
C
C     N          is the number of elements in ELEMS.
C
C     LOC        is the location (within the active group of the pod)
C                at which the new elements are to be inserted. The new
C                elements are inserted in front of the element currently
C                at this location.
C
C     POD        on input, is a pod.
C
C$ Detailed_Output
C
C     POD        on output, is the same pod, the active group of
C                which contains the new elements.
C
C$ Parameters
C
C     None.
C
C$ Exceptions
C
C      1) If N is not positive, the pod is not changed.
C
C     2) If there is insufficient room in the pod to insert all
C        of the new elements, the pod is not changed, and the error
C        SPICE(TOOMANYPEAS) is signalled.
C
C     3) If the location specified for location is not in the range
C        [1,NC+1], where NC is the number of elements in the active
C        group of the pod, the pod is not changed, and the error
C        SPICE(BADPODLOCATION) is signalled.
C
C$ Files
C
C     None.
C
C$ Particulars
C
C     This routine allows you to insert elements into the active
C     group of a pod without having to worry about checking for
C     overflow beforehand, or updating the cardinality afterwards.
C
C$ Examples
C
C     Elements can be inserted into the active group of a pod
C     by hand,
C
C        CALL PODOND ( POD, OFFSET, NUMBER )
C        END = OFFSET + NUMBER
C
C        CALL INSLAD ( ELEMS, N, OFFSET + LOC, POD(1), CUREND )
C        CALL SCARDD ( CUREND,                 POD            )
C
C     However, this is tedious, and it gets worse when you have to
C     check for possible overflow. PODIE accomplishes the same thing,
C
C        CALL PODIED ( ELEMS, N, LOC, POD )
C
C     more simply, and with error-handling built in.
C
C$ Restrictions
C
C     1) In any pod, only the active group should be accessed,
C        and its location should always be determined by PODBE
C        or PODON. Never assume that the active group begins
C        at POD(1).
C
C$ Literature_References
C
C     None.
C
C$ Author_and_Institution
C
C     W.L. Taber     (JPL)
C     I.M. Underwood (JPL)
C
C$ Version
C
C-    Beta Version 1.0.0, 15-JUL-1989 (WLT) (IMU)
C
C
C-&

C
C     SPICELIB functions
C
      LOGICAL               RETURN
      INTEGER               SIZED
      INTEGER               CARDD
C
C     Local variables
C
      INTEGER               NUMBER
      INTEGER               OFFSET
      INTEGER               END

C
C     Standard SPICE error handling.
C
      IF ( RETURN () ) THEN
         RETURN
      ELSE
         CALL CHKIN ( 'PODIED' )
      END IF

C
C     Three things can go `wrong':
C
C        1) No items to insert.
C
C        2) Too many items to insert.
C
C        3) No place to insert them.
C
      CALL PODOND ( POD, OFFSET, NUMBER )

      IF ( N .LT. 1 ) THEN

         CALL CHKOUT ( 'PODIED' )
         RETURN

      ELSE IF ( CARDD(POD) + N .GT. SIZED(POD) ) THEN

         CALL SETMSG ( 'Cannot fit # elements into # spaces.' )
         CALL ERRINT ( '#', N                                 )
         CALL ERRINT ( '#', SIZED(POD) - CARDD(POD)           )
         CALL SIGERR ( 'SPICE(TOOMANYPEAS)'                   )
         CALL CHKOUT ( 'PODIED'                               )
         RETURN

      ELSE IF ( LOC .LT. 1  .OR.  LOC .GT. NUMBER + 1 ) THEN

         CALL SETMSG ( 'Location (#) must be in the range [1,#].' )
         CALL ERRINT ( '#',  LOC                                  )
         CALL ERRINT ( '#',  NUMBER + 1                           )
         CALL SIGERR ( 'SPICE(BADPODLOCATION)'                    )
         CALL CHKOUT ( 'PODIED'                                   )
         RETURN

      END IF

C
C     In theory, we are home free. The rest looks just like the
C     code in $Examples, above.
C
      END = OFFSET + NUMBER

      CALL INSLAD ( ELEMS, N, OFFSET + LOC, POD(1), END )
      CALL SCARDD ( END,                    POD         )

      CALL CHKOUT ( 'PODIED' )
      RETURN
      END

max_stars_repo_name: agforero/FTFramework

max_issues_repo_name: agforero/fortran-testing-framework

max_forks_repo_name: agforero/fortran-testing-framework