text:
      SUBROUTINE AB08NW( EQUIL, N, M, P, A, LDA, B, LDB, C, LDC, D, LDD,
     $                   NFZ, NRANK, NIZ, DINFZ, NKROR, NINFE, NKROL,
     $                   INFZ, KRONR, INFE, KRONL, E, LDE, 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 extract from the system pencil
C
C                       ( A-lambda*I B )
C           S(lambda) = (              )
C                       (      C     D )
C
C     a regular pencil Af-lambda*Ef which has the finite Smith zeros of
C     S(lambda) as generalized eigenvalues. The routine also computes
C     the orders of the infinite Smith zeros and determines the singular
C     and infinite Kronecker structure of the system pencil, i.e., the
C     right and left Kronecker indices, and the multiplicities of the
C     infinite eigenvalues.
C
C     ARGUMENTS
C
C     Mode Parameters
C
C     EQUIL   CHARACTER*1
C             Specifies whether the user wishes to balance the system
C             matrix as follows:
C             = 'S':  Perform balancing (scaling);
C             = 'N':  Do not perform balancing.
C
C     Input/Output Parameters
C
C     N       (input) INTEGER.
C             The order of the square matrix A, the number of rows of
C             the matrix B, and number of columns of the matrix C.
C             N >= 0.
C
C     M       (input) INTEGER.
C             The number of columns of the matrix B.  M >= 0.
C
C     P       (input) INTEGER.
C             The number of rows of the matrix C.  P >= 0.
C
C     A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
C             On entry, the leading N-by-N part of this array must
C             contain the state dynamics matrix A of the system.
C             On exit, the leading NFZ-by-NFZ part of this array
C             contains the matrix Af of the reduced pencil.
C
C     LDA     INTEGER
C             The leading dimension of the array A.  LDA >= MAX(1,N).
C
C     B       (input/output) DOUBLE PRECISION array, dimension (LDB,M)
C             On entry, the leading N-by-M part of this array must
C             contain the input/state matrix B of the system.
C             On exit, this matrix does not contain useful information.
C
C     LDB     INTEGER
C             The leading dimension of the array B.  LDB >= 1, and
C             LDB >= MAX(1,N), 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 of the system.
C             On exit, this matrix does not contain useful information.
C
C     LDC     INTEGER
C             The leading dimension of the array C.  LDC >= MAX(1,P).
C
C     D       (input) DOUBLE PRECISION array, dimension (LDD,M)
C             The leading P-by-M part of this array must contain the
C             direct transmission matrix D of the system.
C
C     LDD     INTEGER
C             The leading dimension of the array D.  LDD >= MAX(1,P).
C
C     NFZ     (output) INTEGER
C             The number of finite zeros.
C
C     NRANK   (output) INTEGER
C             The normal rank of the system pencil.
C
C     NIZ     (output) INTEGER
C             The number of infinite zeros.
C
C     DINFZ   (output) INTEGER
C             The maximal multiplicity of infinite Smith zeros.
C
C     NKROR   (output) INTEGER
C             The number of right Kronecker indices.
C
C     NINFE   (output) INTEGER
C             The number of elementary infinite blocks.
C
C     NKROL   (output) INTEGER
C             The number of left Kronecker indices.
C
C     INFZ    (output) INTEGER array, dimension (N+1)
C             The leading DINFZ elements of INFZ contain information
C             on the infinite elementary divisors as follows:
C             the system has INFZ(i) infinite elementary divisors in
C             the Smith form of degree i, where i = 1,2,...,DINFZ.
C
C     KRONR   (output) INTEGER array, dimension (N+1)
C             The leading NKROR elements of this array contain the
C             right Kronecker (column) indices.
C
C     INFE    (output) INTEGER array, dimension (N+1)
C             The leading NINFE elements of INFE contain the
C             multiplicities of infinite eigenvalues.
C
C     KRONL   (output) INTEGER array, dimension (N+1)
C             The leading NKROL elements of this array contain the
C             left Kronecker (row) indices.
C
C     E       (output) DOUBLE PRECISION array, dimension (LDE,N)
C             The leading NFZ-by-NFZ part of this array contains the
C             matrix Ef of the reduced pencil.
C
C     LDE     INTEGER
C             The leading dimension of the array E.  LDE >= MAX(1,N).
C
C     Tolerances
C
C     TOL     DOUBLE PRECISION
C             A tolerance used in rank decisions to determine the
C             effective rank, which is defined as the order of the
C             largest leading (or trailing) triangular submatrix in the
C             QR (or RQ) factorization with column (or row) pivoting
C             whose estimated condition number is less than 1/TOL.
C             If the user sets TOL <= 0, then an implicitly computed,
C             default tolerance  TOLDEF = MAX(N+P,N+M)**2*EPS,  is used
C             instead, where EPS is the machine precision (see LAPACK
C             Library routine DLAMCH).  TOL < 1.
C
C     Workspace
C
C     IWORK   INTEGER array, dimension (MAX(M,P))
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 >= 1, if MAX(N,M,P) = 0; otherwise,
C             LDWORK >= MAX( MIN(P,M) + M + MAX(2*M,N) - 1,
C                            MIN(P,N) + MAX(N + MAX(P,M), 3*P - 1 ) ) +
C                            MAX(P+N,M+N)*MAX(P+N,M+N).
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 extracts from the system matrix of a state space
C     system, (A-lambda*I,B,C,D), a regular pencil Af-lambda*Ef, which
C     has the finite zeros of the system as generalized eigenvalues.
C     The procedure has the following main computational steps:
C
C        (a) construct the (N+P)-by-(M+N) system pencil
C
C             S(lambda) = (B  A)-lambda*( 0  I );
C                         (D  C)        ( 0  0 )
C
C        (b) reduce S(lambda) to S1(lambda) with the same finite zeros
C            and right Kronecker structure, but with D of full row rank;
C
C        (c) reduce the pencil S1(lambda) to S2(lambda) with the same
C            finite zeros and with D square invertible;
C
C        (d) perform a unitary transformation on the columns of
C            S2(lambda) = (A-lambda*I   B), in order to reduce it to
C                         (     C       D)
C
C            (Af-lambda*Ef   X), with Y and Ef square invertible;
C            (     0         Y)
C
C        (e) compute the right and left Kronecker indices of the system
C            matrix, which, together with the multiplicities of the
C            finite and infinite eigenvalues, constitute the complete
C            set of structural invariants under strict equivalence
C            transformations of a linear system.
C
C     REFERENCES
C
C     [1] Svaricek, F.
C         Computation of the Structural Invariants of Linear
C         Multivariable Systems with an Extended Version of the
C         Program ZEROS.
C         System & Control Letters, 6, pp. 261-266, 1985.
C
C     [2] Emami-Naeini, A. and Van Dooren, P.
C         Computation of Zeros of Linear Multivariable Systems.
C         Automatica, 18, pp. 415-430, 1982.
C
C     NUMERICAL ASPECTS
C
C     The algorithm is backward stable (see [2] and [1]).
C
C     FURTHER COMMENTS
C
C     In order to compute the finite Smith zeros of the system
C     explicitly, a call to this routine may be followed by a call to
C     the LAPACK Library routine DGGEV.
C
C     CONTRIBUTOR
C
C     A. Varga, German Aerospace Center, DLR Oberpfaffenhofen.
C     V. Sima, Katholieke Univ. Leuven, Belgium, May 1999.
C
C     REVISIONS
C
C     A. Varga, May 2003, German Aerospace Center, DLR Oberpfaffenhofen.
C     V. Sima, Dec. 2016, Jan. 2017, Apr. 2017.
C
C     KEYWORDS
C
C     Generalized eigenvalue problem, Kronecker indices, multivariable
C     system, orthogonal transformation, structural invariant.
C
C     ******************************************************************
C
C     .. Parameters ..
      DOUBLE PRECISION  ONE, ZERO
      PARAMETER         ( ONE = 1.0D0, ZERO = 0.0D0 )
C     .. Scalar Arguments ..
      CHARACTER         EQUIL
      INTEGER           DINFZ, INFO, LDA, LDB, LDC, LDD, LDE, LDWORK,
     $                  M, N, NFZ, NINFE, NIZ, NKROL, NKROR, NRANK, P
      DOUBLE PRECISION  TOL
C     .. Array Arguments ..
      INTEGER           INFE(*),  INFZ(*),  IWORK(*), KRONL(*), KRONR(*)
      DOUBLE PRECISION  A(LDA,*), B(LDB,*), C(LDC,*), D(LDD,*),
     $                  DWORK(*), E(LDE,*)
C     .. Local Scalars ..
      LOGICAL           LEQUIL, LQUERY, QRET
      INTEGER           I, I0, I1, II, ITAU, J, JWORK, KABCD, LABCD2,
     $                  LDABCD, MM, MPM, MPN, MU, NN, NSINFE, NU, NU1,
     $                  PP, WRKOPT
      DOUBLE PRECISION  MAXRED, SVLMAX, TOLER
C     .. External Functions ..
      LOGICAL           LSAME
      DOUBLE PRECISION  DLAMCH, DLANGE
      EXTERNAL          DLAMCH, DLANGE, LSAME
C     .. External Subroutines ..
      EXTERNAL          AB08NY, DLACPY, DLASET, DORMRZ, DTZRZF, MA02BD,
     $                  TB01ID, TB01XD, XERBLA
C     .. Intrinsic Functions ..
      INTRINSIC         DBLE, INT, MAX, MIN
C     .. Executable Statements ..
C
      INFO   = 0
      LDABCD = N + MAX( P, M )
      LABCD2 = LDABCD*LDABCD
      LEQUIL = LSAME( EQUIL, 'S' )
C
C     Test the input scalar arguments.
C
      IF( .NOT.LEQUIL .AND. .NOT.LSAME( EQUIL, 'N' ) ) THEN
         INFO = -1
      ELSE IF( N.LT.0 ) THEN
         INFO = -2
      ELSE IF( M.LT.0 ) THEN
         INFO = -3
      ELSE IF( P.LT.0 ) THEN
         INFO = -4
      ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
         INFO = -6
      ELSE IF( LDB.LT.1 .OR. ( M.GT.0 .AND. LDB.LT.N ) ) THEN
         INFO = -8
      ELSE IF( LDC.LT.MAX( 1, P ) ) THEN
         INFO = -10
      ELSE IF( LDD.LT.MAX( 1, P ) ) THEN
         INFO = -12
      ELSE IF( LDE.LT.MAX( 1, N ) ) THEN
         INFO = -25
      ELSE IF( TOL.GE.ONE ) THEN
         INFO = -26
      ELSE
         MPN  = MIN( P, N )
         MPM  = MIN( P, M )
         QRET = MAX( N, M, P ).EQ.0
         LQUERY = ( LDWORK.EQ.-1 )
         IF( QRET ) THEN
            JWORK = 1
         ELSE
            JWORK = MAX( MPM + M + MAX( 2*M, N ) - 1,
     $                   MPN + MAX( LDABCD, 3*P  - 1 ) ) + LABCD2
         END IF
         IF( LQUERY ) THEN
            IF( QRET ) THEN
               WRKOPT = 1
            ELSE
               SVLMAX = ZERO
               NIZ    = 0
               CALL AB08NY( .TRUE., N, M, P, SVLMAX, DWORK, LDABCD, NIZ,
     $                      NU, MU, DINFZ, NKROL, INFZ, KRONL, TOL,
     $                      IWORK, DWORK, -1, INFO )
               WRKOPT = MAX( JWORK, LABCD2 + INT( DWORK(1) ) )
               CALL AB08NY( .FALSE., N, M, M, SVLMAX, DWORK, LDABCD,
     $                      NIZ, NU, MU, I1, NKROR, IWORK, KRONR, TOL,
     $                      IWORK, DWORK, -1, INFO )
               WRKOPT = MAX( WRKOPT, LABCD2 + INT( DWORK(1) ) )
               CALL DTZRZF( MPM, N+MPM, DWORK, LDABCD, DWORK, DWORK, -1,
     $                      INFO )
               WRKOPT = MAX( WRKOPT, LABCD2 + MPM + INT( DWORK(1) ) )
               CALL DORMRZ( 'Right', 'Transpose', N, N+MPM, MPM, N,
     $                      DWORK, LDABCD, DWORK, DWORK, LDABCD, DWORK,
     $                      -1, INFO )
               WRKOPT = MAX( WRKOPT, LABCD2 + MPM + INT( DWORK(1) ) )
            END IF
         ELSE IF( LDWORK.LT.JWORK ) THEN
            INFO = -29
         END IF
      END IF
C
      IF( INFO.NE.0 ) THEN
C
C        Error return.
C
         CALL XERBLA( 'AB08NW', -INFO )
         RETURN
      ELSE IF( LQUERY ) THEN
         DWORK(1) = WRKOPT
         RETURN
      END IF
C
      NIZ   = 0
      NKROL = 0
      NKROR = 0
      NINFE = 0
C
C     Quick return if possible.
C
      IF( QRET ) THEN
         DINFZ = 0
         NFZ   = 0
         NRANK = 0
         DWORK(1) = ONE
         RETURN
      END IF
C
C     (Note: Comments in the code beginning "Workspace:" describe the
C     minimal amount of real workspace needed at that point in the
C     code, as well as the preferred amount for good performance.)
C
      WRKOPT = 1
      KABCD  = 1
      JWORK  = KABCD + LABCD2
C
C     If required, balance the system pencil.
C     Workspace: need   N.
C
      IF( LEQUIL ) THEN
         MAXRED = ZERO
         CALL TB01ID( 'A', N, M, P, MAXRED, A, LDA, B, LDB, C, LDC,
     $                DWORK, INFO )
         WRKOPT = N
      END IF
C
C     Construct the system pencil
C
C                      ( B A-lambda*I )
C         S(lambda) =  (              )
C                      ( D     C      )
C
C     of dimension (N+P)-by-(M+N).
C
      NN = N
      MM = M
      PP = P
C
      CALL DLACPY( 'Full', N, M, B, LDB, DWORK(KABCD),          LDABCD )
      CALL DLACPY( 'Full', P, M, D, LDD, DWORK(KABCD+N),        LDABCD )
      CALL DLACPY( 'Full', N, N, A, LDA, DWORK(KABCD+LDABCD*M), LDABCD )
      CALL DLACPY( 'Full', P, N, C, LDC, DWORK(KABCD+LDABCD*M+N),
     $              LDABCD )
C
C     If required, set tolerance.
C
      TOLER = TOL
      IF( TOLER.LE.ZERO ) THEN
         TOLER = DBLE( LABCD2 ) * DLAMCH( 'Precision' )
      END IF
      SVLMAX = DLANGE( 'Frobenius', NN+PP, NN+MM, DWORK(KABCD), LDABCD,
     $                 DWORK(JWORK) )
C
C     Extract the reduced pencil S1(lambda)
C
C             ( Bc  Ac-lambda*I ),
C             ( Dc      Cc      )
C
C     having the same finite Smith zeros as the system pencil S(lambda),
C     but with Dc, a MU-by-MM full row rank left upper-trapezoidal
C     matrix, with the first MU columns in an upper triangular form.
C
C     Workspace: need   MAX( MIN(P,M) + M + MAX(2*M,N) - 1,
C                            MIN(P,N) + MAX(N + MAX(P,M), 3*P - 1 ) ).
C                prefer larger.
C     Int.work.  need   MAX(M,P).
C
      CALL AB08NY( .TRUE., NN, MM, PP, SVLMAX, DWORK(KABCD), LDABCD,
     $             NIZ, NU, MU, DINFZ, NKROL, INFZ, KRONL, TOLER, IWORK,
     $             DWORK(JWORK), LDWORK-JWORK+1, INFO )
C
      WRKOPT = MAX( WRKOPT, INT( DWORK(JWORK) ) + JWORK - 1 )
C
C     Set the number of simple (non-dynamic) infinite eigenvalues, and
C     the normal rank of the system pencil.
C
      NSINFE = MU
      NRANK  = NN + MU
C
C     Pertranspose the system.
C
      CALL TB01XD( 'D', NU, MM, MM, MAX( 0, NU-1 ), MAX( 0, NU-1 ),
     $              DWORK(KABCD+LDABCD*MM),    LDABCD,
     $              DWORK(KABCD),              LDABCD,
     $              DWORK(KABCD+LDABCD*MM+NU), LDABCD,
     $              DWORK(KABCD+NU),           LDABCD, INFO )
      CALL MA02BD( 'Right', NU+MM, MM, DWORK(KABCD),    LDABCD )
      CALL MA02BD( 'Left',  MM, NU+MM, DWORK(KABCD+NU), LDABCD )
C
      IF( MU.NE.MM ) THEN
         NN = NU
         PP = MM
         MM = MU
         KABCD = KABCD + ( PP - MM )*LDABCD
C
C        Extract the reduced pencil S2(lambda) to
C
C             ( Br  Ar-lambda*I ),
C             ( Dr      Cr      )
C
C        having the same finite Smith zeros as the pencil S(lambda),
C        but with Dr, a MU-by-MU invertible upper-triangular matrix.
C
C        Workspace: need  MAX( MIN(P,M) + M + MAX(2*M,N) - 1,
C                              MIN(P,N) + MAX(N + MAX(P,M), 3*P-1 ) ).
C                   prefer larger.
C
         CALL AB08NY( .FALSE., NN, MM, PP, SVLMAX, DWORK(KABCD), LDABCD,
     $                I0, NU, MU, I1, NKROR, IWORK, KRONR, TOLER,
     $                IWORK, DWORK(JWORK), LDWORK-JWORK+1, INFO )
C
         WRKOPT = MAX( WRKOPT, INT( DWORK(JWORK) ) + JWORK - 1 )
      END IF
C
      IF( MIN( NU, MU ).NE.0 ) THEN
C
C        Perform a unitary transformation on the columns of
C                     ( Br Ar-lambda*I  ),
C                     ( Dr     Cr       )
C        in order to reduce it to
C                     ( Af-lambda*Ef  * ),
C                     (      0        Y )
C        with Y and Ef square invertible.
C
         NU1   = NU + KABCD
         I1    = NU + MU
         ITAU  = JWORK
         JWORK = ITAU + MU
C
C        Workspace: need   2*MIN(P,M);
C                   prefer MIN(P,M) + MIN(P,M)*NB.
C
         CALL DTZRZF( MU, I1, DWORK(NU1), LDABCD, DWORK(ITAU),
     $                DWORK(JWORK), LDWORK-JWORK+1, INFO )
         WRKOPT = MAX( WRKOPT, INT( DWORK(JWORK) ) + JWORK - 1 )
C
C        Compute and save Af.
C
C        Workspace: need   MIN(P,M) + N;
C                   prefer MIN(P,M) + N*NB.
C
         CALL DORMRZ( 'Right', 'Transpose', NU, I1, MU, NU,
     $                DWORK(NU1), LDABCD, DWORK(ITAU), DWORK(KABCD),
     $                LDABCD, DWORK(JWORK), LDWORK-JWORK+1, INFO )
         WRKOPT = MAX( WRKOPT, INT( DWORK(JWORK) ) + JWORK - 1 )
C
         CALL DLACPY( 'Full', NU, NU, DWORK(KABCD+MU*LDABCD), LDABCD, A,
     $                LDA )
C
C        Compute and save Ef.
C
         CALL DLASET( 'Full', NU, MU, ZERO, ZERO, DWORK(KABCD), LDABCD )
         CALL DLASET( 'Full', NU, NU, ZERO, ONE, DWORK(KABCD+MU*LDABCD),
     $                LDABCD )
         CALL DORMRZ( 'Right', 'Transpose', NU, I1, MU, NU,
     $                DWORK(NU1), LDABCD, DWORK(ITAU), DWORK(KABCD),
     $                LDABCD, DWORK(JWORK), LDWORK-JWORK+1, INFO )
         WRKOPT = MAX( WRKOPT, INT( DWORK(JWORK) ) + JWORK - 1 )
C
         CALL DLACPY( 'Full', NU, NU, DWORK(KABCD+MU*LDABCD), LDABCD, E,
     $                LDE )
      ELSE
C
C        Save Af and set Ef.
C
         CALL DLACPY( 'Full', NU, NU, DWORK(KABCD+MU*LDABCD), LDABCD, A,
     $                LDA )
         CALL DLASET( 'Full', NU, NU, ZERO, ONE, E, LDE )
      END IF
C
      NFZ = NU
C
C     Set right Kronecker indices (column indices).
C
      DO 30 I = 1, NKROR
         IWORK(I) = KRONR(I)
   30 CONTINUE
C
      J = 0
C
      DO 50 I = 1, NKROR
         DO 40 II = J + 1, J + IWORK(I)
            KRONR(II) = I - 1
   40    CONTINUE
C
         J = J + IWORK(I)
   50 CONTINUE
C
      NKROR = J
C
C     Set left Kronecker indices (row indices).
C
      DO 60 I = 1, NKROL
         IWORK(I) = KRONL(I)
   60 CONTINUE
C
      J = 0
C
      DO 80 I = 1, NKROL
         DO 70 II = J + 1, J + IWORK(I)
            KRONL(II) = I - 1
   70    CONTINUE
C
         J = J + IWORK(I)
   80 CONTINUE
C
      NKROL = J
C
C     Determine the number of simple infinite blocks as the difference
C     between the order of Dr and the number of infinite blocks of order
C     greater than one.
C
      DO 90 I = 1, DINFZ
         NINFE = NINFE + INFZ(I)
   90 CONTINUE
C
      NINFE = NSINFE - NINFE
C
      DO 100 I = 1, NINFE
         INFE(I) = 1
  100 CONTINUE
C
C     Set the structure of infinite eigenvalues.
C
      DO 120 I = 1, DINFZ
C
         DO 110 II = NINFE + 1, NINFE + INFZ(I)
            INFE(II) = I + 1
  110    CONTINUE
C
         NINFE = NINFE + INFZ(I)
  120 CONTINUE
C
      DWORK(1) = WRKOPT
      RETURN
C *** Last line of AB08NW ***
      END

max_stars_repo_name: bnavigator/SLICOT-Reference

max_issues_repo_name: RJHKnight/slicotr

max_forks_repo_name: RJHKnight/slicotr


text:
        subroutine AsBMFG (y,       x,       shpb,    shglb,
     &                     ienb,    materb,  iBCB,    BCB,
     &                     res,     rmes,    EGmass)
c
c----------------------------------------------------------------------
c
c This routine computes and assembles the data corresponding to the
c  boundary elements.
c
c Zdenek Johan, Winter 1991.  (Fortran 90)
c----------------------------------------------------------------------
c
        include "common.h"
c
        dimension y(nshg,ndofl),           x(numnp,nsd),
     &            shpb(nshl,ngaussb),
     &            shglb(nsd,nshl,ngaussb),
     &            ienb(npro,nshl),          materb(npro),
     &            iBCB(npro,ndiBCB),        BCB(npro,nshlb,ndBCB),
     &            res(nshg,nflow),         rmes(nshg,nflow)
c
        dimension ycl(npro,nshl,ndofl),  xlb(npro,nenl,nsd),
     &            rl(npro,nshl,nflow),
     &            rml(npro,nshl,nflow),
     &            EGmass(npro, nshl, nshl)
c
        dimension sgn(npro,nshl)
c
c.... create the matrix of mode signs for the hierarchic basis
c     functions.
c
        if (ipord .gt. 1) then
           call getsgn(ienb,sgn)
        endif
c
c.... gather the variables
c

        call localy(y,      ycl,     ienb,   ndofl,  'gather  ')
        call localx(x,      xlb,    ienb,   nsd,    'gather  ')
c

        !get the boundary element residuals

        rl  = zero
        rml = zero
c
!  pass the memory location of ycl to both yl and ycl in e3b.  This may
 !  seem dangerous since yl in e3b is :,nflow and ycl is :,ndof but they
 !  do not write to yl (out of bounds at least), only use the data there
 !  so both will access data
 !  properly from this location.
c
        call e3b  (ycl,     ycl,     iBCB,    BCB,     shpb,    shglb,
     &             xlb,     rl,      rml,     sgn,     EGmass)

        !assemble the residual and the modified residual
        call local(res,    rl,     ienb,   nflow,  'scatter ')
        if (Navier .eq. 1 .and. ires.ne.1 )
     &    call local(rmes,   rml,    ienb,   nflow,  'scatter ')

        !end
        return
        end
c
c
c
        subroutine AsBMFGSclr (y,       x,       shpb,    shglb,
     &                         ienb,    materb,  iBCB,
     &                         BCB,     rest,    rmest)
c
c----------------------------------------------------------------------
c
c This routine computes and assembles the data corresponding to the
c  boundary elements.
c
c Zdenek Johan, Winter 1991.  (Fortran 90)
c----------------------------------------------------------------------
c
        include "common.h"
c
        dimension y(nshg,ndofl),           x(numnp,nsd),
     &            shpb(nshl,maxsh),
     &            shglb(nsd,nshl,maxsh),
     &            ienb(npro,nshl),       materb(npro),
     &            iBCB(npro,ndiBCB),   BCB(npro,nshlb,ndBCB),
     &            rest(nshg),         rmest(nshg)
c
        dimension ycl(npro,nshl,ndofl),   xlb(npro,nenl,nsd),
     &            rtl(npro,nshl),
     &            rmtl(npro,nshl)
c
c.... create the matrix of mode signs for the hierarchic basis
c     functions.
c
        if (ipord .gt. 1) then
           call getsgn(ienb,sgn)
        endif
c
c.... gather the variables
c
        call localy (y,      ycl,     ienb,   ndofl,  'gather  ')
        call localx (x,      xlb,    ienb,   nsd,    'gather  ')
c
c.... get the boundary element residuals
c
        rtl  = zero
        rmtl = zero
c
c.... 3D
c
            call e3bSclr (ycl,    iBCB,    BCB,
     &                    shpb,  shglb,   sgn,
     &                    xlb,   rtl,     rmtl)
c
c.... assemble the residual and the modified residual
c

        call local (rest,    rtl,     ienb,   1,  'scatter ')

c
        if (Navier .eq. 1)
     &  call local (rmest,   rmtl,    ienb,   1,  'scatter ')
c
c.... end
c
        return
        end

ext: f

max_stars_repo_name: polmes/phasta

max_issues_repo_name: polmes/phasta

max_forks_repo_name: polmes/phasta


text:
      subroutine AMACH(MODE, I, I1, R1, D1)
c Copyright (c) 1996 California Institute of Technology, Pasadena, CA.
c ALL RIGHTS RESERVED.
c Based on Government Sponsored Research NAS7-03001.
c>> 2003-04-25 AMACH  Krogh   Added comment for the Absoft compiler.
c>> 2001-03-13 AMACH  Krogh   Fixed HOW for case when not set in m77job.
c>> 1997-04-16 AMACH  Krogh   Remove blank lines -- was confusing m77con
c>> 1996-03-30 AMACH  Krogh   Added external statement.
c>> 1994-10-26 AMACH  Krogh   Changes to use M77CON
c>> 1994-09-23 AMACH  Snyder  Add VAX G parameters
c>> 1994-06-21 AMACH  Snyder  Compute only round-off and u-flow at first
c>> 1994-05-25 AMACH  Snyder  Added an option to compute at run time.
c>> 1992-04-07 AMACH  Oken    Removed ^Z at EOF (error found by VAX comp
c>> 1992-02-20 AMACH  Snyder  Added Cray-J90 stuff, q.v.
c>> 1990-06-11 AMACH  Snyder  Added Apollo DN-10000 stuff, q.v.
c>> 1990-12-14 AMACH  Lawson  Changed to eliminate ENTRY statements.
c>> 1990-08-21 AMACH  Krogh   No test was getting done for bad machine.
c>> 1990-02-28 AMACH  Krogh   Correct missing DOUBLE PRECISION AMSUB1
c>> 1989-08-14 AMACH  Krogh   Parameterized everything -- Massive change
c>> 1989-03-30 AMACH  Snyder  Correct missing "/" line 921
c>> 1989-01-30 AMACH  Snyder  Incorporate more constants from NETLIB.
C>> 1988-05-19 AMACH  Lawson  Initial code.
c File AMACH.FOR contains user-callable functions I1MACH, D1MACH, and
c R1MACH, plus second-level subroutines AMACH, AMTEST, and AMSUB1.
c Appropriate lines must be switched between comment and non-comment
c status when this code is moved to a different computer system.
c     These changes can be done with any text editor, however the "c++"
c lines permit automation of the change using the M77CON processor.
c Note that when the M77CON processor activates a line it shifts
c Columns 2-72 to 1-71 and puts a blank in Column 72.  When it inactiv-
c ates a line it shifts Columns 1-71 to 2-72 and puts a C in Column 1.
c     The possible choices using M77CON (don't include parenthetical
c     comments) are:
c      c++ CURRENT HAS SYS = IEEE
c      c++ CURRENT HAS SYS = ALPHA_D3
c      c++ CURRENT HAS SYS = AMDAHL
c      c++ CURRENT HAS SYS = APOLLO_10000
c      c++ CURRENT HAS SYS = BUR1700
c      c++ CURRENT HAS SYS = BUR5700
c      c++ CURRENT HAS SYS = BUR67_7700
c      c++ CURRENT HAS SYS = CDC60_7000
c      c++ CURRENT HAS SYS = CONVEXC_1
c      c++ CURRENT HAS SYS = CRAY1
c      c++ CURRENT HAS SYS = CRAY1_SD (Sngl prec.arith. used for dble.)
c      c++ CURRENT HAS SYS = CRAY1_64 (64 bit integers)
c      c++ CURRENT HAS SYS = CRAY1_SD_64 (64 bit int, SP used for DP)
c      c++ CURRENT HAS SYS = CRAY_T3D
c      c++ CURRENT HAS SYS = CRAY_J90
c      c++ CURRENT HAS SYS = CRAY_J90_SD (Sngl prec. used for dble.)
c      c++ CURRENT HAS SYS = DG_S2000
c      c++ CURRENT HAS SYS = HARRIS220
c      c++ CURRENT HAS SYS = HON600_6000
c      c++ CURRENT HAS SYS = HON_DPS_8_70
c      c++ CURRENT HAS SYS = HP700Q
c      c++ CURRENT HAS SYS = IBM360_370
c      c++ CURRENT HAS SYS = INTERDATA_8_32
c      c++ CURRENT HAS SYS = PDP10_KA
c      c++ CURRENT HAS SYS = PDP10_KB
c      c++ CURRENT HAS SYS = PDP11
c      c++ CURRENT HAS SYS = PRIME50
c      c++ CURRENT HAS SYS = SEQ_BAL_8000
c      c++ CURRENT HAS SYS = UNIVAC
c      c++ CURRENT HAS SYS = VAX
c      c++ CURRENT HAS SYS = VAX_G
c     The current choice is:
c++ CURRENT HAS SYS = IEEE
c
c     One can also select whether floating point constants are created
c     by the compiler or created at run time.  The choices using M77CON
c     are:
c      c++ CURRENT HAS HOW = COMPILER
c      c++ CURRENT HAS HOW = RUN
c
c If this is not set in m77job, it is set as follows:
c++((SYS==VAX)|(SYS==ALPHA_D3)) Default HOW = RUN
c++ Default HOW = COMPILER
c
c     The current choice is:
c++ CURRENT HAS HOW = COMPILER
c
c     If the constants are created at run time, and they fail the run-
c     time check for reasonableness, they are re-created assuming IEEE.
c     If they still fail, the program stops.
c
C  I/O UNIT NUMBERS:
C
C    IM1 = I1MACH( 1) = THE STANDARD INPUT UNIT.
C    IM2 = I1MACH( 2) = THE STANDARD OUTPUT UNIT.
C    IM3 = I1MACH( 3) = THE STANDARD PUNCH UNIT.
C    IM4 = I1MACH( 4) = THE STANDARD ERROR MESSAGE UNIT.
C
C  WORDS:
C
C    IM5 = I1MACH( 5) = THE NUMBER OF BITS PER INTEGER STORAGE UNIT.
C    IM6 = I1MACH( 6) = THE NUMBER OF CHARACTERS/INTEGER STORAGE UNIT.
C
C  INTEGERS:
C
C    ASSUME INTEGERS ARE REPRESENTED IN THE S-DIGIT, BASE-A FORM
C
C               SIGN ( X(S-1)*A**(S-1) + ... + X(1)*A + X(0) )
C
C               WHERE 0 .LE. X(I) .LT. A FOR I=0,...,S-1.
C
C    IM7 = I1MACH( 7) = A, THE BASE.
C    IM8 = I1MACH( 8) = S, THE NUMBER OF BASE-A DIGITS.
C    IM9 = I1MACH( 9) = A**S - 1, THE LARGEST MAGNITUDE.
C
C  FLOATING-POINT NUMBERS:
C
C    ASSUME FLOATING-POINT NUMBERS ARE REPRESENTED IN THE T-DIGIT,
C    BASE-B FORM
C
C               SIGN (B**E)*( (X(1)/B) + ... + (X(T)/B**T) )
C
C               WHERE 0 .LE. X(I) .LT. B FOR I=1,...,T,
C               0 .LT. X(1), AND EMIN .LE. E .LE. EMAX.
C
C    IM10 = I1MACH(10) = B, THE BASE.
C
C  SINGLE-PRECISION:
C
C    IM11 = I1MACH(11) = T, THE NUMBER OF BASE-B DIGITS.
C    IM12 = I1MACH(12) = EMIN, THE SMALLEST EXPONENT E.
C    IM13 = I1MACH(13) = EMAX, THE LARGEST EXPONENT E.
C
C  DOUBLE-PRECISION:
C
C    IM14 = I1MACH(14) = T, THE NUMBER OF BASE-B DIGITS.
C    IM15 = I1MACH(15) = EMIN, THE SMALLEST EXPONENT E.
C    IM16 = I1MACH(16) = EMAX, THE LARGEST EXPONENT E.
C
C  CONVERSION FROM FUNCTIONAL TO STRUCTURAL FLOATING POINT CONSTANTS
C
C    IM17 = CONSTANT SUCH THAT IM14 + IM17 = ACTUAL NUMBER OF BASE-B
C           DIGITS IN DOUBLE PRECISION, USED FOR CHECKING THAT CORRECT
C           VERSION OF THIS PROGRAM IS INSTALLED.  (SEE DEFINITION OF
C           DM6, AND THE USE OF DM6 IN CALLING AMTEST.)
C
C  TO ALTER THIS FUNCTION FOR A PARTICULAR ENVIRONMENT,
C  THE DESIRED SET OF PARAMETER STATEMENTS SHOULD BE ACTIVATED BY
C  REMOVING THE C FROM COLUMN 1.  ALSO, THE VALUES OF
C  IM1 - IM4 SHOULD BE CHECKED FOR CONSISTENCY
C  WITH THE LOCAL OPERATING SYSTEM.
c     -----------------------------------------------------------------
c     Original design and code due to P. A. Fox, A. D. Hall, and
c     N. L. Schryer, Bell Laboratories.  See ACM TOMS, 4,(1978),177-188.
c     Adapted to Univac 1100 by Kris Stewart, JPL, 7/30/81.
c     Adapted for the JPL MATH77 library by C. L. Lawson and F. T. Krogh
c     Sept, 1987.
c     1989-08-14 AMACH  Krogh   Parameterized everything. Major changes.
C     1990 Dec. CLL reorganized code to avoid using ENTRY statements
c     for functions of different types.  Also added save statements.
c     -----------------------------------------------------------------
c     On the first call to this function, tests are done to verify that
c     IM10 and IM14 are not grossly wrong for the host environment.
c     This gives some protection against using the wrong version of this
c     subprogram.
c     -----------------------------------------------------------------
      integer MODE, I, I1
      real R1
      double precision D1, TEST
c
      integer IMACH(17)
      integer IM1, IM2, IM3, IM4, IM5, IM6, IM7, IM8, IM9, IM10, IM11,
     1            IM12, IM13, IM14, IM15, IM16, IM17
c++ Code for HOW=RUN is INACTIVE
C      integer IEEE
C      integer ID1, ID2, ID3, ID4, ID5, ID6, ID7, ID8, ID10, ID11,
C     1   ID12, ID13, ID14, ID15, ID16, ID17
c++ Code for (HOW=RUN) | SYS=IEEE is ACTIVE
      integer IE1, IE2, IE3, IE4, IE5, IE6, IE7, IE8, IE10, IE11,
     1   IE12, IE13, IE14, IE15, IE16, IE17
c++ end
      real             RMACH(5), RM1, RM2, RM3, RM4, RM5,
     1                 RMA, RMB, RBASE
      double precision DMACH(5), DM1, DM2, DM3, DM4, DM5, DM6,
     1                 DMA, DMB, DBASE
      save TEST, IMACH, RMACH, DMACH
C     -----------------------------------------------------------------
C     Machine constants for IEEE standard binary floating-point
c     processors.  This includes PC's and work-stations using the
c     Intel 8087, 80287, 80387, ... processors or the
c     Motorola 68881, 68882, ... processors.
c     Note:  We are setting the "most negative exponent" (IMACH(12) and
c     IMACH(15)) to be the exponent of the smallest normalized number.
c     An IEEE processor actually handles smaller numbers before
c     underflowing, however these "unnormalized" numbers have
c     diminished precision.
c
c++ Code for (HOW=RUN) | SYS=IEEE is ACTIVE
c     Parameters for IEEE when generating at run time:
      PARAMETER (IE1 =5, IE2 =6, IE3 =7, IE4 =6)
      PARAMETER (IE5 =32, IE6 =4, IE7 =2, IE8 =31)
      PARAMETER (IE10 =2, IE11 =24, IE12 =-125, IE13 =128)
      PARAMETER (IE14 =53, IE15 =-1021, IE16 =1024, IE17=0)
c++ Code for SYS = IEEE is ACTIVE
      PARAMETER (IM1 = IE1, IM2 = IE2, IM3 = IE3, IM4 = IE4)
      PARAMETER (IM5 = IE5, IM6 = IE6, IM7 = IE7, IM8 = IE8)
      PARAMETER (IM10 = IE10, IM11 = IE11, IM12 = IE12, IM13 = IE13)
      PARAMETER (IM14 = IE14, IM15 = IE15, IM16 = IE16, IM17 = IE17)
C     -----------------------------------------------------------------
c++ Code for SYS = ALPHA_D3 is INACTIVE
Cc     MACHINE CONSTANTS for the VAX/VMS F and D-3 format for Alpha
Cc
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
C      PARAMETER (IM10 =2, IM11 =24, IM12 =-127, IM13 =127)
C      PARAMETER (IM14 =53, IM15 =-127, IM16 =127, IM17=0)
c++ end
C     -----------------------------------------------------------------
c++ Code for HOW = RUN is INACTIVE
Cc     MACHINE CONSTANTS for the VAX/VMS F and D-3 format for Alpha
Cc
C      PARAMETER (ID1 =5, ID2 =6, ID3 =7, ID4 =6)
C      PARAMETER (ID5 =32, ID6 =4, ID7 =2, ID8 =31)
C      PARAMETER (ID10 =2, ID11 =24, ID12 =-127, ID13 =127)
C      PARAMETER (ID14 =53, ID15 =-127, ID16 =127, ID17=0)
c++ end
C     -----------------------------------------------------------------
c++ Code for SYS = AMDAHL is INACTIVE
CC     MACHINE CONSTANTS FOR AMDAHL MACHINES.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
C      PARAMETER (IM10 =16, IM11 =6, IM12 =-64, IM13 =63)
C      PARAMETER (IM14 =14, IM15 =-64, IM16 =63, IM17=0)
C      -----------------------------------------------------------------
c++ Code for SYS = APOLLO_10000 is INACTIVE
cc     MACHINE CONSTANTS FOR APOLLO DN_10000 MACHINES.
cc     The only difference from IEEE is IM13.  This difference has
cc     nothing to do with the arithmetic or representation used by the
cc     machine.  It is caused by a bug in the compiler:  The right-hand
cc     side of RM2 (below) is apparently evaluated in double precision.
cc     When the compiler is ready to store the resulting value into its
cc     internal data structures, it compares it to an incorrect value
cc     of the overflow limit.  It appears the incorrect value has the
cc     correct exponent, but the fraction is 1.5 instead of 2-2**(-p),
cc     where p is the precision in bits.  You can get the correct result
cc     by changing IM13 to 128, changing RM2 from a parameter to a
cc     variable, and changing the parameter statement that assigns a
cc     value to RM2 into an ordinary assignment statement.
CC
c      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
c      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
c      PARAMETER (IM10 =2, IM11 =24, IM12 =-125, IM13 =127)
c      PARAMETER (IM14 =53, IM15 =-1021, IM16 =1024, IM17 =0)
CC     -----------------------------------------------------------------
c++ Code for SYS = BUR1700 is INACTIVE
CC     MACHINE CONSTANTS FOR THE BURROUGHS 1700 SYSTEM.
CC
C      PARAMETER (IM1 =7, IM2 =2, IM3 =2, IM4 =2)
C      PARAMETER (IM5 =36, IM6 =4, IM7 =2, IM8 =33)
C      PARAMETER (IM10 =2, IM11 =24, IM12 =-256, IM13 =255)
C      PARAMETER (IM14 =60, IM15 =-256, IM16 =255, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = BUR5700 is INACTIVE
CC     MACHINE CONSTANTS FOR THE BURROUGHS 5700 SYSTEM.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =48, IM6 =6, IM7 =2, IM8 =39)
C      PARAMETER (IM10 =8, IM11 =13, IM12 =-50, IM13 =76)
C      PARAMETER (IM14 =26, IM15 =-50, IM16 =76, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = BUR67_7700 is INACTIVE
CC     MACHINE CONSTANTS FOR THE BURROUGHS 6700/7700 SYSTEMS.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =48, IM6 =6, IM7 =2, IM8 =39)
C      PARAMETER (IM10 =8, IM11 =13, IM12 =-50, IM13 =76)
C      PARAMETER (IM14 =26, IM15 =-32754, IM16 =32780, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = CDC60_7000 is INACTIVE
CC     MACHINE CONSTANTS FOR THE CDC 6000/7000 SERIES.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =60, IM6 =10, IM7 =2, IM8 =48)
C      PARAMETER (IM10 =2, IM11 =47, IM12 =-929, IM13 =1070)
C      PARAMETER (IM14 =94, IM15 =-929, IM16 =1069, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = CONVEXC_1 is INACTIVE
CC     MACHINE CONSTANTS FOR CONVEX C-1.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
C      PARAMETER (IM10 =2, IM11 =24, IM12 =-128, IM13 =127)
C      PARAMETER (IM14 =53, IM15 =-1024, IM16 =1023, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = CRAY1 is INACTIVE
CC     MACHINE CONSTANTS FOR THE CRAY 1, XMP, 2, AND 3.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =102, IM4 =6)
C      PARAMETER (IM5 =64, IM6 =8, IM7 =2, IM8 =46)
C      PARAMETER (IM10 =2, IM11 =47, IM12 =-8189, IM13 =8190)
C      PARAMETER (IM14 =94, IM15 =-8099, IM16 =8190, IM17=2)
CC     -----------------------------------------------------------------
c++ Code for SYS = CRAY_T3D is INACTIVE
cc     Machine constants for Cray T3D.  IEEE double for both precisions.
c      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
c      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
c      PARAMETER (IM10 =2, IM11 =53, IM12 =-1021, IM13 =1024)
c      PARAMETER (IM14 =53, IM15 =-1021, IM16 =1024, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = CRAY_J90 is INACTIVE
CC     MACHINE CONSTANTS FOR THE CRAY J90
CC     Cray claims the overflow exponent (IM13 and IM16) is 8189, and
CC     the underflow exponent (IM12 and IM15) is -8189, but these values
CC     don't seem to work in cf77:  the underflow limit underflows, and
CC     the overflow limit overflows when using Cray's values.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =102, IM4 =6)
C      PARAMETER (IM5 =64, IM6 =8, IM7 =2, IM8 =46)
C      PARAMETER (IM10 =2, IM11 =47, IM12 =-8188, IM13 =8189)
C      PARAMETER (IM14 =94, IM15 =-8188, IM16 =8189, IM17=2)
CC     -----------------------------------------------------------------
c++ Code for SYS = CRAY_J90_SD is INACTIVE
CC     MACHINE CONSTANTS FOR THE CRAY J90
CC     Cray claims the overflow exponent (IM13 and IM16) is 8189, and
CC     the underflow exponent (IM12 and IM15) is -8189, but these
CC     values don't seem to work in cf77:  the underflow limit under-
CC     flows, and the overflow limit overflows when using Cray's values.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =102, IM4 =6)
C      PARAMETER (IM5 =64, IM6 =8, IM7 =2, IM8 =46)
C      PARAMETER (IM10 =2, IM11 =47, IM12 =-8188, IM13 =8189)
C      PARAMETER (IM14 =47, IM15 =-8188, IM16 =8189, IM17=1)
CC     -----------------------------------------------------------------
c++ Code for SYS = CRAY1_SD is INACTIVE
CC     MACHINE CONSTANTS FOR THE CRAY 1, XMP, 2, AND 3, WHEN DOUBLE
CC     PRECISION IS TO USE SINGLE PRECISION ARITHMETIC.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =102, IM4 =6)
C      PARAMETER (IM5 =64, IM6 =8, IM7 =2, IM8 =46)
C      PARAMETER (IM10 =2, IM11 =47, IM12 =-8189, IM13 =8190)
C      PARAMETER (IM14 =47, IM15 =-8189, IM16 =8190, IM17=1)
CC     -----------------------------------------------------------------
c++ Code for SYS = CRAY1_64 is INACTIVE
CC     MACHINE CONSTANTS FOR THE CRAY 1, XMP, 2, AND 3.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =102, IM4 =6)
C      PARAMETER (IM5 =64, IM6 =8, IM7 =2, IM8 =63)
C      PARAMETER (IM10 =2, IM11 =47, IM12 =-8189, IM13 =8190)
C      PARAMETER (IM14 =94, IM15 =-8099, IM16 =8190, IM17=2)
CC     -----------------------------------------------------------------
c++ Code for SYS = CRAY1_SD_64 is INACTIVE
CC     MACHINE CONSTANTS FOR THE CRAY 1, XMP, 2, AND 3, WHEN DOUBLE
CC     PRECISION IS TO USE SINGLE PRECISION ARITHMETIC.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =102, IM4 =6)
C      PARAMETER (IM5 =64, IM6 =8, IM7 =2, IM8 =63)
C      PARAMETER (IM10 =2, IM11 =47, IM12 =-8189, IM13 =8190)
C      PARAMETER (IM14 =47, IM15 =-8189, IM16 =8190, IM17=1)
CC     -----------------------------------------------------------------
c++ Code for SYS = DG_S2000 is INACTIVE
CC     MACHINE CONSTANTS FOR THE DATA GENERAL ECLIPSE S/200
CC
C      PARAMETER (IM1 =11, IM2 =12, IM3 =8, IM4 =10)
C      PARAMETER (IM5 =16, IM6 =2, IM7 =2, IM8 =15)
C      PARAMETER (IM10 =16, IM11 =6, IM12 =-64, IM13 =63)
C      PARAMETER (IM14 =14, IM15 =-64, IM16 =63, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = HARRIS220 is INACTIVE
CC     MACHINE CONSTANTS FOR THE HARRIS 220, SLASH 6, SLASH 7.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =0, IM4 =6)
C      PARAMETER (IM5 =24, IM6 =3, IM7 =2, IM8 =23)
C      PARAMETER (IM10 =2, IM11 =23, IM12 =-127, IM13 =127)
C      PARAMETER (IM14 =38, IM15 =-127, IM16 =127, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = HON600_6000 is INACTIVE
CC     MACHINE CONSTANTS FOR THE HONEYWELL 600/6000 SERIES.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =43, IM4 =6)
C      PARAMETER (IM5 =36, IM6 =6, IM7 =2, IM8 =35)
C      PARAMETER (IM10 =2, IM11 =27, IM12 =-127, IM13 =127)
C      PARAMETER (IM14 =63, IM15 =-127, IM16 =127, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = HON_DPS_8_70 is INACTIVE
CC     MACHINE CONSTANTS FOR THE HONEYWELL DPS 8/70 SERIES.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =43, IM4 =6)
C      PARAMETER (IM5 =36, IM6 =4, IM7 =2, IM8 =35)
C      PARAMETER (IM10 =2, IM11 =27, IM12 =-127, IM13 =127)
C      PARAMETER (IM14 =63, IM15 =-127, IM16 =127, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = HP700Q is INACTIVE
cc     Machine constants for HP-700 using the +autodblpad option,
cc     which automatically increases DOUBLE PRECISION to REAL*16, and
cc     REAL to DOUBLE PRECISION.
c      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
c      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
c      PARAMETER (IM10 =2, IM11 =53, IM12 =-1021, IM13 =1024)
c      PARAMETER (IM14 = 113, IM15 = -16381, IM16 = 16384, IM17 = 0)
CC     -----------------------------------------------------------------
c++ Code for SYS = IBM360_370 is INACTIVE
CC     MACHINE CONSTANTS FOR THE IBM 360/370 SERIES,
CC     THE XEROX SIGMA 5/7/9 AND THE SEL SYSTEMS 85/86.
Cc
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
C      PARAMETER (IM10 =16, IM11 =6, IM12 =-64, IM13 =63)
C      PARAMETER (IM14 =14, IM15 =-64, IM16 =63, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = INTERDATA_8_32 is INACTIVE
CC     MACHINE CONSTANTS FOR THE INTERDATA 8/32
CC     WITH THE UNIX SYSTEM FORTRAN 77 COMPILER.
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =6, IM4 =6)
C      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
C      PARAMETER (IM10 =16, IM11 =6, IM12 =-64, IM13 =62)
C      PARAMETER (IM14 =14, IM15 =-64, IM16 =62, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = PDP10_KA is INACTIVE
CC     MACHINE CONSTANTS FOR THE PDP-10 (KA PROCESSOR).
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =36, IM6 =5, IM7 =2, IM8 =35)
C      PARAMETER (IM10 =2, IM11 =27, IM12 =-128, IM13 =127)
C      PARAMETER (IM14 =54, IM15 =-101, IM16 =127, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = PDP10_KB is INACTIVE
CC     MACHINE CONSTANTS FOR THE PDP-10 (KI PROCESSOR).
CC
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =36, IM6 =5, IM7 =2, IM8 =35)
C      PARAMETER (IM10 =2, IM11 =27, IM12 =-128, IM13 =127)
C      PARAMETER (IM14 =62, IM15 =-128, IM16 =127, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = PDP11 is INACTIVE
CC     MACHINE CONSTANTS FOR PDP-11 FORTRAN'S SUPPORTING
CC     16-BIT INTEGER ARITHMETIC.
Cc
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =16, IM6 =2, IM7 =2, IM8 =15)
C      PARAMETER (IM10 =2, IM11 =24, IM12 =-127, IM13 =127)
C      PARAMETER (IM14 =56, IM15 =-127, IM16 =127, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = PRIME50 is INACTIVE
CC     MACHINE CONSTANTS FOR THE PRIME 50 SERIES SYSTEMS
CC     WITH 32-BIT INTEGERS AND 64V MODE INSTRUCTIONS,
CC     SUPPLIED BY IGOR BRAY.
Cc
C      PARAMETER (IM1 =1, IM2 =1, IM3 =2, IM4 =1)
C      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
C      PARAMETER (IM10 =2, IM11 =23, IM12 =-127, IM13 =127)
C      PARAMETER (IM14 =47, IM15 =-32895, IM16 =32637, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = SEQ_BAL_8000 is INACTIVE
CC     MACHINE CONSTANTS FOR THE SEQUENT BALANCE 8000.
CC
C      PARAMETER (IM1 =0, IM2 =0, IM3 =7, IM4 =0)
C      PARAMETER (IM5 =32, IM6 =1, IM7 =2, IM8 =31)
C      PARAMETER (IM10 =2, IM11 =24, IM12 =-125, IM13 =128)
C      PARAMETER (IM14 =53, IM15 =-1021, IM16 =1024, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = UNIVAC is INACTIVE
CC     MACHINE CONSTANTS FOR THE UNIVAC 1100 SERIES.
CC
CC     NOTE THAT THE PUNCH UNIT, I1MACH(3), HAS BEEN SET TO 1
CC     WHICH IS APPROPRIATE FOR THE UNIVAC-FTN SYSTEM.
CC     IF YOU HAVE THE UNIVAC-FOR SYSTEM, SET IT TO 7.
CC     IM6 = 4 for FTN (4 chars per word), 6 for FOR (6 chars per word).
Cc
C      PARAMETER (IM1 =5, IM2 =6, IM3 =1, IM4 =6)
C      PARAMETER (IM5 =36, IM6 =4, IM7 =2, IM8 =35)
C      PARAMETER (IM10 =2, IM11 =27, IM12 =-128, IM13 =127)
C      PARAMETER (IM14 =60, IM15 =-1024, IM16 =1023, IM17=0)
CC     -----------------------------------------------------------------
c++ Code for SYS = VAX is INACTIVE
Cc     MACHINE CONSTANTS for the VAX/VMS F and D formats
Cc     and for PDP-11 FORTRAN SUPPORTING 32-BIT INTEGER ARITHMETIC.
Cc
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
C      PARAMETER (IM10 =2, IM11 =24, IM12 =-127, IM13 =127)
C      PARAMETER (IM14 =56, IM15 =-127, IM16 =127, IM17=0)
c++ end
C     -----------------------------------------------------------------
c++ Code for SYS = VAX_G is INACTIVE
Cc     MACHINE CONSTANTS for the VAX/VMS F and G formats
Cc     and for PDP-11 FORTRAN SUPPORTING 32-BIT INTEGER ARITHMETIC.
Cc
C      PARAMETER (IM1 =5, IM2 =6, IM3 =7, IM4 =6)
C      PARAMETER (IM5 =32, IM6 =4, IM7 =2, IM8 =31)
C      PARAMETER (IM10 =2, IM11 =24, IM12 =-127, IM13 =127)
C      PARAMETER (IM14 =53, IM15 =-1023, IM16 =1023, IM17=0)
c++ end
C     -----------------------------------------------------------------
C
C
C Real parameters
C
C  RM1 = R1MACH(1) = B**(EMIN-1), The smallest positive number, i.e.,
c                    the underflow limit.
C  RM2 = R1MACH(2) = B**EMAX*(1 - B**(-T)), The largest number, i.e.,
c                    the overflow limit.
C  RM3 = R1MACH(3) = B**(-T), The smallest relative spacing, i.e., the
c                    difference between 1.0 and the next smaller number.
C  RM4 = R1MACH(4) = B**(1-T), The largest relative spacing, i.e., the
c                     difference between 1.0 and the next larger number.
C  RM5 = R1MACH(5) = LOG10(B).  When B = 2 this value is
c              Log10(2) = 0.30102_99956_63981_19521_37388_94724
C
C Parameter RMA and RMB are selected so that for values of the base =
C 2, 8, 16, 10, RMA has the values 1, 3, 4, 0, and RMB has the values 0,
C 0, 0, 1.  These values are used in computing RM5.
C $$$$ Note that if other bases are to be supported, the calculation of
C $$$$ RMA and RMB will have to be generalized.
C
c++   Code for HOW = COMPILER is ACTIVE
      PARAMETER (IM9 = 2 * (2**(IM8-1) - 1) + 1)
      PARAMETER (RMA = ((IM10 - 10) * (-3 + ((IM10 - 2) * (-77 +
     1    12 * (IM10 - 8))) / 14)) / 24)
      PARAMETER (RMB = ((IM10 - 2) * (IM10 - 8) * (16 - IM10)) / 96)
      PARAMETER (RBASE = IM10)
C
C     Weird subterfuges below are NECESSARY to compute DM1 and DM2 on
C     some systems.  DON'T SIMPLIFY THEM.  We compute RM1 and RM2 using
C     these subterfuges so it will be clear we're computing the REAL
C     and DOUBLE PRECISION characteristics in the same way.
      PARAMETER (RM1 = (RBASE**(IM12/2)) * (RBASE**(IM12-IM12/2-1)))
      PARAMETER (RM2 = RBASE**(IM13-IM11) * ((RBASE**IM11 - RBASE)
     1               + (RBASE - 1.0E0)))
      PARAMETER (RM3 = RBASE**(-IM11))
      PARAMETER (RM4 = RBASE**(1-IM11))
c     PARAMETER (RM5 = RMA*0.30102 99956 63981 19521 37388 94724E0+RMB)
      PARAMETER (RM5 = RMA*0.301029995663981195213738894724E0+RMB)
C
C Double precision parameters -- (Defined like the real ones.)
C
      PARAMETER (DMA = ((IM10 - 10) * (-3 + ((IM10 - 2) * (-77 +
     1    12 * (IM10 - 8))) / 14)) / 24)
      PARAMETER (DMB = ((IM10 - 2) * (IM10 - 8) * (16 - IM10)) / 96)
      PARAMETER (DBASE = IM10)
C
C     Weird subterfuges below are NECESSARY to compute DM1 and DM2 on
C     some systems.  DON'T SIMPLIFY THEM.
C      If you are using the Absoft Pro Fortran for Windows v8.0
C      Uncomment the following line and comment the one after.
C      PARAMETER (DM1=.222507385850721D-307)
      PARAMETER (DM1 = (DBASE**(IM15/2)) * (DBASE**(IM15-IM15/2-1)))
      PARAMETER (DM2 = DBASE**(IM16-IM14) * ((DBASE**IM14 - DBASE)
     1               + (DBASE - 1.0D0)))
      PARAMETER (DM3 = DBASE**(-IM14))
      PARAMETER (DM4 = DBASE**(1-IM14))
c     PARAMETER (DM5 = DMA *
c    1 0.30102 99956 63981 19521 37388 94724 49302 67681 89881 46211 D0
c    2 + DMB)
c
      PARAMETER (DM5 = DMA*
     1 0.30102999566398119521373889472449302676818988146211D0 + DMB)
C DM6 and TEST are used in checking that the correct constants have
C been selected.
      PARAMETER (DM6 = DBASE**(-IM14-IM17))
c++   END
      data TEST / 0.D0 /
C
c     DATA IMACH / IM1, IM2, IM3, IM4, IM5, IM6, IM7, IM8, IM9, IM10,
c    1   IM11, IM12, IM13, IM14, IM15, IM16 /
c     DATA RMACH / RM1, RM2, RM3, RM4, RM5 /
c     DATA DMACH / DM1, DM2, DM3, DM4, DM5 /
C     -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
C
      if (TEST .eq. 0.0D0) then
C         IM9 = 2 * (2**(IM8-1) - 1) + 1
         IMACH(1) = IM1
         IMACH(2) = IM2
         IMACH(3) = IM3
         IMACH(4) = IM4
         IMACH(5) = IM5
         IMACH(6) = IM6
         IMACH(7) = IM7
         IMACH(8) = IM8
         IMACH(10) = IM10
         IMACH(11) = IM11
         IMACH(12) = IM12
         IMACH(13) = IM13
         IMACH(14) = IM14
         IMACH(15) = IM15
         IMACH(16) = IM16
         IMACH(17) = IM17
c++   Code for HOW = RUN is INACTIVE
C         IEEE = 0
C100      continue
C      DBASE = IMACH(10)
CC
CC     Weird subterfuge below is NECESSARY to compute DM1 on
CC     some systems.  DON'T SIMPLIFY IT.
C      DM1=(DBASE**(IMACH(15)/2)) * (DBASE**(IMACH(15)-IMACH(15)/2-1))
CC DM6 and TEST are used in checking that the correct constants have
CC been selected.
C      DM6 = DBASE**(-IMACH(14)-IMACH(17))
c++   end
         CALL AMTEST (TEST, DM6)
         if (dm1 .eq. 0.0d0 .or. test .eq. 0.0d0) then
c++   Code for HOW = RUN is INACTIVE
C           if (IEEE .eq. 0) then
C              IEEE = 1
C              IMACH(1) = IE1
C              IMACH(2) = IE2
C              IMACH(3) = IE3
C              IMACH(4) = IE4
C              IMACH(5) = IE5
C              IMACH(6) = IE6
C              IMACH(7) = IE7
C              IMACH(8) = IE8
C              IMACH(10) = IE10
C              IMACH(11) = IE11
C              IMACH(12) = IE12
C              IMACH(13) = IE13
C              IMACH(14) = IE14
C              IMACH(15) = IE15
C              IMACH(16) = IE16
C              IMACH(17) = IE17
C              go to 100
C           end if
C           if (IEEE .eq. 1) then
C              IEEE = 2
C              IMACH(1) = ID1
C              IMACH(2) = ID2
C              IMACH(3) = ID3
C              IMACH(4) = ID4
C              IMACH(5) = ID5
C              IMACH(6) = ID6
C              IMACH(7) = ID7
C              IMACH(8) = ID8
C              IMACH(10) = ID10
C              IMACH(11) = ID11
C              IMACH(12) = ID12
C              IMACH(13) = ID13
C              IMACH(14) = ID14
C              IMACH(15) = ID15
C              IMACH(16) = ID16
C              IMACH(17) = ID17
C              go to 100
C           end if
c++   END
            print*,'AMACH has bad parameters for current environment.'
            stop
         end if
c++   Code for HOW = RUN is INACTIVE
C         IM9 = 2 * (2**(IMACH(8)-1) - 1) + 1
C         RMA = ((IMACH(10) - 10) * (-3 + ((IMACH(10) - 2) * (-77 +
C     1       12 * (IMACH(10) - 8))) / 14)) / 24
C         RMB = ((IMACH(10)-2) * (IMACH(10)-8) * (16-IMACH(10)))/96
C         RBASE = IMACH(10)
CC
CC        Weird subterfuges below are NECESSARY to compute DM1 and DM2
CC        on some systems.  DON'T SIMPLIFY THEM.  We compute RM1 and
CC        RM2 using these subterfuges so it will be clear we're
CC        computing the REAL and DOUBLE PRECISION characteristics in
Cc        the same way.
C         RM1=(RBASE**(IMACH(12)/2))*(RBASE**(IMACH(12)-IMACH(12)/2-1))
C         RM2 = RBASE**(IMACH(13)-IMACH(11))*((RBASE**IMACH(11) - RBASE)
C     1                  + (RBASE - 1.0E0))
C         RM3 = RBASE**(-IMACH(11))
C         RM4 = RBASE**(1-IMACH(11))
Cc        RM5 = RMA*0.30102 99956 63981 19521 37388 94724E0+RMB
C         RM5 = RMA*0.301029995663981195213738894724E0+RMB
CC
CC Double precision parameters -- (Defined like the real ones.)
CC
C         DMA = ((IMACH(10) - 10) * (-3 + ((IMACH(10) - 2) * (-77 +
C     1       12 * (IMACH(10) - 8))) / 14)) / 24
C         DMB = ((IMACH(10)-2) * (IMACH(10)-8) * (16-IMACH(10)))/96
CC
CC        Weird subterfuge below is NECESSARY to compute DM2 on
CC        some systems.  DON'T SIMPLIFY IT.
C         DM2 = DBASE**(IMACH(16)-IMACH(14))*((DBASE**IMACH(14) - DBASE)
C     1                  + (DBASE - 1.0D0))
C         DM3 = DBASE**(-IMACH(14))
C         DM4 = DBASE**(1-IMACH(14))
Cc        DM5 = DMA*0.30102 99956 63981 19521 37388 94724D0+DMB
C         DM5 = DMA*0.301029995663981195213738894724D0+DMB
c++   END
         IMACH(9) = IM9
         RMACH(1) = RM1
         RMACH(2) = RM2
         RMACH(3) = RM3
         RMACH(4) = RM4
         RMACH(5) = RM5
         DMACH(1) = DM1
         DMACH(2) = DM2
         DMACH(3) = DM3
         DMACH(4) = DM4
         DMACH(5) = DM5
      ENDIF
C
      if (MODE .eq. 0) then
         I1=IMACH(I)
      else if (MODE .eq. 1) then
         R1=RMACH(I)
c                                  Here we assume MODE = 2.
      else
         D1=DMACH(I)
      endif
      return
      end
c     ==================================================================
      integer function I1MACH(I)
      integer I, I1
      real R1
      double precision D1
      IF (I .LT. 1  .OR.  I .GT. 16) THEN
         PRINT*,'I1MACH.. Bad argument: I =',I
         STOP 'I1MACH error'
      END IF
      call AMACH (0, I, I1, R1, D1)
      I1MACH = I1
      return
      end
c     ==================================================================
c
      real function R1MACH(I)
      integer I, I1
      real R1
      double precision D1
      IF (I .lt. 1  .or.  I .gt. 5) THEN
         print*,'R1MACH.. Bad argument: I = ',I
         stop 'R1MACH error'
      END IF
      call AMACH (1, I, I1, R1, D1)
      R1MACH = R1
      RETURN
      end
c     ==================================================================
c
      double precision function D1MACH(I)
      integer I, I1
      real R1
      double precision D1
      IF (I .lt. 1  .or.  I .gt. 5) THEN
         print*,'D1MACH.. Bad argument: I = ',I
         stop 'D1MACH error'
      END IF
      call AMACH (2, I, I1, R1, D1)
      D1MACH = D1
      RETURN
      END
c     ==================================================================
c
      SUBROUTINE AMTEST (TEST, D6)
c Verifies that D6 is an appropriate value for DM6.
c Returns TEST = D6 + D6 - 1, .ne. 0 if D6 is an appropriate value for
c DM6, else returns TEST = 0.  The caller uses TEST = 0 as a signal to
c try again with IEEE settings (unless that's already been done).
      external AMSUB1
      DOUBLE PRECISION AMSUB1, D6, TEST
      TEST = AMSUB1(1.D0 + D6)
C
C The comparison with 1.875E0*D6 in the line below is to guard
C against the possibility that TEST is > 0 as a result of rounding
C up in the addition of D6 to 1.
C
      IF ((TEST .eq. 0.D0) .or. (TEST .gt. 1.875D0*D6)) THEN
         TEST = (D6 + D6) + 1.D0
         IF (AMSUB1(TEST) .ne. 0.D0) RETURN
      END IF
      test = 0.0d0
      END
c     ==================================================================
c
      DOUBLE PRECISION FUNCTION AMSUB1 (TEST1)
      DOUBLE PRECISION TEST1
C     Returns the value of TEST1 - 1.
      AMSUB1 = TEST1 - 1.0D0
      RETURN
      END


text:
      subroutine strdcl
      common /cdefio/ bp, buf (4096)
      integer bp
      integer buf
      common /cfname/ fcname (30)
      integer fcname
      common /cfor/ fordep, forstk (200)
      integer fordep
      integer forstk
      common /cgoto/ xfer
      integer xfer
      common /clabel/ label, retlab, memflg, col, logic0
      integer label
      integer retlab
      integer memflg
      integer col
      integer logic0
      common /cline/ dbgout, dbglev, level, linect (5), infile (5), fnam
     *p, fnames ( 150)
      integer dbgout
      integer dbglev
      integer level
      integer linect
      integer infile
      integer fnamp
      integer fnames
      common /cmacro/ cp, ep, evalst (500), deftbl
      integer cp
      integer ep
      integer evalst
      integer deftbl
      common /coutln/ outp, outbuf (74)
      integer outp
      integer outbuf
      common /csbuf/ sbp, sbuf(2048), smem(240)
      integer sbp
      integer sbuf
      integer smem
      common /cswtch/ swtop, swlast, swstak(1000), swvnum, swvlev, swvst
     *k(10), swinrg
      integer swtop
      integer swlast
      integer swstak
      integer swvnum
      integer swvlev
      integer swvstk
      integer swinrg
      common /ckword/ rkwtbl
      integer rkwtbl
      common /clname/ fkwtbl, namtbl, gentbl, errtbl, xpptbl
      integer fkwtbl
      integer namtbl
      integer gentbl
      integer errtbl
      integer xpptbl
      common /erchek/ ername, body, esp, errstk(30)
      integer ername
      integer body
      integer esp
      integer errstk
      integer mem( 60000)
      common/cdsmem/mem
      integer t, token (100), dchar (100)
      integer gnbtok
      integer i, j, k, n, len
      integer length, ctoi, lex
      integer char(11)
      integer dat(6)
      integer eoss(3)
      data char(1)/105/,char(2)/110/,char(3)/116/,char(4)/101/,char(5)/1
     *03/,char(6)/101/,char(7)/114/,char(8)/42/,char(9)/50/,char(10)/47/
     *,char(11)/-2/
      data dat(1)/100/,dat(2)/97/,dat(3)/116/,dat(4)/97/,dat(5)/32/,dat(
     *6)/-2/
      data eoss(1)/48/,eoss(2)/47/,eoss(3)/-2/
      t = gnbtok (token, 100)
      if (.not.(t .ne. -9))goto 23000
      call synerr (21Hmissing string token.)
23000 continue
      call squash (token)
      call outtab
      call pbstr (char)
23002 continue
      t = gnbtok (dchar, 100)
      if (.not.(t .eq. 47))goto 23005
      goto 23004
23005 continue
      call outstr (dchar)
23003 goto 23002
23004 continue
      call outch (32)
      call outstr (token)
      call addstr (token, sbuf, sbp, 2048)
      call addchr (-2, sbuf, sbp, 2048)
      if (.not.(gnbtok (token, 100) .ne. 40))goto 23007
      len = length (token) + 1
      if (.not.(token (1) .eq. 39 .or. token (1) .eq. 34))goto 23009
      len = len - 2
23009 continue
      goto 23008
23007 continue
      t = gnbtok (token, 100)
      i = 1
      len = ctoi (token, i)
      if (.not.(token (i) .ne. -2))goto 23011
      call synerr (20Hinvalid string size.)
23011 continue
      if (.not.(gnbtok (token, 100) .ne. 41))goto 23013
      call synerr (20Hmissing right paren.)
      goto 23014
23013 continue
      t = gnbtok (token, 100)
23014 continue
23008 continue
      call outch (40)
      call outnum (len)
      call outch (41)
      call outdon
      if (.not.(token (1) .eq. 39 .or. token (1) .eq. 34))goto 23015
      len = length (token)
      token (len) = -2
      call addstr (token (2), sbuf, sbp, 2048)
      goto 23016
23015 continue
      call addstr (token, sbuf, sbp, 2048)
23016 continue
      call addchr (-2, sbuf, sbp, 2048)
      t = lex (token)
      call pbstr (token)
      if (.not.(t .ne. -75))goto 23017
      i = 1
23019 if (.not.(i .lt. sbp))goto 23021
      call outtab
      call outstr (dat)
      k = 1
      j = i + length (sbuf (i)) + 1
23022 continue
      if (.not.(k .gt. 1))goto 23025
      call outch (44)
23025 continue
      call outstr (sbuf (i))
      call outch (40)
      call outnum (k)
      call outch (41)
      call outch (47)
      if (.not.(sbuf (j) .eq. -2))goto 23027
      goto 23024
23027 continue
      n = sbuf (j)
      call outnum (n)
      call outch (47)
      k = k + 1
23023 j = j + 1
      goto 23022
23024 continue
      call pbstr (eoss)
23029 continue
      t = gnbtok (token, 100)
      call outstr (token)
23030 if (.not.(t .eq. 47))goto 23029
23031 continue
      call outdon
23020 i = j + 1
      goto 23019
23021 continue
      sbp = 1
23017 continue
      return
      end
c     logic0  logical_column

max_stars_repo_name: ysBach/irafdocgen

max_issues_repo_name: kirxkirx/iraf

max_forks_repo_name: kirxkirx/iraf