text:
      subroutine rdgrid(ndat0s,ndbmax,ndbptg,ndbptl,qend,qerr,q500nd,
     $ udbval,xdbval)
c
c     This subroutine reads data from the standard log K temperature
c     grid. Other kinds of data are read from such grids.
c
c     This suboutine is called by:
c
c       EQPT/rdpar.f
c       EQPT/pcraq.f
c       EQPT/pcrsg.f
c
c-----------------------------------------------------------------------
c
c     Principal input:
c
c       ndbptg = the number of points on the grid
c       ndbptl = the maximum number of points on a single line
c
c     Principal output:
c
c       qend   = logical flag, = .true. if end-of-file was encountered
c       qerr   = logical flag, = .true. if a read format error was
c                  encountered
c       udbval = string array, the equivalent of xdbval
c       xdbval = array containing the data read from the grid
c
c-----------------------------------------------------------------------
c
      implicit none
c
c-----------------------------------------------------------------------
c
c     Calling sequence variable declarations.
c
      integer ndbmax
c
      integer ndat0s
c
      integer ndbptg,ndbptl
c
      logical qend,qerr,q500nd
c
      character(len=16) udbval(ndbmax)
c
      real(8) xdbval(ndbmax)
c
c-----------------------------------------------------------------------
c
c     Local variable declarations.
c
      integer i,ii,j,jj,k,n
c
      character(len=80) uline,ulbufa,ulbufb
      character(len=16) ux16
c
c-----------------------------------------------------------------------
c
      qend = .false.
      qerr = .false.
c
c     Note: the following coding generally assumes that a blank
c     space trails all but the last number on a line. If a blank
c     space is not found, it is assumed that a number ends after
c     four decimal places.
c
      n = 0
      do i = 1,ndbptg/ndbptl
        read (ndat0s,'(a)',end=990,err=995) uline
        ulbufa = uline
        do k = 1,ndbptl
          call lejust(ulbufa)
          ii = index(ulbufa,' ')
          jj = index(ulbufa,'.')
          if (jj .eq. 0) then
            jj = 80
          else
            jj = jj + 5
          endif
          ii = min(ii,jj)
          if (ii .gt. 1) then
            ux16 = ulbufa(1:ii - 1)
            call locase(ux16)
            udbval(n + k) = ux16
            if (ux16(1:7) .eq. 'no_data') then
              xdbval(n + k) = 9999999.
            else
              read (ux16,'(f10.4)',err=995) xdbval(n + k)
            endif
          else
            udbval(n + k) = ' '
            xdbval(n + k) = 0.
          endif
          if (k .lt. ndbptl) then
            ulbufb = ulbufa(ii:80)
            ulbufa = ulbufb
          endif
        enddo
c       read (uline,udbfmt,err=995) (xdbval(n + k), k = 1,ndbptl)
c       Note: udbfmt = '( (5x,6f10.4) )' with the "6" repeat
c       count replaced by ndbptl.
        n = n + ndbptl
      enddo
c
      j = mod(ndbptg,ndbptl)
      if (j .gt. 0) then
        read (ndat0s,'(a)',end=990,err=995) uline
        ulbufa = uline
        do k = 1,j
          call lejust(ulbufa)
          ii = index(ulbufa,' ')
          jj = index(ulbufa,'.')
          if (jj .eq. 0) then
            jj = 80
          else
            jj = jj + 5
          endif
          ii = min(ii,jj)
          if (ii .gt. 1) then
            ux16 = ulbufa(1:ii - 1)
            call locase(ux16)
            udbval(n + k) = ux16
            if (ux16(1:7) .eq. 'no_data') then
              xdbval(n + k) = 9999999.
            else
              read (ux16,'(f10.4)',err=995) xdbval(n + k)
            endif
          else
            udbval(n + k) = ' '
            xdbval(n + k) = 0.
          endif
          if (k .lt. j) then
            ulbufb = ulbufa(ii:80)
            ulbufa = ulbufb
          endif
        enddo
c       read (uline,udbfmt,err=995) (xdbval(n + k), k = 1,j)
        n = n + j
      endif
c
      if (q500nd) then
c
c       Apply the following filter. Treat values of "500.0000" as
c       "no data". This accommodates older EQ3/6 data files which
c       use "500.0000" to mean no data.
c
        do i = 1,ndbptg
          if (xdbval(i).gt.499.9999 .and. xdbval(i).lt.500.0001) then
            udbval(i) = 'no_data'
            xdbval(i) = 9999999.
          endif
        enddo
      endif
c
      go to 999
c
c* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
c
  990 qend = .true.
      go to 999
c
c* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
c
  995 qerr = .true.
      go to 999
c
c* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
c
  999 continue
      end

ext: f


text:
      SUBROUTINE DT_FLAHAD(Id,If1,If2,If3)

C***********************************************************************
C sampling of FLAvor composition for HADrons/photons                   *
C              ID         BAMJET-id of hadron                          *
C              IF1,2,3    flavor content                               *
C                         (u,d,s: 1,2,3;  au,ad,as: -1,-1,-3)          *
C Note:  -  u,d numbering as in BAMJET                                 *
C        -  All ID are accepted. ID above 39 taken from PHOJET.        *
C This version dated 12.03.96 is written by S. Roesler                 *
C***********************************************************************

      IMPLICIT NONE
      DOUBLE PRECISION DT_RNDM , one
      INTEGER Id , idcpc , idpdg , IDT_IPDGHA , If1 , If2 , If3 ,
     &        IPHO_PDG2ID , ipoqua , iq , ix , jsel , k
      SAVE

C auxiliary common for reggeon exchange (DTUNUC 1.x)
      INCLUDE 'inc/dtquar'
      INCLUDE 'inc/popar2'

      DIMENSION jsel(3,6)
      DATA jsel/1 , 2 , 3 , 2 , 3 , 1 , 3 , 1 , 2 , 1 , 3 , 2 , 2 , 1 ,
     &     3 , 3 , 2 , 1/
      DIMENSION ipoqua(3)

      one = 1.0D0
      IF ( Id.EQ.7 ) THEN
C photon (charge dependent flavour sampling)
         k = INT(DT_RNDM(one)*6.D0+1.D0)
         IF ( k.LE.4 ) THEN
            If1 = 2
            If2 = -2
         ELSE IF ( k.EQ.5 ) THEN
            If1 = 1
            If2 = -1
         ELSE
            If1 = 3
            If2 = -3
         END IF
         IF ( DT_RNDM(one).LT.0.5D0 ) THEN
            k = If1
            If1 = If2
            If2 = k
         END IF
         If3 = 0
      ELSE
C hadron
         ix = INT(1.0D0+5.99999D0*DT_RNDM(one))
         IF ( Id.LE.39 ) THEN
            If1 = MQUark(jsel(1,ix),Id)
            If2 = MQUark(jsel(2,ix),Id)
            If3 = MQUark(jsel(3,ix),Id)

C Use PHOJET quark tables if BAMJET id exceeds 39 (limit of MQUARK)
         ELSE
            idpdg = IDT_IPDGHA(Id)
            ! WRITE(6,*) 'DT_FLAHAD: ID, IDPDG', ID, IDPDG
            idcpc = IPHO_PDG2ID(idpdg)
            DO iq = 1 , 3
               ipoqua(iq) = IQ_list(iq,ABS(idcpc))
               ipoqua(iq) = SIGN(ipoqua(iq),idcpc)

C Swap u and d quarks to convert from PHOJET to DPMJET
               IF ( ABS(ipoqua(iq)).EQ.1 ) THEN
                  ipoqua(iq) = SIGN(2,ipoqua(iq))
               ELSE IF ( ABS(ipoqua(iq)).EQ.2 ) THEN
                  ipoqua(iq) = SIGN(1,ipoqua(iq))
               END IF
            END DO
            ! WRITE(6,*) 'DT_FLAHAD: ID, IPOQUA', ID, IPOQUA
            If1 = ipoqua(jsel(1,ix))
            If2 = ipoqua(jsel(2,ix))
            If3 = ipoqua(jsel(3,ix))
         END IF
         ! WRITE(6,*) 'DT_FLAHAD: ID, IF1, IF2, IF3', ID,IF1,IF2,IF3
         IF ( (If1.EQ.0) .AND. (If3.NE.0) ) THEN
            If1 = If3
            If3 = 0
         ELSE IF ( (If2.EQ.0) .AND. (If3.NE.0) ) THEN
            If2 = If3
            If3 = 0
         END IF
      END IF

      END SUBROUTINE

max_stars_repo_name: pzhristov/DPMJET

max_issues_repo_name: pzhristov/DPMJET

max_forks_repo_name: pzhristov/DPMJET


text:
      LOGICAL FUNCTION GREATR(A,B,NWDS)
      INTEGER NWDS
      DOUBLE PRECISION A(NWDS), B(NWDS)
      GREATR= A(1) .GT. B(1)
      RETURN
      END

max_stars_repo_name: jnnccc/ring-gravity

max_issues_repo_name: jnnccc/ring-gravity

max_forks_repo_name: jnnccc/ring-gravity


text:
      SUBROUTINE CLAGS2( UPPER, A1, A2, A3, B1, B2, B3, CSU, SNU, CSV,
     $                   SNV, CSQ, SNQ )
*
*  -- LAPACK auxiliary routine (version 3.0) --
*     Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
*     Courant Institute, Argonne National Lab, and Rice University
*     September 30, 1994
*
*     .. Scalar Arguments ..
      LOGICAL            UPPER
      REAL               A1, A3, B1, B3, CSQ, CSU, CSV
      COMPLEX            A2, B2, SNQ, SNU, SNV
*     ..
*
*  Purpose
*  =======
*
*  CLAGS2 computes 2-by-2 unitary matrices U, V and Q, such
*  that if ( UPPER ) then
*
*            U'*A*Q = U'*( A1 A2 )*Q = ( x  0  )
*                        ( 0  A3 )     ( x  x  )
*  and
*            V'*B*Q = V'*( B1 B2 )*Q = ( x  0  )
*                        ( 0  B3 )     ( x  x  )
*
*  or if ( .NOT.UPPER ) then
*
*            U'*A*Q = U'*( A1 0  )*Q = ( x  x  )
*                        ( A2 A3 )     ( 0  x  )
*  and
*            V'*B*Q = V'*( B1 0  )*Q = ( x  x  )
*                        ( B2 B3 )     ( 0  x  )
*  where
*
*    U = (     CSU      SNU ), V = (     CSV     SNV ),
*        ( -CONJG(SNU)  CSU )      ( -CONJG(SNV) CSV )
*
*    Q = (     CSQ      SNQ )
*        ( -CONJG(SNQ)  CSQ )
*
*  Z' denotes the conjugate transpose of Z.
*
*  The rows of the transformed A and B are parallel. Moreover, if the
*  input 2-by-2 matrix A is not zero, then the transformed (1,1) entry
*  of A is not zero. If the input matrices A and B are both not zero,
*  then the transformed (2,2) element of B is not zero, except when the
*  first rows of input A and B are parallel and the second rows are
*  zero.
*
*  Arguments
*  =========
*
*  UPPER   (input) LOGICAL
*          = .TRUE.: the input matrices A and B are upper triangular.
*          = .FALSE.: the input matrices A and B are lower triangular.
*
*  A1      (input) REAL
*  A2      (input) COMPLEX
*  A3      (input) REAL
*          On entry, A1, A2 and A3 are elements of the input 2-by-2
*          upper (lower) triangular matrix A.
*
*  B1      (input) REAL
*  B2      (input) COMPLEX
*  B3      (input) REAL
*          On entry, B1, B2 and B3 are elements of the input 2-by-2
*          upper (lower) triangular matrix B.
*
*  CSU     (output) REAL
*  SNU     (output) COMPLEX
*          The desired unitary matrix U.
*
*  CSV     (output) REAL
*  SNV     (output) COMPLEX
*          The desired unitary matrix V.
*
*  CSQ     (output) REAL
*  SNQ     (output) COMPLEX
*          The desired unitary matrix Q.
*
*  =====================================================================
*
*     .. Parameters ..
      REAL               ZERO, ONE
      PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0 )
*     ..
*     .. Local Scalars ..
      REAL               A, AUA11, AUA12, AUA21, AUA22, AVB11, AVB12,
     $                   AVB21, AVB22, CSL, CSR, D, FB, FC, S1, S2, SNL,
     $                   SNR, UA11R, UA22R, VB11R, VB22R
      COMPLEX            B, C, D1, R, T, UA11, UA12, UA21, UA22, VB11,
     $                   VB12, VB21, VB22
*     ..
*     .. External Subroutines ..
      EXTERNAL           CLARTG, SLASV2
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          ABS, AIMAG, CMPLX, CONJG, REAL
*     ..
*     .. Statement Functions ..
      REAL               ABS1
*     ..
*     .. Statement Function definitions ..
      ABS1( T ) = ABS( REAL( T ) ) + ABS( AIMAG( T ) )
*     ..
*     .. Executable Statements ..
*
      IF( UPPER ) THEN
*
*        Input matrices A and B are upper triangular matrices
*
*        Form matrix C = A*adj(B) = ( a b )
*                                   ( 0 d )
*
         A = A1*B3
         D = A3*B1
         B = A2*B1 - A1*B2
         FB = ABS( B )
*
*        Transform complex 2-by-2 matrix C to real matrix by unitary
*        diagonal matrix diag(1,D1).
*
         D1 = ONE
         IF( FB.NE.ZERO )
     $      D1 = B / FB
*
*        The SVD of real 2 by 2 triangular C
*
*         ( CSL -SNL )*( A B )*(  CSR  SNR ) = ( R 0 )
*         ( SNL  CSL ) ( 0 D ) ( -SNR  CSR )   ( 0 T )
*
         CALL SLASV2( A, FB, D, S1, S2, SNR, CSR, SNL, CSL )
*
         IF( ABS( CSL ).GE.ABS( SNL ) .OR. ABS( CSR ).GE.ABS( SNR ) )
     $        THEN
*
*           Compute the (1,1) and (1,2) elements of U'*A and V'*B,
*           and (1,2) element of |U|'*|A| and |V|'*|B|.
*
            UA11R = CSL*A1
            UA12 = CSL*A2 + D1*SNL*A3
*
            VB11R = CSR*B1
            VB12 = CSR*B2 + D1*SNR*B3
*
            AUA12 = ABS( CSL )*ABS1( A2 ) + ABS( SNL )*ABS( A3 )
            AVB12 = ABS( CSR )*ABS1( B2 ) + ABS( SNR )*ABS( B3 )
*
*           zero (1,2) elements of U'*A and V'*B
*
            IF( ( ABS( UA11R )+ABS1( UA12 ) ).EQ.ZERO ) THEN
               CALL CLARTG( -CMPLX( VB11R ), CONJG( VB12 ), CSQ, SNQ,
     $                      R )
            ELSE IF( ( ABS( VB11R )+ABS1( VB12 ) ).EQ.ZERO ) THEN
               CALL CLARTG( -CMPLX( UA11R ), CONJG( UA12 ), CSQ, SNQ,
     $                      R )
            ELSE IF( AUA12 / ( ABS( UA11R )+ABS1( UA12 ) ).LE.AVB12 /
     $               ( ABS( VB11R )+ABS1( VB12 ) ) ) THEN
               CALL CLARTG( -CMPLX( UA11R ), CONJG( UA12 ), CSQ, SNQ,
     $                      R )
            ELSE
               CALL CLARTG( -CMPLX( VB11R ), CONJG( VB12 ), CSQ, SNQ,
     $                      R )
            END IF
*
            CSU = CSL
            SNU = -D1*SNL
            CSV = CSR
            SNV = -D1*SNR
*
         ELSE
*
*           Compute the (2,1) and (2,2) elements of U'*A and V'*B,
*           and (2,2) element of |U|'*|A| and |V|'*|B|.
*
            UA21 = -CONJG( D1 )*SNL*A1
            UA22 = -CONJG( D1 )*SNL*A2 + CSL*A3
*
            VB21 = -CONJG( D1 )*SNR*B1
            VB22 = -CONJG( D1 )*SNR*B2 + CSR*B3
*
            AUA22 = ABS( SNL )*ABS1( A2 ) + ABS( CSL )*ABS( A3 )
            AVB22 = ABS( SNR )*ABS1( B2 ) + ABS( CSR )*ABS( B3 )
*
*           zero (2,2) elements of U'*A and V'*B, and then swap.
*
            IF( ( ABS1( UA21 )+ABS1( UA22 ) ).EQ.ZERO ) THEN
               CALL CLARTG( -CONJG( VB21 ), CONJG( VB22 ), CSQ, SNQ, R )
            ELSE IF( ( ABS1( VB21 )+ABS( VB22 ) ).EQ.ZERO ) THEN
               CALL CLARTG( -CONJG( UA21 ), CONJG( UA22 ), CSQ, SNQ, R )
            ELSE IF( AUA22 / ( ABS1( UA21 )+ABS1( UA22 ) ).LE.AVB22 /
     $               ( ABS1( VB21 )+ABS1( VB22 ) ) ) THEN
               CALL CLARTG( -CONJG( UA21 ), CONJG( UA22 ), CSQ, SNQ, R )
            ELSE
               CALL CLARTG( -CONJG( VB21 ), CONJG( VB22 ), CSQ, SNQ, R )
            END IF
*
            CSU = SNL
            SNU = D1*CSL
            CSV = SNR
            SNV = D1*CSR
*
         END IF
*
      ELSE
*
*        Input matrices A and B are lower triangular matrices
*
*        Form matrix C = A*adj(B) = ( a 0 )
*                                   ( c d )
*
         A = A1*B3
         D = A3*B1
         C = A2*B3 - A3*B2
         FC = ABS( C )
*
*        Transform complex 2-by-2 matrix C to real matrix by unitary
*        diagonal matrix diag(d1,1).
*
         D1 = ONE
         IF( FC.NE.ZERO )
     $      D1 = C / FC
*
*        The SVD of real 2 by 2 triangular C
*
*         ( CSL -SNL )*( A 0 )*(  CSR  SNR ) = ( R 0 )
*         ( SNL  CSL ) ( C D ) ( -SNR  CSR )   ( 0 T )
*
         CALL SLASV2( A, FC, D, S1, S2, SNR, CSR, SNL, CSL )
*
         IF( ABS( CSR ).GE.ABS( SNR ) .OR. ABS( CSL ).GE.ABS( SNL ) )
     $        THEN
*
*           Compute the (2,1) and (2,2) elements of U'*A and V'*B,
*           and (2,1) element of |U|'*|A| and |V|'*|B|.
*
            UA21 = -D1*SNR*A1 + CSR*A2
            UA22R = CSR*A3
*
            VB21 = -D1*SNL*B1 + CSL*B2
            VB22R = CSL*B3
*
            AUA21 = ABS( SNR )*ABS( A1 ) + ABS( CSR )*ABS1( A2 )
            AVB21 = ABS( SNL )*ABS( B1 ) + ABS( CSL )*ABS1( B2 )
*
*           zero (2,1) elements of U'*A and V'*B.
*
            IF( ( ABS1( UA21 )+ABS( UA22R ) ).EQ.ZERO ) THEN
               CALL CLARTG( CMPLX( VB22R ), VB21, CSQ, SNQ, R )
            ELSE IF( ( ABS1( VB21 )+ABS( VB22R ) ).EQ.ZERO ) THEN
               CALL CLARTG( CMPLX( UA22R ), UA21, CSQ, SNQ, R )
            ELSE IF( AUA21 / ( ABS1( UA21 )+ABS( UA22R ) ).LE.AVB21 /
     $               ( ABS1( VB21 )+ABS( VB22R ) ) ) THEN
               CALL CLARTG( CMPLX( UA22R ), UA21, CSQ, SNQ, R )
            ELSE
               CALL CLARTG( CMPLX( VB22R ), VB21, CSQ, SNQ, R )
            END IF
*
            CSU = CSR
            SNU = -CONJG( D1 )*SNR
            CSV = CSL
            SNV = -CONJG( D1 )*SNL
*
         ELSE
*
*           Compute the (1,1) and (1,2) elements of U'*A and V'*B,
*           and (1,1) element of |U|'*|A| and |V|'*|B|.
*
            UA11 = CSR*A1 + CONJG( D1 )*SNR*A2
            UA12 = CONJG( D1 )*SNR*A3
*
            VB11 = CSL*B1 + CONJG( D1 )*SNL*B2
            VB12 = CONJG( D1 )*SNL*B3
*
            AUA11 = ABS( CSR )*ABS( A1 ) + ABS( SNR )*ABS1( A2 )
            AVB11 = ABS( CSL )*ABS( B1 ) + ABS( SNL )*ABS1( B2 )
*
*           zero (1,1) elements of U'*A and V'*B, and then swap.
*
            IF( ( ABS1( UA11 )+ABS1( UA12 ) ).EQ.ZERO ) THEN
               CALL CLARTG( VB12, VB11, CSQ, SNQ, R )
            ELSE IF( ( ABS1( VB11 )+ABS1( VB12 ) ).EQ.ZERO ) THEN
               CALL CLARTG( UA12, UA11, CSQ, SNQ, R )
            ELSE IF( AUA11 / ( ABS1( UA11 )+ABS1( UA12 ) ).LE.AVB11 /
     $               ( ABS1( VB11 )+ABS1( VB12 ) ) ) THEN
               CALL CLARTG( UA12, UA11, CSQ, SNQ, R )
            ELSE
               CALL CLARTG( VB12, VB11, CSQ, SNQ, R )
            END IF
*
            CSU = SNR
            SNU = CONJG( D1 )*CSR
            CSV = SNL
            SNV = CONJG( D1 )*CSL
*
         END IF
*
      END IF
*
      RETURN
*
*     End of CLAGS2
*
      END

max_stars_repo_name: tenomoto/ncl

max_issues_repo_name: tenomoto/ncl

max_forks_repo_name: tenomoto/ncl


text:
c
c
c     =====================================================
      subroutine rpt2(ixy,maxm,meqn,mwaves,mbc,mx,
     &                  ql,qr,aux1,aux2,aux3,
     &                  ilr,asdq,bmasdq,bpasdq)
c     =====================================================
      use geoclaw_module

      implicit none
c
c     # Riemann solver in the transverse direction using an einfeldt
c     Jacobian.

c-----------------------last modified 1/10/05----------------------

      integer ixy,maxm,meqn,mwaves,mbc,mx,ilr

      double precision  ql(1-mbc:maxm+mbc, meqn)
      double precision  qr(1-mbc:maxm+mbc, meqn)
      double precision  asdq(1-mbc:maxm+mbc, meqn)
      double precision  bmasdq(1-mbc:maxm+mbc, meqn)
      double precision  bpasdq(1-mbc:maxm+mbc, meqn)
      double precision  aux1(1-mbc:maxm+mbc,*)
      double precision  aux2(1-mbc:maxm+mbc,*)
      double precision  aux3(1-mbc:maxm+mbc,*)

      double precision  s(3)
      double precision  r(3,3)
      double precision  beta(3)
      double precision  g,tol,abs_tol
      double precision  hl,hr,hul,hur,hvl,hvr,vl,vr,ul,ur,bl,br
      double precision  uhat,vhat,hhat,roe1,roe3,s1,s2,s3,s1l,s3r
      double precision  delf1,delf2,delf3,dxdcd,dxdcu
      double precision  dxdcm,dxdcp,topo1,topo3,eta

      integer i,m,mw,mu,mv

      g=grav
      tol=drytolerance
      abs_tol=drytolerance

      if (ixy.eq.1) then
	     mu = 2
	     mv = 3
      else
	     mu = 3
	     mv = 2
      endif

      do i=2-mbc,mx+mbc

         hl=qr(i-1,1)
         hr=ql(i,1)
         hul=qr(i-1,mu)
         hur=ql(i,mu)
         hvl=qr(i-1,mv)
         hvr=ql(i,mv)

c===========determine velocity from momentum===========================
       if (hl.lt.abs_tol) then
          hl=0.d0
          ul=0.d0
          vl=0.d0
       else
          ul=hul/hl
          vl=hvl/hl
       endif

       if (hr.lt.abs_tol) then
          hr=0.d0
          ur=0.d0
          vr=0.d0
       else
          ur=hur/hr
          vr=hvr/hr
       endif

       do mw=1,mwaves
          s(mw)=0.d0
          beta(mw)=0.d0
          do m=1,meqn
             r(m,mw)=0.d0
          enddo
       enddo
      dxdcp = 1.d0
      dxdcm = 1.d0

       if (hl.le.drytolerance.and.hr.le.drytolerance) go to 90

*      !check and see if cell that transverse waves are going in is high and dry
       if (ilr.eq.1) then
            eta = qr(i-1,1) + aux2(i-1,1)
            topo1 = aux1(i-1,1)
            topo3 = aux3(i-1,1)
c            s1 = vl-sqrt(g*hl)
c            s3 = vl+sqrt(g*hl)
c            s2 = 0.5d0*(s1+s3)
       else
            eta = ql(i,1) + aux2(i,1)
            topo1 = aux1(i,1)
            topo3 = aux3(i,1)
c            s1 = vr-sqrt(g*hr)
c            s3 = vr+sqrt(g*hr)
c            s2 = 0.5d0*(s1+s3)
       endif
       if (eta.lt.max(topo1,topo3)) go to 90

      if (icoordsys.eq.2) then
         if (ixy.eq.2) then
            dxdcp=(Rearth*pi/180.d0)
            dxdcm = dxdcp
         else
            if (ilr.eq.1) then
               dxdcp = Rearth*pi*cos(aux3(i-1,3))/180.d0
               dxdcm = Rearth*pi*cos(aux1(i-1,3))/180.d0
            else
               dxdcp = Rearth*pi*cos(aux3(i,3))/180.d0
               dxdcm = Rearth*pi*cos(aux1(i,3))/180.d0
            endif
         endif
      endif

c=====Determine some speeds necessary for the Jacobian=================
            vhat=(vr*dsqrt(hr))/(dsqrt(hr)+dsqrt(hl)) +
     &        (vl*dsqrt(hl))/(dsqrt(hr)+dsqrt(hl))

            uhat=(ur*dsqrt(hr))/(dsqrt(hr)+dsqrt(hl)) +
     &        (ul*dsqrt(hl))/(dsqrt(hr)+dsqrt(hl))
            hhat=(hr+hl)/2.d0

            roe1=vhat-dsqrt(g*hhat)
            roe3=vhat+dsqrt(g*hhat)

            s1l=vl-dsqrt(g*hl)
            s3r=vr+dsqrt(g*hr)

            s1=dmin1(roe1,s1l)
            s3=dmax1(roe3,s3r)

            s2=0.5d0*(s1+s3)

           s(1)=s1
           s(2)=s2
           s(3)=s3
c=======================Determine asdq decomposition (beta)============
         delf1=asdq(i,1)
         delf2=asdq(i,mu)
         delf3=asdq(i,mv)

         beta(1) = (s3*delf1/(s3-s1))-(delf3/(s3-s1))
         beta(2) = -s2*delf1 + delf2
         beta(3) = (delf3/(s3-s1))-(s1*delf1/(s3-s1))
c======================End =================================================

c=====================Set-up eigenvectors===================================
         r(1,1) = 1.d0
         r(2,1) = s2
         r(3,1) = s1

         r(1,2) = 0.d0
         r(2,2) = 1.d0
         r(3,2) = 0.d0

         r(1,3) = 1.d0
         r(2,3) = s2
         r(3,3) = s3
c============================================================================
90      continue
c============= compute fluctuations==========================================

            do  m=1,meqn
               bmasdq(i,m)=0.0d0
               bpasdq(i,m)=0.0d0
            enddo
            do  mw=1,3
               if (s(mw).lt.0.d0) then
                 bmasdq(i,1) =bmasdq(i,1) + dxdcm*s(mw)*beta(mw)*r(1,mw)
                 bmasdq(i,mu)=bmasdq(i,mu)+ dxdcm*s(mw)*beta(mw)*r(2,mw)
                 bmasdq(i,mv)=bmasdq(i,mv)+ dxdcm*s(mw)*beta(mw)*r(3,mw)
               elseif (s(mw).gt.0.d0) then
                 bpasdq(i,1) =bpasdq(i,1) + dxdcp*s(mw)*beta(mw)*r(1,mw)
                 bpasdq(i,mu)=bpasdq(i,mu)+ dxdcp*s(mw)*beta(mw)*r(2,mw)
                 bpasdq(i,mv)=bpasdq(i,mv)+ dxdcp*s(mw)*beta(mw)*r(3,mw)
               endif
            enddo
c========================================================================
         enddo
c

c
      return
      end

max_issues_repo_name: che-wenchao/D-Claw

max_forks_repo_name: che-wenchao/D-Claw


text:
! THIS VERSION: GALAHAD 2.1 - 22/03/2007 AT 09:00 GMT.

!-*-*-*-*-*-*-*-  G A L A H A D   C P U _ T I M E   M O D U L E  *-*-*-*-*-*-*-

!  Nick Gould, for GALAHAD productions
!  Copyright reserved
!  January 27th 1995

    MODULE GALAHAD_CPU_time

      IMPLICIT NONE

      PRIVATE
      PUBLIC :: CPU_TIME

    CONTAINS

!-*-*-*-*-*-  G A L A H A D   C P U _ T I M E   S U B R O U T I N E   -*-*-*-*-

      SUBROUTINE CPU_TIME( time )

!  Return the cpu time (to be replaced in 1995 by a Fortran intrinsic)

!NAG  USE F90_UNIX

!-----------------------------------------------
!   D u m m y   A r g u m e n t s
!-----------------------------------------------

      REAL ( KIND = KIND( 1.0E0 ) ), INTENT( OUT ) :: time

! For all machines using the NAG f90 compiler

!NAG  TYPE( TMS ) BUFFER
!NAG  INTEGER :: i
!NAG  i = TIMES( BUFFER )
!NAG  time = REAL( BUFFER%UTIME ) / REAL( CLOCK_TICKS_PER_SECOND( ) )

! For CRAY (UNICOS)

!CRAY CALL SECOND( time )

! For general Unix etime (Linux, DEC/Compaq OSF/Ultrix/Tru64, HP, SUN, SGI ...)

!GEN  REAL ( KIND = KIND( 1.0E0 ) ) :: ETIME, DUMMY( 2 )
!GEN  EXTERNAL ETIME
!GEN  time = ETIME( DUMMY )

! For DEC (VAX D arithmetic)

!VAXD INCLUDE '($JPIDEF)'
!VAXD INTEGER  i
!VAXD CALL LIB$GETJPI( JPI$_CPUTIM, , , i )
!VAXD time = I / 100.

! For DEC (VAX G arithmetic)

!VAXG INCLUDE '($JPIDEF)'
!VAXG INTEGER          i
!VAXG CALL LIB$GETJPI( JPI$_CPUTIM, , , i )
!VAXG time = I / 100.

! For IBM (AIX)

!RS6  INTEGER, EXTERNAL :: MCLOCK
!RS6  time = REAL( MCLOCK( ) ) / 100.0

! For IBM (CMS)

!IBM  INTEGER :: i
!IBM  REAL ( KIND = KIND = KIND( 1.0D0 ) ) :: dum8
!IBM  CALL CPUTIME( dum8, i )
!IBM  time = dum8 / 1000000.0

! For PC (Waterloo Fortran Compiler)

!WFC  INCLUDE 'FSUBLIB.FI'
!WFC  INTEGER * 2  hrs, mins, secs, hsecs
!WFC  CALL GETTIM( hrs, mins, secs, hsecs )
!WFC  time = 3600. * hrs + 60. * mins + secs + hsecs / 100.

! For PC (Salford Fortran Compiler)

!SALF CALL CLOCK@(DUM)
!SALF time = DUM

      RETURN

!  End of subroutine CPU_TIME

      END SUBROUTINE CPU_TIME

!  End of module GALAHAD_CPU_time

    END MODULE GALAHAD_CPU_time

max_forks_repo_name: dynamics-of-stellar-systems/dynamite


text:
Address(Braddock Court) is a residential Culdesacs culdesac in South Davis. There is greenbelt access at its southern end.

 Intersecting Streets

  Albany Avenue and across the intersection Benbow Court

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


text:
!
! ODPS_BIN2NC
!
! Program to convert netCDF format ODPS files to the CRTM Binary
! format.
!
!
! CREATION HISTORY:
!       Written by:     Yong Chen 27-Feb-2009
!                       Yong.Chen@noaa.gov
!

PROGRAM ODPS_BIN2NC

  ! ------------------
  ! Environment set up
  ! ------------------
  ! Module usage
  USE File_Utility   , ONLY: File_Exists
  USE Message_Handler, ONLY: SUCCESS, FAILURE, WARNING, INFORMATION, &
                             Program_Message, Display_Message
  USE ODPS_Define    , ONLY: ODPS_type, &
                             Destroy_ODPS, Equal_ODPS
  USE ODPS_Binary_IO , ONLY: Read_ODPS_Binary, Write_ODPS_Binary
  USE ODPS_netCDF_IO , ONLY: Read_ODPS_netCDF, Write_ODPS_netCDF
  ! Disable implicit typing
  IMPLICIT NONE

  ! ----------
  ! Parameters
  ! ----------
  CHARACTER(*), PARAMETER :: PROGRAM_NAME = 'ODPS_BIN2NC'
  CHARACTER(*), PARAMETER :: PROGRAM_RCS_ID = &

  ! ---------
  ! Variables
  ! ---------
  INTEGER :: Error_Status
  CHARACTER(256) :: NC_Filename
  CHARACTER(256) :: BIN_Filename
  TYPE(ODPS_type) :: ODPS
  TYPE(ODPS_type) :: ODPS_Test

  ! Output prgram header
  CALL Program_Message( PROGRAM_NAME, &
                        'Program to convert netCDF format ODPS files to '//&
                        'their CRTM Binary format.', &
                        '$Revision$' )

  ! Get the input and output filenames
  WRITE( *, FMT     = '( /5x, "Enter the INPUT Binary ODPS file: " )', &
            ADVANCE = 'NO' )
  READ( *, '( a )' ) BIN_Filename
  BIN_Filename = ADJUSTL( BIN_FileNAME )
  IF ( .NOT. File_Exists( TRIM( BIN_Filename ) ) ) THEN
    CALL Display_Message( PROGRAM_NAME, &
                          'File '//TRIM( BIN_Filename )//' not found.', &
                          FAILURE )
    STOP
  END IF

  WRITE( *, FMT     = '( /5x, "Enter the OUTPUT netCDF ODPS file: " )', &
            ADVANCE = 'NO' )
  READ( *, '( a )' ) NC_Filename
  NC_Filename = ADJUSTL( NC_Filename )

  ! Check that the BIN file isn't accidentally overwritten
  IF ( TRIM( BIN_Filename ) == TRIM( NC_Filename ) ) THEN
    CALL Display_Message( PROGRAM_NAME, &
                          'Output filename is the same as the input filename!', &
                          FAILURE )
    STOP
  END IF

  ! Read the input binary file
  WRITE( *, '( /5x, "Reading Binary ODPS data ..." )' )
  Error_Status = Read_ODPS_Binary( BIN_Filename, ODPS )
  IF ( Error_Status /= SUCCESS ) THEN
    CALL Display_Message( PROGRAM_NAME, &
                          'Error reading Binary ODPS file '//&
                          TRIM( BIN_Filename ), &
                          Error_Status )
    STOP
  END IF

  ! Write the netCDF file
  WRITE( *, '( /5x, "Writeing netCDF ODPS data ..." )' )
  Error_Status = Write_ODPS_netCDF( NC_Filename, ODPS )
  IF ( Error_Status /= SUCCESS ) THEN
    CALL Display_Message( PROGRAM_NAME, &
                          'Error writing netCDF ODPS file '//&
                          TRIM( NC_Filename ), &
                          Error_Status )
    STOP
  END IF

  ! Test read the netCDF data file
  WRITE( *, '( /5x, "Test reading the netCDF ODPS data file ..." )' )
  Error_Status = Read_ODPS_netCDF( NC_Filename, ODPS_Test )
  IF ( Error_Status /= SUCCESS ) THEN
    CALL Display_Message( PROGRAM_NAME, &
                          'Error reading netCDF ODPS file '//&
                          TRIM( NC_Filename ), &
                          Error_Status )
    STOP
  END IF

  ! Compare the two structures
  WRITE( *, '( /5x, "Comparing the Binary and netCDF ODPS structures ..." )' )
  Error_Status = Equal_ODPS( ODPS_Test, ODPS, Check_All=1  )
  IF ( Error_Status /= SUCCESS ) THEN
    CALL Display_Message( PROGRAM_NAME, &
                          'Differences found in Binary and netCDF'//&
                          'file ODPS structure comparison.', &
                          Error_Status )
  ELSE
    CALL Display_Message( PROGRAM_NAME, &
                          'Binary and netCDF file ODPS structures are equal.', &
                          INFORMATION )
  END IF

  ! Destroy the structures
  Error_Status = Destroy_ODPS( ODPS )
  IF ( Error_Status /= SUCCESS ) THEN
    CALL Display_Message( PROGRAM_NAME, &
                          'Error destroying ODPS structure.', &
                          WARNING )
  END IF

  Error_Status = Destroy_ODPS( ODPS_Test )
  IF ( Error_Status /= SUCCESS ) THEN
    CALL Display_Message( PROGRAM_NAME, &
                          'Error destroying ODPS_Test structure.', &
                          WARNING )
  END IF

END PROGRAM ODPS_BIN2NC

max_stars_repo_name: hsbadr/crtm

max_issues_repo_name: hsbadr/crtm

max_forks_repo_name: hsbadr/crtm


text:
        subroutine read_nc(filename,along,alat,s0,t0,p0,gamma0,a0,n0,
     &                                          along0,alat0,iocean0, ierr)
ccc
ccc
ccc
ccc     DESCRIPTION :   Read variables from the netcdf labelled data file
ccc
ccc     PRECISION:      Double precision
ccc
ccc     INPUT :
ccc            filename    full path to nc file
ccc             along           longitude of record
ccc                     alat            latitude of record
ccc
ccc     OUTPUT :        s0(nz,2,2)      array of cast salinities
ccc                     t0(nz,2,2)      array of cast in situ temperatures
ccc                     p0(nz)          array of cast pressures
ccc                     gamma0(nz,2,2)  array of cast gamma values
ccc                     a0(nz,2,2)      array of cast a values
ccc                     n0(2,2)         length of casts
ccc                     along0(2)       array of cast longitudes
ccc                     alat0(2)        array of cast latitudes
ccc                     iocean0(2,2)    array of cast oceans
ccc
ccc     UNITS :         salinity        psu (IPSS-78)
ccc                     temperature     degrees C (IPTS-68)
ccc                     pressure        db
ccc                     gamma           kg m-3
ccc
ccc
ccc     AUTHOR :        David Jackett
ccc
ccc     CREATED :       July 1993
ccc
ccc     REVISION :      1.3             15/11/94
ccc
ccc
ccc
        parameter(nx=90,ny=45,nz=33,ndx=4,ndy=4)

        implicit double precision (a-h,o-z)

        integer*4 id_gnc,id_lon,id_lat,id_p,id_n
        integer*4 id_ocean,id_s,id_t,id_gamma,id_a
        integer*4 start(3),count(3),n0(2,2),n0_t(2)

        integer*4 iocean(nx,ny),iocean0(2,2)

        real*4 along_s(nx),alat_s(ny)
        real*4 s0_s(nz,2,2),t0_s(nz,2,2),p0_s(nz),gamma0_s(nz,2,2),a0_s(nz,2,2)
        real*4 s0_t(nz,2),t0_t(nz,2),gamma0_t(nz,2),a0_t(nz,2)

        dimension along_d(nx),alat_d(ny),along0(2),alat0(2)
        dimension s0(nz,2,2),t0(nz,2,2),p0(nz),gamma0(nz,2,2),a0(nz,2,2)

        character*(*) filename

        save along_d, alat_d, i0, j0, iocean, p0_s

        data start/1,1,1/, count/nz,nx,ny/, nw/0/, i0/1/, j0/1/

        ilong(alng) = int(alng/ndx + 1)
        jlat(alt) = int((88+alt)/ndy + 1)

      ierr = 0

ccc
ccc             only read when you have to
ccc

        dx = along-along_d(i0)
        dy = alat-alat_d(j0)

        if(dx.lt.0.0.or.dx.ge.4.0.or.dy.lt.0.0.or.dy.ge.4.0.or.
     &                                      (i0.eq.1.and.j0.eq.1)) then

ccc
ccc             open the netcdf file and get the variable id's
ccc

C         filename = PWD//'/gamma.nc'
          ierr = nf_open(filename, 0, id_gnc)
            if (ierr.ne.0) then
              print *, ierr
                ierr = 13
                return
            end if
          ! Error handling above looks wrong...
          if(i0.eq.1.and.j0.eq.1) then

            ! These look unnecessary, duplicating data statement.
            start(2) = 1
            count(2) = nx
            start(3) = 1
            count(3) = ny

            ierr = nf_inq_varid(id_gnc, 'lon', id_lon)
            ierr = nf_get_var_double(id_gnc, id_lon, along_d)
            ierr = nf_inq_varid(id_gnc, 'lat', id_lat)
            ierr = nf_get_var_double(id_gnc, id_lat, alat_d)
            ierr = nf_inq_varid(id_gnc, 'pressure', id_p)
            ierr = nf_get_var_double(id_gnc, id_p, p0)

            ierr = nf_inq_varid(id_gnc, 'iocean', id_ocean)
            ierr = nf_get_var_int(id_gnc, id_ocean, iocean)

          end if

          ierr = nf_inq_varid(id_gnc, 'n', id_n)
          ierr = nf_inq_varid(id_gnc, 's', id_s)
          ierr = nf_inq_varid(id_gnc, 't', id_t)
          ierr = nf_inq_varid(id_gnc, 'gamma', id_gamma)
          ierr = nf_inq_varid(id_gnc, 'a', id_a)

ccc
ccc             read the appropriate records
ccc

          i0 = ilong(along)
          j0 = jlat(alat)

          along0(1) = along_d(i0)
          alat0(1) = alat_d(j0)
          alat0(2) = alat0(1)+ndy

          if(i0.eq.nx+1) i0 = 1

          if(i0.lt.nx) then

            along0(2) = along0(1)+ndx

            start(2) = i0
            count(2) = 2
            start(3) = j0
            count(3) = 2

            ierr = nf_get_vara_int(id_gnc, id_n, start(2), count(2), n0)
            ierr = nf_get_vara_double(id_gnc, id_s, start, count, s0)
            ierr = nf_get_vara_double(id_gnc, id_t, start, count, t0)
            ierr = nf_get_vara_double(id_gnc, id_gamma, start, count, gamma0)
            ierr = nf_get_vara_double(id_gnc, id_a, start, count, a0)

          elseif(i0.eq.nx) then

            start(2) = i0
            count(2) = 1
            start(3) = j0
            count(3) = 2

            ierr = nf_get_vara_int(id_gnc, id_n, start(2), count(2), n0_t)
            ierr = nf_get_vara_real(id_gnc, id_s, start, count, s0_t)
            ierr = nf_get_vara_real(id_gnc, id_t, start, count, t0_t)
            ierr = nf_get_vara_real(id_gnc, id_gamma, start, count, gamma0_t)
            ierr = nf_get_vara_real(id_gnc, id_a, start, count, a0_t)

            do j = 1,2
              n0(1,j) = n0_t(j)
              do k = 1,nz
                s0_s(k,1,j) = s0_t(k,j)
                t0_s(k,1,j) = t0_t(k,j)
                gamma0_s(k,1,j) = gamma0_t(k,j)
                a0_s(k,1,j) = a0_t(k,j)
              end do
            end do

            along0(2) = 0.0

            start(2) = 1
            count(2) = 1
            start(3) = j0
            count(3) = 2

            ierr = nf_get_vara_int(id_gnc, id_n, start(2), count(2), n0_t)
            ierr = nf_get_vara_real(id_gnc, id_s, start, count, s0_t)
            ierr = nf_get_vara_real(id_gnc, id_t, start, count, t0_t)
            ierr = nf_get_vara_real(id_gnc, id_gamma, start, count, gamma0_t)
            ierr = nf_get_vara_real(id_gnc, id_a, start, count, a0_t)

            do j = 1,2
              n0(2,j) = n0_t(j)
              do k = 1,nz
                s0_s(k,2,j) = s0_t(k,j)
                t0_s(k,2,j) = t0_t(k,j)
                gamma0_s(k,2,j) = gamma0_t(k,j)
                a0_s(k,2,j) = a0_t(k,j)
              end do
            end do

          end if

          ! call ncclos(id_gnc,ierr)
          ierr = nf_close(id_gnc)

ccc
ccc             get the ocean information
ccc

          do j = 1,2
          do i = 1,2
            iocean0(i,j) = iocean(ilong(along0(i)),jlat(alat0(j)))
          end do
          end do

        end if

        return
        end


text:
	SUBROUTINE SNHPLT  ( line, mkcolr, u, v,  npts,   utag, vtag,
     +			     tag,  ntag,   lvert, iret )
C************************************************************************
C* SNHPLT								*
C*									*
C* This subroutine plots the hodograph for SNHODO.			*
C*									*
C* SNHPLT  ( LINE,   MKCOLR, U, V, NPTS, UTAG, VTAG, TAG, NTAG, LVERT,	*
C*           IRET )							*
C*									*
C* Input parameters:							*
C*	LINE		CHAR*		Line color, type, width		*
C*	MKCOLR		INTEGER		Marker color			*
C*	U    (NPTS)	REAL		U for line			*
C*	V    (NPTS)	REAL		V for line			*
C*	NPTS		INTEGER		Number of line points		*
C*	UTAG (NTAG)	REAL		U for tag			*
C*	VTAG (NTAG)	REAL		V for tag			*
C*	TAG  (NTAG)	REAL		Vertical coordinate tags	*
C*	NTAG		INTEGER		Number of tags			*
C*	LVERT		INTEGER		Vertical coordinate type	*
C*									*
C* Output parameters:							*
C*	IRET		INTEGER		Return code			*
C*					  0 = normal return		*
C**									*
C* Log:									*
C* G. Huffman/USRA	 8/89	Adapted from SNPPLT			*
C* S. Schotz/GSC	 7/90	Remove motion parameter			*
C* S. Schotz/GSC	 8/90	Fix marker color, plot meters not KM	*
C* S. Schotz/GSC	10/90	Allow five characters for height text	*
C* K. Brill/NMC		07/91	Cleaned up				*
C************************************************************************
	INCLUDE		'GEMPRM.PRM'
C*
	CHARACTER*(*)	line
	REAL		u (*), v (*), utag (*), vtag (*), tag (*)
C*
	CHARACTER	ctag*8
	INTEGER		larr (3)
C----------------------------------------------------------------------
	iret  = 0
C
C*	Decode the line color, type, width.
C
	CALL ST_ILST  ( line, '/', 1, 3, larr, n, ier )
C
C*	Check that the color is non-zero.
C
	IF  ( larr (1) .eq. 0 )  RETURN
C
C*	Save current line characteristics, draw line, draw markers.
C
	CALL GQLINE  ( i1, i2, i3, i4, ier )
	CALL GSCOLR  ( larr (1), ier )
	CALL GSLINE  ( larr (2), 0, larr (3), 0, ier )
	CALL GLINE   ( 'M', npts, u, v, ier )
C
	IF  ( mkcolr .gt. 0 )  THEN
            CALL GSCOLR ( mkcolr, ier )
            CALL GMARK  ( 'M', ntag, utag, vtag, ier )
            CALL GSCOLR ( larr (1), ier )
        END IF
C
C*	Reset the line type.
C
	CALL GSLINE  ( i1, 0, i3, 0, ier )
C
C*	Convert each coordinate value to character (in whole numbers).
C
	IF  ( lvert .eq. 3 )  THEN
	    DO  i = 1, ntag
		CALL ST_RLCH  (  tag (i), 0, ctag, ier )
		CALL ST_LSTR  ( ctag, len, ier )
		IF  ( ctag ( len:len ) .eq. '.' )  THEN
		    ctag ( len:len ) = ' '
		    len = len - 1
		END IF
		CALL GTEXT  ( 'M', utag (i), vtag (i), ctag, 0.,
     +			      ((-2) * len),    0,        ier )
	    END DO
	  ELSE
	    DO  i = 1, ntag
		CALL ST_RLCH  ( tag (i),             0, ctag, ier )
		CALL ST_LSTR  ( ctag, len, ier )
		IF  ( ctag ( len:len ) .eq. '.' )  THEN
		    ctag ( len:len ) = ' '
		    len = len - 1
		END IF
		CALL GTEXT  ( 'M', utag (i), vtag (i), ctag, 0.,
     +			      ((-2) * len),    0,        ier )
	    END DO
	END IF
C*
	RETURN
	END

max_stars_repo_name: oxelson/gempak

max_issues_repo_name: oxelson/gempak

max_forks_repo_name: oxelson/gempak


        !COMPILER-GENERATED INTERFACE MODULE: Fri May 22 16:58:30 2020
        ! This source file is for reference only and may not completely
        ! represent the generated interface used by the compiler.
        MODULE COMPUTE_CURL_ONE_ELEMENT__genmod
          INTERFACE
            SUBROUTINE COMPUTE_CURL_ONE_ELEMENT(ISPEC)
              INTEGER(KIND=4) :: ISPEC
            END SUBROUTINE COMPUTE_CURL_ONE_ELEMENT
          END INTERFACE
        END MODULE COMPUTE_CURL_ONE_ELEMENT__genmod


text:
program ddot_timing
    use, intrinsic :: iso_fortran_env, only : error_unit, DP => REAL64
    implicit none
    integer :: n, runs, run, i, istat
    real(kind=DP), dimension(:), allocatable :: u, v
    real(kind=DP) :: scale_u, scale_v, r, r_ddot
    real(kind=DP), external :: ddot
    real :: start_time, end_time, ddot_time, loop_time

    call get_arguments(runs, n)
    allocate (u(n), v(n), stat=istat)
    if (istat /= 0) then
        write (unit=error_unit, fmt='(A)') 'error: can not allocate vectors'
        stop 2
    end if

    do run = 1, runs
        scale_u = 1.17*run
        scale_v = 0.17*run
        do i = 1, n
            u(i) = scale_u*sqrt(real(i))/n
            v(i) = scale_v*sqrt(real(i))/n
        end do
        call cpu_time(start_time)
        r_ddot = ddot(n, u, 1, v, 1)
        call cpu_time(end_time)
        ddot_time = end_time - start_time
        call cpu_time(start_time)
        r = 0.0_DP
        do i = 1, n
            r = r + u(i)*v(i)
        end do
        call cpu_time(end_time)
        loop_time = end_time - start_time
        print '(A, F10.6, A, F10.6)', 'ddot time: ', ddot_time, &
                                      ', loop time: ', loop_time
    end do

    deallocate (u, v)

contains

    subroutine get_arguments(runs, n)
        use, intrinsic :: iso_fortran_env, only : error_unit
        implicit none
        integer, intent(out) :: runs, n
        integer :: istat
        character(len=1024) :: buffer, msg

        runs = 1
        if (command_argument_count() >= 1) then
            call get_command_argument(1, buffer)
            read (buffer, fmt=*, iostat=istat, iomsg=msg) runs
            if (istat /= 0) then
                write (unit=error_unit, fmt='(2A)') 'error: ', trim(msg)
                stop 1
            end if
        end if
        n = 50000000
        if (command_argument_count() >= 2) then
            call get_command_argument(2, buffer)
            read (buffer, fmt=*) n
            if (istat /= 0) then
                write (unit=error_unit, fmt='(2A)') 'error: ', trim(msg)
                stop 1
            end if
        end if
    end subroutine get_arguments

end program ddot_timing

max_stars_repo_name: caguerra/Fortran-MOOC

max_issues_repo_name: gjbex/CMake-intro

max_forks_repo_name: gjbex/CMake-intro


text:
LOGICAL :: a, b, c
a = .TRUE.
a = .TRUE._4
b = .FALSE.
b = .FALSE._1
a = .NOT. b
a = b .AND. c
a = a .AND. b .AND. c
a = .NOT. a .AND. b .AND. .NOT. c
END PROGRAM

max_stars_repo_name: ouankou/rose

max_issues_repo_name: sujankh/rose-matlab

max_forks_repo_name: sujankh/rose-matlab


text:
The somewhat elusive Pomology Department studies the production of tree fruit and nut crops, and is the only such department in the country.  The department belongs to the College of Agricultural and Environmental Sciences.  It contributes to academic teaching in other related University Departments departments (like plant biology, crop science and management and biotechnology), but offers no majors major, minors minor or graduate groups graduate degree of its own.  It does, however, offer a http://caes.ucdavis.edu/intprog/eduopp/postgrcert/postgraduatecertificateinpomology postgraduate certificate in pomology through UC Davis Extension.

The department also has a reputation of giving away loads of tasty fruit and nut samples on Picnic Day.


text:
c###selmod.for
      SUBROUTINE SELMOD
C--------------------------------
C
C     THIS ROUTINE SELECTS THE SAMPLE AREA TO BE USED IN RAY PATH
C     CALCULATIONS.
C
      COMMON/FRQ/FREA(13),FREL(29),FREQ,JMODE,ITXRCP(2)
      COMMON /RON /CLAT(5), CLONG(5), GLAT(5), RD(5), FI(3,5), YI(3,5),
     1HI(3,5), HPRIM(30,5), HTRUE(30,5), FVERT(30,5),KM,KFX, AFAC(30,5),
     2HTR(50,3), FNSQ(50,3)
      IF(KFX.EQ.1)THEN
C.......ONLY ONE AVAILABLE
        JMODE =1
      ELSE IF(KFX.EQ.2)THEN
C.......ONLY TWO WITH THE SAME F (SET EQUAL IN IONSET)
C.......SO CHECK E LAYER
        IF(FI(1,1).GT.FI(1,2))THEN
C.........RECEIVER END CONTROLS
          JMODE =2
        ELSE
C.........TRANSMITTER END CONTROLS
          JMODE =1
        ENDIF
      ELSE
C.......THREE SAMPLES BUT CHECK ENDS ONLY
C.......CHECK  F2 LAYERS
        DELFI=FI(3,1)-FI(3,3)
        IF( ABS(DELFI).LT.0.01 )THEN
C.........F2 LAYER IS SAME SO CHECK E LAYER
          IF(FI(1,1).GT.FI(1,3))THEN
C...........RECEIVER END CONTROLS
            JMODE =3
          ELSE
C...........TRANSMITTER END CONTROLS
            JMODE =1
          ENDIF
        ELSE
C.........F2 LAYER CONTROLS
          IF(DELFI.GT.0.)THEN
C...........RECEIVER END CONTROLS
            JMODE =3
          ELSE
C...........TRANSMITTER END CONTROLS
            JMODE =1
          ENDIF
        ENDIF
      ENDIF
      RETURN
      END
C--------------------------------

ext: for

max_stars_repo_name: jawatson/voacapl

max_issues_repo_name: jawatson/voacapl

max_forks_repo_name: jawatson/voacapl


!*==slaruv.f90  processed by SPAG 7.51RB at 20:08 on  3 Mar 2022
!> \brief \b SLARUV returns a vector of n random real numbers from a uniform distribution.
!
!  =========== DOCUMENTATION ===========
!
! Online html documentation available at
!            http://www.netlib.org/lapack/explore-html/
!
!> \htmlonly
!> Download SLARUV + dependencies
!>
!> [TGZ]
!>
!> [ZIP]
!>
!> [TXT]
!> \endhtmlonly
!
!  Definition:
!  ===========
!
!       SUBROUTINE SLARUV( ISEED, N, X )
!
!       .. Scalar Arguments ..
!       INTEGER            N
!       ..
!       .. Array Arguments ..
!       INTEGER            ISEED( 4 )
!       REAL               X( N )
!       ..
!
!
!> \par Purpose:
!  =============
!>
!> \verbatim
!>
!> SLARUV returns a vector of n random real numbers from a uniform (0,1)
!> distribution (n <= 128).
!>
!> This is an auxiliary routine called by SLARNV and CLARNV.
!> \endverbatim
!
!  Arguments:
!  ==========
!
!> \param[in,out] ISEED
!> \verbatim
!>          ISEED is INTEGER array, dimension (4)
!>          On entry, the seed of the random number generator; the array
!>          elements must be between 0 and 4095, and ISEED(4) must be
!>          odd.
!>          On exit, the seed is updated.
!> \endverbatim
!>
!> \param[in] N
!> \verbatim
!>          N is INTEGER
!>          The number of random numbers to be generated. N <= 128.
!> \endverbatim
!>
!> \param[out] X
!> \verbatim
!>          X is REAL array, dimension (N)
!>          The generated random numbers.
!> \endverbatim
!
!  Authors:
!  ========
!
!> \author Univ. of Tennessee
!> \author Univ. of California Berkeley
!> \author Univ. of Colorado Denver
!> \author NAG Ltd.
!
!> \date December 2016
!
!> \ingroup OTHERauxiliary
!
!> \par Further Details:
!  =====================
!>
!> \verbatim
!>
!>  This routine uses a multiplicative congruential method with modulus
!>  2**48 and multiplier 33952834046453 (see G.S.Fishman,
!>  'Multiplicative congruential random number generators with modulus
!>  2**b: an exhaustive analysis for b = 32 and a partial analysis for
!>  b = 48', Math. Comp. 189, pp 331-344, 1990).
!>
!>  48-bit integers are stored in 4 integer array elements with 12 bits
!>  per element. Hence the routine is portable across machines with
!>  integers of 32 bits or more.
!> \endverbatim
!>
!  =====================================================================
      SUBROUTINE SLARUV(Iseed,N,X)
      IMPLICIT NONE
!*--SLARUV99
!
!  -- LAPACK auxiliary routine (version 3.7.0) --
!  -- LAPACK is a software package provided by Univ. of Tennessee,    --
!  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
!     December 2016
!
!     .. Scalar Arguments ..
      INTEGER N
!     ..
!     .. Array Arguments ..
      INTEGER Iseed(4)
      REAL X(N)
!     ..
!
!  =====================================================================
!
!     .. Parameters ..
      REAL ONE
      PARAMETER (ONE=1.0E0)
      INTEGER LV , IPW2
      REAL R
      PARAMETER (LV=128,IPW2=4096,R=ONE/IPW2)
!     ..
!     .. Local Scalars ..
      INTEGER i , i1 , i2 , i3 , i4 , it1 , it2 , it3 , it4 , j
!     ..
!     .. Local Arrays ..
      INTEGER mm(LV,4)
!     ..
!     .. Intrinsic Functions ..
      INTRINSIC MIN , MOD , REAL
!     ..
!     .. Data statements ..
      DATA (mm(1,j),j=1,4)/494 , 322 , 2508 , 2549/
      DATA (mm(2,j),j=1,4)/2637 , 789 , 3754 , 1145/
      DATA (mm(3,j),j=1,4)/255 , 1440 , 1766 , 2253/
      DATA (mm(4,j),j=1,4)/2008 , 752 , 3572 , 305/
      DATA (mm(5,j),j=1,4)/1253 , 2859 , 2893 , 3301/
      DATA (mm(6,j),j=1,4)/3344 , 123 , 307 , 1065/
      DATA (mm(7,j),j=1,4)/4084 , 1848 , 1297 , 3133/
      DATA (mm(8,j),j=1,4)/1739 , 643 , 3966 , 2913/
      DATA (mm(9,j),j=1,4)/3143 , 2405 , 758 , 3285/
      DATA (mm(10,j),j=1,4)/3468 , 2638 , 2598 , 1241/
      DATA (mm(11,j),j=1,4)/688 , 2344 , 3406 , 1197/
      DATA (mm(12,j),j=1,4)/1657 , 46 , 2922 , 3729/
      DATA (mm(13,j),j=1,4)/1238 , 3814 , 1038 , 2501/
      DATA (mm(14,j),j=1,4)/3166 , 913 , 2934 , 1673/
      DATA (mm(15,j),j=1,4)/1292 , 3649 , 2091 , 541/
      DATA (mm(16,j),j=1,4)/3422 , 339 , 2451 , 2753/
      DATA (mm(17,j),j=1,4)/1270 , 3808 , 1580 , 949/
      DATA (mm(18,j),j=1,4)/2016 , 822 , 1958 , 2361/
      DATA (mm(19,j),j=1,4)/154 , 2832 , 2055 , 1165/
      DATA (mm(20,j),j=1,4)/2862 , 3078 , 1507 , 4081/
      DATA (mm(21,j),j=1,4)/697 , 3633 , 1078 , 2725/
      DATA (mm(22,j),j=1,4)/1706 , 2970 , 3273 , 3305/
      DATA (mm(23,j),j=1,4)/491 , 637 , 17 , 3069/
      DATA (mm(24,j),j=1,4)/931 , 2249 , 854 , 3617/
      DATA (mm(25,j),j=1,4)/1444 , 2081 , 2916 , 3733/
      DATA (mm(26,j),j=1,4)/444 , 4019 , 3971 , 409/
      DATA (mm(27,j),j=1,4)/3577 , 1478 , 2889 , 2157/
      DATA (mm(28,j),j=1,4)/3944 , 242 , 3831 , 1361/
      DATA (mm(29,j),j=1,4)/2184 , 481 , 2621 , 3973/
      DATA (mm(30,j),j=1,4)/1661 , 2075 , 1541 , 1865/
      DATA (mm(31,j),j=1,4)/3482 , 4058 , 893 , 2525/
      DATA (mm(32,j),j=1,4)/657 , 622 , 736 , 1409/
      DATA (mm(33,j),j=1,4)/3023 , 3376 , 3992 , 3445/
      DATA (mm(34,j),j=1,4)/3618 , 812 , 787 , 3577/
      DATA (mm(35,j),j=1,4)/1267 , 234 , 2125 , 77/
      DATA (mm(36,j),j=1,4)/1828 , 641 , 2364 , 3761/
      DATA (mm(37,j),j=1,4)/164 , 4005 , 2460 , 2149/
      DATA (mm(38,j),j=1,4)/3798 , 1122 , 257 , 1449/
      DATA (mm(39,j),j=1,4)/3087 , 3135 , 1574 , 3005/
      DATA (mm(40,j),j=1,4)/2400 , 2640 , 3912 , 225/
      DATA (mm(41,j),j=1,4)/2870 , 2302 , 1216 , 85/
      DATA (mm(42,j),j=1,4)/3876 , 40 , 3248 , 3673/
      DATA (mm(43,j),j=1,4)/1905 , 1832 , 3401 , 3117/
      DATA (mm(44,j),j=1,4)/1593 , 2247 , 2124 , 3089/
      DATA (mm(45,j),j=1,4)/1797 , 2034 , 2762 , 1349/
      DATA (mm(46,j),j=1,4)/1234 , 2637 , 149 , 2057/
      DATA (mm(47,j),j=1,4)/3460 , 1287 , 2245 , 413/
      DATA (mm(48,j),j=1,4)/328 , 1691 , 166 , 65/
      DATA (mm(49,j),j=1,4)/2861 , 496 , 466 , 1845/
      DATA (mm(50,j),j=1,4)/1950 , 1597 , 4018 , 697/
      DATA (mm(51,j),j=1,4)/617 , 2394 , 1399 , 3085/
      DATA (mm(52,j),j=1,4)/2070 , 2584 , 190 , 3441/
      DATA (mm(53,j),j=1,4)/3331 , 1843 , 2879 , 1573/
      DATA (mm(54,j),j=1,4)/769 , 336 , 153 , 3689/
      DATA (mm(55,j),j=1,4)/1558 , 1472 , 2320 , 2941/
      DATA (mm(56,j),j=1,4)/2412 , 2407 , 18 , 929/
      DATA (mm(57,j),j=1,4)/2800 , 433 , 712 , 533/
      DATA (mm(58,j),j=1,4)/189 , 2096 , 2159 , 2841/
      DATA (mm(59,j),j=1,4)/287 , 1761 , 2318 , 4077/
      DATA (mm(60,j),j=1,4)/2045 , 2810 , 2091 , 721/
      DATA (mm(61,j),j=1,4)/1227 , 566 , 3443 , 2821/
      DATA (mm(62,j),j=1,4)/2838 , 442 , 1510 , 2249/
      DATA (mm(63,j),j=1,4)/209 , 41 , 449 , 2397/
      DATA (mm(64,j),j=1,4)/2770 , 1238 , 1956 , 2817/
      DATA (mm(65,j),j=1,4)/3654 , 1086 , 2201 , 245/
      DATA (mm(66,j),j=1,4)/3993 , 603 , 3137 , 1913/
      DATA (mm(67,j),j=1,4)/192 , 840 , 3399 , 1997/
      DATA (mm(68,j),j=1,4)/2253 , 3168 , 1321 , 3121/
      DATA (mm(69,j),j=1,4)/3491 , 1499 , 2271 , 997/
      DATA (mm(70,j),j=1,4)/2889 , 1084 , 3667 , 1833/
      DATA (mm(71,j),j=1,4)/2857 , 3438 , 2703 , 2877/
      DATA (mm(72,j),j=1,4)/2094 , 2408 , 629 , 1633/
      DATA (mm(73,j),j=1,4)/1818 , 1589 , 2365 , 981/
      DATA (mm(74,j),j=1,4)/688 , 2391 , 2431 , 2009/
      DATA (mm(75,j),j=1,4)/1407 , 288 , 1113 , 941/
      DATA (mm(76,j),j=1,4)/634 , 26 , 3922 , 2449/
      DATA (mm(77,j),j=1,4)/3231 , 512 , 2554 , 197/
      DATA (mm(78,j),j=1,4)/815 , 1456 , 184 , 2441/
      DATA (mm(79,j),j=1,4)/3524 , 171 , 2099 , 285/
      DATA (mm(80,j),j=1,4)/1914 , 1677 , 3228 , 1473/
      DATA (mm(81,j),j=1,4)/516 , 2657 , 4012 , 2741/
      DATA (mm(82,j),j=1,4)/164 , 2270 , 1921 , 3129/
      DATA (mm(83,j),j=1,4)/303 , 2587 , 3452 , 909/
      DATA (mm(84,j),j=1,4)/2144 , 2961 , 3901 , 2801/
      DATA (mm(85,j),j=1,4)/3480 , 1970 , 572 , 421/
      DATA (mm(86,j),j=1,4)/119 , 1817 , 3309 , 4073/
      DATA (mm(87,j),j=1,4)/3357 , 676 , 3171 , 2813/
      DATA (mm(88,j),j=1,4)/837 , 1410 , 817 , 2337/
      DATA (mm(89,j),j=1,4)/2826 , 3723 , 3039 , 1429/
      DATA (mm(90,j),j=1,4)/2332 , 2803 , 1696 , 1177/
      DATA (mm(91,j),j=1,4)/2089 , 3185 , 1256 , 1901/
      DATA (mm(92,j),j=1,4)/3780 , 184 , 3715 , 81/
      DATA (mm(93,j),j=1,4)/1700 , 663 , 2077 , 1669/
      DATA (mm(94,j),j=1,4)/3712 , 499 , 3019 , 2633/
      DATA (mm(95,j),j=1,4)/150 , 3784 , 1497 , 2269/
      DATA (mm(96,j),j=1,4)/2000 , 1631 , 1101 , 129/
      DATA (mm(97,j),j=1,4)/3375 , 1925 , 717 , 1141/
      DATA (mm(98,j),j=1,4)/1621 , 3912 , 51 , 249/
      DATA (mm(99,j),j=1,4)/3090 , 1398 , 981 , 3917/
      DATA (mm(100,j),j=1,4)/3765 , 1349 , 1978 , 2481/
      DATA (mm(101,j),j=1,4)/1149 , 1441 , 1813 , 3941/
      DATA (mm(102,j),j=1,4)/3146 , 2224 , 3881 , 2217/
      DATA (mm(103,j),j=1,4)/33 , 2411 , 76 , 2749/
      DATA (mm(104,j),j=1,4)/3082 , 1907 , 3846 , 3041/
      DATA (mm(105,j),j=1,4)/2741 , 3192 , 3694 , 1877/
      DATA (mm(106,j),j=1,4)/359 , 2786 , 1682 , 345/
      DATA (mm(107,j),j=1,4)/3316 , 382 , 124 , 2861/
      DATA (mm(108,j),j=1,4)/1749 , 37 , 1660 , 1809/
      DATA (mm(109,j),j=1,4)/185 , 759 , 3997 , 3141/
      DATA (mm(110,j),j=1,4)/2784 , 2948 , 479 , 2825/
      DATA (mm(111,j),j=1,4)/2202 , 1862 , 1141 , 157/
      DATA (mm(112,j),j=1,4)/2199 , 3802 , 886 , 2881/
      DATA (mm(113,j),j=1,4)/1364 , 2423 , 3514 , 3637/
      DATA (mm(114,j),j=1,4)/1244 , 2051 , 1301 , 1465/
      DATA (mm(115,j),j=1,4)/2020 , 2295 , 3604 , 2829/
      DATA (mm(116,j),j=1,4)/3160 , 1332 , 1888 , 2161/
      DATA (mm(117,j),j=1,4)/2785 , 1832 , 1836 , 3365/
      DATA (mm(118,j),j=1,4)/2772 , 2405 , 1990 , 361/
      DATA (mm(119,j),j=1,4)/1217 , 3638 , 2058 , 2685/
      DATA (mm(120,j),j=1,4)/1822 , 3661 , 692 , 3745/
      DATA (mm(121,j),j=1,4)/1245 , 327 , 1194 , 2325/
      DATA (mm(122,j),j=1,4)/2252 , 3660 , 20 , 3609/
      DATA (mm(123,j),j=1,4)/3904 , 716 , 3285 , 3821/
      DATA (mm(124,j),j=1,4)/2774 , 1842 , 2046 , 3537/
      DATA (mm(125,j),j=1,4)/997 , 3987 , 2107 , 517/
      DATA (mm(126,j),j=1,4)/2573 , 1368 , 3508 , 3017/
      DATA (mm(127,j),j=1,4)/1148 , 1848 , 3525 , 2141/
      DATA (mm(128,j),j=1,4)/545 , 2366 , 3801 , 1537/
!     ..
!     .. Executable Statements ..
!
      i1 = Iseed(1)
      i2 = Iseed(2)
      i3 = Iseed(3)
      i4 = Iseed(4)
!
      DO i = 1 , MIN(N,LV)
         DO
!
!
!        Multiply the seed by i-th power of the multiplier modulo 2**48
!
            it4 = i4*mm(i,4)
            it3 = it4/IPW2
            it4 = it4 - IPW2*it3
            it3 = it3 + i3*mm(i,4) + i4*mm(i,3)
            it2 = it3/IPW2
            it3 = it3 - IPW2*it2
            it2 = it2 + i2*mm(i,4) + i3*mm(i,3) + i4*mm(i,2)
            it1 = it2/IPW2
            it2 = it2 - IPW2*it1
            it1 = it1 + i1*mm(i,4) + i2*mm(i,3) + i3*mm(i,2)            &
     &            + i4*mm(i,1)
            it1 = MOD(it1,IPW2)
!
!        Convert 48-bit integer to a real number in the interval (0,1)
!
            X(i) = R*(REAL(it1)+R*(REAL(it2)+R*(REAL(it3)+R*REAL(it4))))
!
            IF ( X(i)==1.0 ) THEN
!           If a real number has n bits of precision, and the first
!           n bits of the 48-bit integer above happen to be all 1 (which
!           will occur about once every 2**n calls), then X( I ) will
!           be rounded to exactly 1.0. In IEEE single precision arithmetic,
!           this will happen relatively often since n = 24.
!           Since X( I ) is not supposed to return exactly 0.0 or 1.0,
!           the statistically correct thing to do in this situation is
!           simply to iterate again.
!           N.B. the case X( I ) = 0.0 should not be possible.
               i1 = i1 + 2
               i2 = i2 + 2
               i3 = i3 + 2
               i4 = i4 + 2
               CYCLE
            ENDIF
            EXIT
         ENDDO
!
      ENDDO
!
!     Return final value of seed
!
      Iseed(1) = it1
      Iseed(2) = it2
      Iseed(3) = it3
      Iseed(4) = it4
!
!     End of SLARUV
!
      END SUBROUTINE SLARUV


text:
!>  Partitioning a graph with non-equal number of edges per node.
!
program test_PartGraphRecursive2

    use metis_interface, only: idx_t, METIS_PartGraphRecursive, METIS_SetDefaultOptions, &
        METIS_NOPTIONS, METIS_OPTION_NUMBERING, METIS_OK

    integer(idx_t) :: n, m
    integer(idx_t), allocatable :: xadj(:), adjncy(:), part(:)
    integer(idx_t) :: options(0:METIS_NOPTIONS-1), ios, ncon, objval, nparts

    write(*,'(A)') "TEST METIS_PartGraphRecursive 2"

    !   1---5
    !   |\ / \
    !   | 3   6
    !   |/ \ / \
    !   2---4---7

    xadj = [1,4,7,11,15,18,21,23]
    adjncy = [5,3,2,1,3,4,5,4,2,1,2,3,6,7,1,3,6,5,4,7,6,4]

    n = size(xadj) - 1
    m = size(adjncy)/2
    write(*,'(A,I0)') "n = ", n
    write(*,'(A,I0)') "m = ", m

    ios = METIS_SetDefaultOptions(options)
    if (ios /= METIS_OK) then
        write(*,*) "METIS_SetDefaultOptions failed with error: ", ios
        error stop 1
    end if
    options(METIS_OPTION_NUMBERING) = 1     ! Fortran-style numbering

    ncon = 1
    nparts = 2
    allocate(part(n))
    ios = METIS_PartGraphRecursive(n,ncon,xadj,adjncy,&
                nparts=nparts,objval=objval,part=part,options=options)
    if (ios /= METIS_OK) then
        write(*,*) "METIS_PartGraphKway failed with error: ", ios
        error stop 1
    end if

    write(*,'(A,I0)') "objval = ", objval
    write(*,'(A,*(I1,:,1X))') "part = ", part

end program

max_stars_repo_name: dylanrubini/fmetis

max_issues_repo_name: dylanrubini/fmetis

max_forks_repo_name: dylanrubini/fmetis


text:
program exp_distr_program
    implicit none
    integer, parameter :: NR_VALUES = 1000
    integer :: i

    do i = 1, NR_VALUES
        print *, get_exp_rand()
    end do

contains

    function get_exp_rand() result(val)
        implicit none
        real :: val
        real, save :: next_val
        logical, save :: is_next_set = .false.
        real :: x, y, radius, factor

        if (is_next_set) then
            is_next_set = .false.
            val = next_val
        else
            do while (.not. is_next_set)
                call random_number(x)
                call random_number(y)
                x = 2.0*x - 1.0
                y = 2.0*y - 1.0
                radius = x**2 + y**2
                if (0.0 < radius .and. radius <= 1.0) exit
            end do
            factor = sqrt(-2.0*log(radius)/radius)
            val = x*factor
            next_val =y*factor
            is_next_set = .true.
        end if
    end function get_exp_rand

end program exp_distr_program

max_issues_repo_name: caguerra/Fortran-MOOC

max_forks_repo_name: caguerra/Fortran-MOOC


!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!  unit testing module for mod_optic
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

module mod_utest
   use mod_optic
   use mod_ioset

   !for 3-mirror test and 1-slab test
   real(R_KD) n_mr(3,3)
   complex(R_KD) E_vec_in(3,2)

contains

    subroutine mcr_test

       select case (nutdbg)

          case(1)
             call ut_flec
          case(2)
             call ut_3_mirror
          case(3)
             call ut_1_slab
      end select
      call StopMe

    end subroutine mcr_test

    subroutine ut_flec

       implicit none
       real(R_KD) dir_i_ut(3), dir_r_ut(3),dir_n_ut(3)
       real(R_KD) r_n_ut

       dir_n_ut=[0d0,0d0,1d0]
       dir_cos=[0d0,0d0,-1d0]

       call GetFlecDir(dir_n_ut)

       write(*,"('dir_r_ut=', 3(f8.4,x))") dir_cos

       dir_cos=[dsin(pi/3d0), 0d0, -dcos(pi/3d0)]
       call GetFlecDir(dir_n_ut)
       write(*,"('dir_r_ut=', 3(f8.4,x))") dir_cos

       r_n_ut=0.8d0
       dir_cos=[dsin(pi/3d0), 0d0, -dcos(pi/3d0)]
       call GetFracDir(dir_n_ut, r_n_ut)
       write(*,"('dir_r_ut=', 3(f8.4,x))") dir_cos

    end subroutine

    subroutine ut_3_mirror

        integer i, j, k
        real(R_KD) rval,r_n
        logical :: is_reflected

        rval=dsqrt(2.0d0)/2d0

        !set-up the mirror normal direction
        n_mr(:,1)=[rval,0d0, -rval]
        n_mr(:,2)=[-rval, rval,0d0]
        n_mr(:,3)=[0d0,-rval, rval]

        !incident light
        dir_cos=[0d0,0d0,1d0]
        E_vec_in(:,1)=[1d0,0d0,0d0]  !x_polarized
        E_vec_in(:,2)=[0d0,1d0,0d0]  !y_polarized

        r_n=0d0

        do i=1, 2  ! two tests: x and y polarized incident lights
           E_vec=E_vec_in(:,i)
           do j=1,3  !3-mirror
              is_reflected=IsFlec(n_mr(:,j),r_n)
              if (.not. is_reflected) then
                 write(disp_msg(1), "('Utest failure: 3-mirror test light not reflected')")
                 call ShowMsg(1)
              endif
           enddo
           write(disp_msg(1),"('3-mirror test Part ', I0)") i
           write(disp_msg(2),"('Incident polarization E=', 3('(',ES12.5,x,ES12.5,')') )") &
                    (realpart(E_vec_in(k,i)), imagpart(E_vec_in(k,i)), k=1,3)
           write(disp_msg(3),"('Outgoing polarization E=', 3('(',ES12.5,x,ES12.5,')') )") &
                    (realpart(E_vec(k)), imagpart(E_vec(k)), k=1,3)
           write(disp_msg(4),"('     ')")
           call ShowMsg(4)
        enddo

    end subroutine ut_3_mirror

    subroutine ut_1_slab

        integer i, j, k,ei
        real(R_KD) rval,r_n
        integer :: counter(3)=0
        logical :: is_reflected, inSlab

        rval=dsqrt(2.0d0)/2d0

        !set-up the mirror normal direction
        n_mr(:,1)=[0d0,0d0, -1d0]
        n_mr(:,2)=[0d0,0d0,1d0]
        n_mr(:,3)=[0d0,-rval, rval]

        !incident light
        E_vec_in(:,1)=[1d0,0d0,0d0]  !x_polarized
        E_vec_in(:,2)=[0d0,1d0,0d0]  !y_polarized

        ei=1
        counter=0
        do i=1, num_p
           E_vec=E_vec_in(:,ei)
           r_n=n_frac(1)/n_frac(2)
           dir_cos=[0d0,rval ,rval]
           inSlab=.false.
           is_reflected=IsFlec(n_mr(:,1),r_n)
           if (is_reflected) then
               counter(3)=counter(3)+1
           else
               inSlab=.true.
           endif
           k=0
           do while (inSlab)
             j=mod(k,2)+1  !j=1: upper surface; j=2: lower surface
             r_n=n_frac(2)/n_frac(1)
             inSlab=IsFlec(n_mr(:,j), r_n)
             if (.not. inSlab) then
                counter(j)=counter(j)+1
             endif
             k=k+1
           enddo
        enddo

        write(disp_msg(1),"('1-slab test: tot_particle=',I0)") num_p
        write(disp_msg(2),"('Incident polarization E=', 3('(',ES12.5,x,ES12.5,')') )") &
                    (realpart(E_vec_in(k,ei)), imagpart(E_vec_in(k,ei)), k=1,3)
        write(disp_msg(3),"('  Reflected   : ', I0 )") counter(3)
        write(disp_msg(4),"('  backscatted : ', I0 )") counter(2)

        write(disp_msg(5),"('  Passed      : ', I0 )") counter(1)
        call ShowMsg(5)
    end subroutine ut_1_slab

end module mod_utest


program main
    implicit none
    integer :: clock, n
    integer, allocatable :: seed(:)
    integer(4) :: i, n_tot = 100000000
    real(8), parameter :: l = 0.d0, r = 1.d0
    real(8) :: x, y
    real(8) :: integral = 0.d0

    call SYSTEM_CLOCK(clock)
    call RANDOM_SEED(size = n)
    allocate(seed(n))
    do i = 1, n
        seed(i) = clock + 37 * i
    end do
    call RANDOM_SEED(put = seed)
    do i = 1, n_tot
        call RANDOM_NUMBER(x)
        x = x * (r - l) + l
        call func(x, y)
        integral = integral + y
    end do
    integral = (r - l) / dble(n_tot) * integral
    write(*,'(f10.5)') integral
end program main

subroutine func(x, y)
    ! the function to be integrated
    implicit none
    real(8), intent(in) :: x
    real(8), intent(out) :: y

    y = x**2
end subroutine func


text:
module kernel_share_mod

use constants_mod

!! ------------------------------------------------------------------------- !!
!! AMR operators

real(dp), save :: rfaceMat(N_NODES,N_NODES,4)
real(dp), save :: rfaceTam(N_NODES,N_NODES,4)

real(dp), save :: cfaceMat(N_NODES,N_NODES,4)
real(dp), save :: cfaceTam(N_NODES,N_NODES,4)

real(dp), save :: refineMat(N_NODES,N_NODES,3,8)
real(dp), save :: coarseMat(N_NODES,N_NODES,3,8)

real(dp), save :: refineTam(N_NODES,N_NODES,3,8)
real(dp), save :: coarseTam(N_NODES,N_NODES,3,8)

!! ------------------------------------------------------------------------- !!

contains

subroutine hook_init_400

    use linalg_mod, only: shiftmatrix

    logical :: ok

    integer :: i,j,k,n,q,f,p,o, d ,qq

    real(dp) :: tempMat(N_NODES,N_NODES)

    !! --------------------------------------------------------------------- !!

    rfaceMat = 0.0_dp

    tempMat = 0.0_dp
    tempMat(1,1) = 1.0_dp
    tempMat(2,1) = 1.0_dp

    do f = 1,4
        do i = 1,N_NODES/2
            rfaceMat(:,:,f) = rfaceMat(:,:,f) + shiftMatrix(2*(i-1),(f-1)*N_NODES/2+(i-1),tempMat)
        end do
    end do

    do f = 1,4
        cfaceMat(:,:,f) = 0.5_dp*transpose(rfaceMat(:,:,f))
    end do

    rfaceMat(:,:,3) = rfaceMat(:,:,1)
    rfaceMat(:,:,4) = rfaceMat(:,:,2)

    rfaceTam(:,:,1) = transpose(rfaceMat(:,:,1))
    rfaceTam(:,:,2) = transpose(rfaceMat(:,:,1))
    rfaceTam(:,:,3) = transpose(rfaceMat(:,:,2))
    rfaceTam(:,:,4) = transpose(rfaceMat(:,:,2))

    cfaceMat(:,:,3) = cfaceMat(:,:,1)
    cfaceMat(:,:,4) = cfaceMat(:,:,2)

    cfaceTam(:,:,1) = transpose(cfaceMat(:,:,1))
    cfaceTam(:,:,2) = transpose(cfaceMat(:,:,1))
    cfaceTam(:,:,3) = transpose(cfaceMat(:,:,2))
    cfaceTam(:,:,4) = transpose(cfaceMat(:,:,2))

    refineMat(:,:,1,1) = rfaceMat(:,:,1)
    refineMat(:,:,2,1) = rfaceMat(:,:,1)
    refineMat(:,:,3,1) = rfaceMat(:,:,1)

    refineMat(:,:,1,2) = rfaceMat(:,:,2)
    refineMat(:,:,2,2) = rfaceMat(:,:,1)
    refineMat(:,:,3,2) = rfaceMat(:,:,1)

    refineMat(:,:,1,3) = rfaceMat(:,:,1)
    refineMat(:,:,2,3) = rfaceMat(:,:,2)
    refineMat(:,:,3,3) = rfaceMat(:,:,1)

    refineMat(:,:,1,4) = rfaceMat(:,:,2)
    refineMat(:,:,2,4) = rfaceMat(:,:,2)
    refineMat(:,:,3,4) = rfaceMat(:,:,1)

    refineMat(:,:,1,5) = rfaceMat(:,:,1)
    refineMat(:,:,2,5) = rfaceMat(:,:,1)
    refineMat(:,:,3,5) = rfaceMat(:,:,2)

    refineMat(:,:,1,6) = rfaceMat(:,:,2)
    refineMat(:,:,2,6) = rfaceMat(:,:,1)
    refineMat(:,:,3,6) = rfaceMat(:,:,2)

    refineMat(:,:,1,7) = rfaceMat(:,:,1)
    refineMat(:,:,2,7) = rfaceMat(:,:,2)
    refineMat(:,:,3,7) = rfaceMat(:,:,2)

    refineMat(:,:,1,8) = rfaceMat(:,:,2)
    refineMat(:,:,2,8) = rfaceMat(:,:,2)
    refineMat(:,:,3,8) = rfaceMat(:,:,2)

    coarseMat(:,:,1,1) = cfaceMat(:,:,1)
    coarseMat(:,:,2,1) = cfaceMat(:,:,1)
    coarseMat(:,:,3,1) = cfaceMat(:,:,1)

    coarseMat(:,:,1,2) = cfaceMat(:,:,2)
    coarseMat(:,:,2,2) = cfaceMat(:,:,1)
    coarseMat(:,:,3,2) = cfaceMat(:,:,1)

    coarseMat(:,:,1,3) = cfaceMat(:,:,1)
    coarseMat(:,:,2,3) = cfaceMat(:,:,2)
    coarseMat(:,:,3,3) = cfaceMat(:,:,1)

    coarseMat(:,:,1,4) = cfaceMat(:,:,2)
    coarseMat(:,:,2,4) = cfaceMat(:,:,2)
    coarseMat(:,:,3,4) = cfaceMat(:,:,1)

    coarseMat(:,:,1,5) = cfaceMat(:,:,1)
    coarseMat(:,:,2,5) = cfaceMat(:,:,1)
    coarseMat(:,:,3,5) = cfaceMat(:,:,2)

    coarseMat(:,:,1,6) = cfaceMat(:,:,2)
    coarseMat(:,:,2,6) = cfaceMat(:,:,1)
    coarseMat(:,:,3,6) = cfaceMat(:,:,2)

    coarseMat(:,:,1,7) = cfaceMat(:,:,1)
    coarseMat(:,:,2,7) = cfaceMat(:,:,2)
    coarseMat(:,:,3,7) = cfaceMat(:,:,2)

    coarseMat(:,:,1,8) = cfaceMat(:,:,2)
    coarseMat(:,:,2,8) = cfaceMat(:,:,2)
    coarseMat(:,:,3,8) = cfaceMat(:,:,2)

    do f = 1,8
    do d = 1,3
        refineTam(:,:,d,f) = transpose(refineMat(:,:,d,f))
        coarseTam(:,:,d,f) = transpose(coarseMat(:,:,d,f))
    end do
    end do

end subroutine

end module

max_stars_repo_name: jmark/nemo

max_issues_repo_name: jmark/nemo

max_forks_repo_name: jmark/nemo


text:
The Tech Exchange has been a popular European high street shop for years and the concept is simple. Preowned tech gear is purchased for cash then refurbished, updated, tested and offered back to the local Sacramento community as a cost effective affordable solution.

Visit our West Sacramento store for a wide selection of tech including laptop computers, desktop computers, computer monitors, digital cameras, smart phones, video games consoles, video games, dj gear, turntables, cdjs, studio gear, pro audio, home audio, cables, accessories and more!

Need a little extra cash? Bring in your gear for evaluation for cash today!
Want to make a trade? Bring in what you have and trade it for something you want!
Need a solid cost effective laptop or desktop computer to surf the web or do homework? We have a huge selection of refurbished computers in stock! All Pcs come with Open Office and Anti Virus! The majority of our laptops are priced between $150 and $250!
Need service or repair for your desktop or laptop computer? Our West Sacramento service department has you covered. Click here for details.
Dont want to spend a fortune on a cable? We have a wide selection of computer and audio cables including HDMI cables for $10!
Have a video game you are bored of? Bring it in and trade it for a different title.
The bottom line is if you dont want to spend a fortune on tech then visit the Exchange!


text:
      subroutine qqb_wtbndk(p,msq)
************************************************************************
*    W+t+b production, with massive b                                  *
*    Matrix element squared and averaged over initial colours and spins*
*--- W+t production (nwz=-1)                                           *
*     q(-p1) + qbar(-p2) --> W + t(p5) + b(p6)                         *
*                            |                                         *
*                            |                                         *
*                            |                                         *
*                            --> e^-(p3) + nubar(p4)                   *
************************************************************************
      implicit none
      include 'constants.f'
      include 'masses.f'
      include 'qcdcouple.f'
      include 'ewcouple.f'
      integer ia,ig,ib,it,j,k
      double precision p(mxpart,4),msq(-nf:nf,-nf:nf),fac,dot,prop,
     . msq_gg
      double complex ampgg_ga(2,2,2,2),ampgg_ag(2,2,2,2)

      do j=-nf,nf
      do k=-nf,nf
      msq(j,k)=0d0
      enddo
      enddo

      fac=gsq**2*gwsq**2

      prop=(two*dot(p,3,4)-wmass**2)**2+(wmass*wwidth)**2

c--- calculate amplitudes for g+g in initial state
      call BBamps_nores(p,1,2,3,4,5,6,ampgg_ag)
      call BBamps_nores(p,2,1,3,4,5,6,ampgg_ga)
c      call BBamps(p,1,2,3,4,5,6,ampgg_ag)
c      call BBamps(p,2,1,3,4,5,6,ampgg_ga)

      msq_gg=0d0
c--- sum over helicities of gluons and massive quarks
      do ib=1,2
      do it=1,2
      do ia=1,2
      do ig=1,2
        msq_gg=msq_gg+xn*cf**2*(
     .        + abs(ampgg_ag(ia,ig,ib,it))**2
     .        + abs(ampgg_ga(ig,ia,ib,it))**2
     .        - one/xn/cf*dble(ampgg_ag(ia,ig,ib,it)
     .                 *dconjg(ampgg_ga(ig,ia,ib,it))))
      enddo
      enddo
      enddo
      enddo

      msq_gg=msq_gg/prop

      do j=-nf,nf,nf
      do k=-nf,nf,nf
      msq(j,k)=0d0
      if ((j .eq. 0) .and. (k .eq. 0)) then
          msq(j,k)=avegg*fac*msq_gg
      endif
      enddo
      enddo

      return
      end

max_stars_repo_name: tmartini/JHUGen

max_issues_repo_name: tmartini/JHUGen

max_forks_repo_name: tmartini/JHUGen