text:
C
C       $Id$
C
      SUBROUTINE SHGETNP3 (PX,PY,PZ,X,Y,Z,NR,LCELL,LNEXT,
     +                    XYZMIN,XYZDEL, NP,DSQ)
      INTEGER NR, LCELL(NR,NR,NR), LNEXT(1), NP
      REAL    PX, PY, PZ, X(1), Y(1), Z(1), XYZMIN(3),
     +        XYZDEL(3), DSQ
C
C***********************************************************
C
C                                               ROBERT RENKA
C                                       Univ. of North Texas
C                                             (817) 565-2767
C
C  Given a set of N nodes and the data structure defined in
C  subroutine SHSTORE3, this subroutine uses the cell method to
C  find the closest unmarked node NP to a specified point P.
C  NP is then marked by setting LNEXT(NP) to -LNEXT(NP).  (A
C  node is marked if and only if the corresponding LNEXT element
C  is negative.  The absolute values of LNEXT elements,
C  however, must be preserved.)  Thus, the closest M nodes to
C  P may be determined by a sequence of M calls to this routine
C  Note that if the nearest neighbor to node K is to be
C  determined (PX = X(K), PY = Y(K), and PZ = Z(K)), then
C  K should be marked before the call to this routine.
C
C  The search is begun in the cell containing (or closest
C  to) P and proceeds outward in box-shaped layers until all
C  cells which contain points within distance R of P have
C  been searched, where R is the distance from P to the first
C  unmarked node encountered (infinite if no unmarked nodes
C  are present).
C
C  On input --
C
C    PX,PY,PZ = Cartesian coordinates of the point P whose
C               nearest unmarked neighbor is to be found.
C
C    X,Y,Z = Arrays of length N, for N .GE. 2, containing
C            the Cartesian coordinates of the nodes.
C
C    NR = Number of rows, columns, and planes in the cell
C         grid.  NR .GE. 1.
C
C    LCELL = NR by NR by NR array of nodal indices associated
C            with cells.
C
C    LNEXT = Array of length N containing next-node indices
C            (or their negatives).
C
C    XYZMIN,XYZDEL = Arrays of length 3 containing minimum
C                    nodal coordinates and cell dimensions,
C                    respectively.  XYZDEL elements must
C                    be positive.
C
C  Input parameters other than LNEXT are not altered by
C  this routine.  With the exception of (PX,PY,PZ) and the
C  signs of LNEXT elements, these parameters should be
C  unaltered from their values on output from subroutine SHSTORE3.
C
C  On output --
C
C    NP = index (for X, Y, and Z) of the nearest unmarked
C         node to P, or 0 if all nodes are marked or NR
C         .LT. 1 or an element of XYZDEL is not positive.
C         LNEXT(NP) .LT. 0 if NP .NE. 0.
C
C    DSQ = Squared Euclidean distance between P and node
C          NP, or 0 if NP = 0.
C
C  Modules required by SHGETNP3 -- none
C
C  Intrinsic functions called by SHGETNP3 -- IABS, IFIX, SQRT
C
C***********************************************************
C
      LOGICAL FIRST
      XP = PX
      YP = PY
      ZP = PZ
      DX = XYZDEL(1)
      DY = XYZDEL(2)
      DZ = XYZDEL(3)
C
C  Test for invalid input parameters.
C
      IF (NR .LT. 1  .OR.  DX .LE. 0.  .OR.  DY .LE. 0.
     +    .OR.  DZ .LE. 0.) GO TO 10
C
C Initialize parameters --
C
C   FIRST = TRUE iff the first unmarked node has yet to be
C           encountered,
C   IMIN,...,KMAX = cell indices defining the range of the
C                   search,
C   DELX,DELY,DELZ = PX-XYZMIN(1), PY-XYZMIN(2), and
C                    PZ-XYZMIN(3),
C   I0,J0,K0 = cell containing or closest to P,
C   I1,...,K2 = cell indices of the layer whose intersection
C               with the range defined by IMIN,...,KMAX is
C               currently being searched.
C
      FIRST = .TRUE.
      IMIN = 1
      IMAX = NR
      JMIN = 1
      JMAX = NR
      KMIN = 1
      KMAX = NR
      DELX = XP - XYZMIN(1)
      DELY = YP - XYZMIN(2)
      DELZ = ZP - XYZMIN(3)
      I0 = IFIX(DELX/DX) + 1
      IF (I0 .LT. 1) I0 = 1
      IF (I0 .GT. NR) I0 = NR
      J0 = IFIX(DELY/DY) + 1
      IF (J0 .LT. 1) J0 = 1
      IF (J0 .GT. NR) J0 = NR
      K0 = IFIX(DELZ/DZ) + 1
      IF (K0 .LT. 1) K0 = 1
      IF (K0 .GT. NR) K0 = NR
      I1 = I0
      I2 = I0
      J1 = J0
      J2 = J0
      K1 = K0
      K2 = K0
C
C  Outer loop on layers, inner loop on layer cells, excluding
C  those outside the range (IMIN,IMAX) x (JMIN,JMAX) x (KMIN,KMAX).
C
    1 DO 7 K = K1,K2
        IF (K .GT. KMAX) GO TO 8
        IF (K .LT. KMIN) GO TO 7
        DO 6 J = J1,J2
          IF (J .GT. JMAX) GO TO 7
          IF (J .LT. JMIN) GO TO 6
          DO 5 I = I1,I2
            IF (I .GT. IMAX) GO TO 6
            IF (I .LT. IMIN) GO TO 5
            IF (K .NE. K1  .AND.  K .NE. K2  .AND.
     +          J .NE. J1  .AND.  J .NE. J2  .AND.
     +          I .NE. I1  .AND.  I .NE. I2) GO TO 5
C
C  Search cell (I,J,K) for unmarked nodes L.
C
            L = LCELL(I,J,K)
            IF (L .EQ. 0) GO TO 5
C
C  Loop on nodes in cell (I,J,K).
C
    2       LN = LNEXT(L)
            IF (LN .LT. 0) GO TO 4
C
C  Node L is not marked.
C
            RSQ = (X(L)-XP)**2 + (Y(L)-YP)**2 + (Z(L)-ZP)**2
            IF (.NOT. FIRST) GO TO 3
C
C  Node L is the first unmarked neighbor of P encountered.
C  Initialize LMIN to the current candidate for NP, and
C  RSMIN to the squared distance from P to LMIN.  IMIN,
C  IMAX, JMIN, JMAX, KMIN, and KMAX are updated to define
C  the smallest rectangle containing a sphere of radius
C  R = SQRT(RSMIN) centered at P, and contained in (1,NR)
C  x (1,NR) x (1,NR) (except that, if P is outside the
C  box defined by the nodes, it is possible that IMIN
C  .GT. NR or IMAX .LT. 1, etc.).  FIRST is reset to FALSE.
C
            LMIN = L
            RSMIN = RSQ
            R = SQRT(RSMIN)
            IMIN = IFIX((DELX-R)/DX) + 1
            IF (IMIN .LT. 1) IMIN = 1
            IMAX = IFIX((DELX+R)/DX) + 1
            IF (IMAX .GT. NR) IMAX = NR
            JMIN = IFIX((DELY-R)/DY) + 1
            IF (JMIN .LT. 1) JMIN = 1
            JMAX = IFIX((DELY+R)/DY) + 1
            IF (JMAX .GT. NR) JMAX = NR
            KMIN = IFIX((DELZ-R)/DZ) + 1
            IF (KMIN .LT. 1) KMIN = 1
            KMAX = IFIX((DELZ+R)/DZ) + 1
            IF (KMAX .GT. NR) KMAX = NR
            FIRST = .FALSE.
            GO TO 4
C
C  Test for node L closer than LMIN to P.
C
    3       IF (RSQ .GE. RSMIN) GO TO 4
C
C  Update LMIN and RSMIN.
C
            LMIN = L
            RSMIN = RSQ
C
C  Test for termination of loop on nodes in cell (I,J,K).
C
    4       IF (IABS(LN) .EQ. L) GO TO 5
            L = IABS(LN)
            GO TO 2
    5       CONTINUE
    6     CONTINUE
    7   CONTINUE
C
C  Test for termination of loop on cell layers.
C
    8 IF (I1 .LE. IMIN  .AND.  I2 .GE. IMAX  .AND.
     .    J1 .LE. JMIN  .AND.  J2 .GE. JMAX  .AND.
     .    K1 .LE. KMIN  .AND.  K2 .GE. KMAX) GO TO 9
      I1 = I1 - 1
      I2 = I2 + 1
      J1 = J1 - 1
      J2 = J2 + 1
      K1 = K1 - 1
      K2 = K2 + 1
      GO TO 1
C
C  Unless no unmarked nodes were encountered, LMIN is the
C  closest unmarked node to P.
C
    9 IF (FIRST) GO TO 10
      NP = LMIN
      DSQ = RSMIN
      LNEXT(LMIN) = -LNEXT(LMIN)
      RETURN
C
C  Error -- NR or XYZDEL is invalid or all nodes are marked.
C
   10 NP = 0
      DSQ = 0.
C
      RETURN
      END

ext: f

max_stars_repo_name: xylar/cdat

max_issues_repo_name: xylar/cdat

max_forks_repo_name: CDAT/uvcdat


text:
c
      Subroutine  FOMCLV(IND,NPOINT,AA,B,OM,RS,CL,CS)
c
C     11 JAN 78   ...JJL...       INITIAL RELEASE
c     15 Jan 93   ...W.P. Lee...  Ported for UNIX Conversion
C
C     SET UP RETICLE DESCRIPTIVE ARRAY A AND CALL FOMCAL
C
C     AA CONTAINS NPOINT ENTRYS OF
C          AA(1,J) LINE VALUE OF JTH RETICLE
C          AA(2,J) SAMPLE VALUE
C          AA(3,J) LATITUDE
C          AA(4,J) LONGITUDE
C
      Real*4  AA(4,20),A(5,20),B(6)
      Real*8  OM(3,3)
      Real*8  RS(3)
      Real*4  CL, CS
c
      Do 100 I = 1, NPOINT
       A(1,I) = AA(1,I)
       A(2,I) = AA(2,I)
       A(3,I) =  B(3)
       A(4,I) = AA(3,I)
       A(5,I) = AA(4,I)
100   Continue
C
      Call FOMCAV(IND,NPOINT,A,B,OM,RS,CL,CS)
c
      Return
      End

max_stars_repo_name: NASA-AMMOS/VICAR

max_issues_repo_name: NASA-AMMOS/VICAR

max_forks_repo_name: NASA-AMMOS/VICAR


text:
        function func result(rst)
        end

max_stars_repo_name: xevolver/xevxml

max_issues_repo_name: xevolver/xevxml

max_forks_repo_name: xevolver/xevxml


text:
! { dg-do compile }
! { dg-options "-std=legacy" }
!
! Testing g77 intrinsics as subroutines
      integer(kind=8) i8
      integer i4
      integer i
      character*80 c

      i8 = time ()
      i4 = time ()
      i8 = time8 ()
      i4 = time8 ()

      i8 = hostnm (c)
      i4 = hostnm (c)
      i = hostnm (c)

      i8 = ierrno ()
      i4 = ierrno ()
      i = ierrno ()

      i8 = kill (i8, i8)
      i8 = kill (i8, i4)
      i8 = kill (i4, i8)
      i8 = kill (i4, i4)
      i4 = kill (i8, i8)
      i4 = kill (i8, i4)
      i4 = kill (i4, i8)
      i4 = kill (i4, i4)

      i8 = link ('foo', 'bar')
      i4 = link ('foo', 'bar')
      i = link ('foo', 'bar')

      i8 = rename ('foo', 'bar')
      i4 = rename ('foo', 'bar')
      i = rename ('foo', 'bar')

      i8 = symlnk ('foo', 'bar')
      i4 = symlnk ('foo', 'bar')
      i = symlnk ('foo', 'bar')

! Cleaning our mess
      call unlink ('bar')

! This should be the last test, unless you want garbage everywhere in
! your filesystem.
      i8 = chdir ('..')
      i4 = chdir ('..')
      i = chdir ('..')

      end

max_issues_repo_name: sujankh/rose-matlab

max_forks_repo_name: sujankh/rose-matlab


text:
    character(len = *), intent(in), optional :: label
    logical :: failure

    if (present(tol)) then
       failure = abs(val - expected) >= tol
    else
       failure = (val /= expected)
    end if
    if (failure) then
       if (present(label)) then
          call f_err_throw(label // ": expected " // yaml_toa(expected) // ", value " // yaml_toa(val))
       else
          call f_err_throw("expected " // yaml_toa(expected) // ", value " // yaml_toa(val))
       end if
    end if

max_stars_repo_name: adegomme/aiida-bigdft-plugin


text:
!!# MODULE <<PAR_MEMORY>>
MODULE PAR_MEMORY

!!## PURPOSE
!! Define the number of bytes per element for
!! the intrinsic types and overhead for pointers and
!! arrays.

!!## DETAILS
!! Used for the <MEMORY> and <MEMORYn> functions.

!!## HISTORY
! 1.[waw] Author   = William A. Wieselquist
!         Modified = 104.2006
!         Contact  = william.wieselquist AT gmail.com

!!## DEFAULT IMPLICIT
IMPLICIT NONE

!!## DEFAULT ACCESS
PUBLIC

!!## EXTERNAL PARAMETERS
!! * bytes per element of intrinsic types
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_R   = 4
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_Rsp = 4
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_Rdp = 8
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_C   = 2*BYTES_PER_ELEMENT_R
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_Csp = 2*BYTES_PER_ELEMENT_Rsp
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_Cdp = 2*BYTES_PER_ELEMENT_Rdp
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_I1  = 1
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_I2  = 2
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_I   = 4
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_I4  = 4
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_I8  = 8
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_L   = 4
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_L1  = 1
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_L2  = 2
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_L4  = 4
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_S   = 4
INTEGER,PARAMETER :: BYTES_PER_ELEMENT_VS  = 4
!!
!! * overhead for containers
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_A0 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_A1 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_A2 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_A3 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_A4 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_A5 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_A6 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_A7 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_P0 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_P1 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_P2 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_P3 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_P4 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_P5 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_P6 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_P7 = 0
INTEGER,PARAMETER :: BYTES_FOR_OVERHEAD_VS = BYTES_PER_ELEMENT_I

END MODULE

max_stars_repo_name: wawiesel/tapack

max_issues_repo_name: wawiesel/tapack

max_forks_repo_name: wawiesel/tapack


text:
      SUBROUTINE VER4_20
      CONTINUE
      RETURN
      END

ext: for


text:
C    @(#)loadarcv.f	20.3 2/13/96
C****************************************************************
C
C   File: loadarcv.f
C   Purpose: IPF subroutine to archive all commands which effect the
C            the base data
C
C   Author: Walt Powell  Date: 24 November 1992
C   Called by:
C
C****************************************************************
C
	subroutine loadarcv
        return
        end

c### arcvfile ###      disable archive logging

c        include 'ipfinc:parametr.inc'
c
c        include 'ipfinc:arcvfile.inc'
c        include 'ipfinc:lfiles.inc'
c
c        character text * 120
c        logical eof
c        integer findex
c
c        eof = .false.
c        rewind inp
c        do while (.not. eof)
c           read (inp, 100, end = 120) text
c  100      format (a)
c           if (text(1:1) .eq. '(' .and.
c     &         findex (text(2:), 'end') .ne. 0) then
c              eof = .true.
c           else
c              last = lastchln (text)
c              if (last .gt. 0) write (arcvfile, 100) text(1:last)
c           endif
c        enddo
c  120   rewind inp
c        return
c        end

max_stars_repo_name: mbheinen/bpa-ipf-tsp

max_issues_repo_name: cuihantao/bpa-ipf-tsp

max_forks_repo_name: mbheinen/bpa-ipf-tsp


text:
       subroutine pebi (x,y,cpb,vol)
!***********************************************************************
!  Copyright, 1993, 2004,  The  Regents of the University of California.
!  This program was prepared by the Regents of the University of
!  California at Los Alamos National Laboratory (the University) under
!  contract No. W-7405-ENG-36 with the U.S. Department of Energy (DOE).
!  All rights in the program are reserved by the DOE and the University.
!  Permission is granted to the public to copy and use this software
!  without charge, provided that this Notice and any statement of
!  authorship are reproduced on all copies. Neither the U.S. Government
!  nor the University makes any warranty, express or implied, or
!  assumes any liability or responsibility for the use of this software.
C***********************************************************************
CD1
CD1  PURPOSE
CD1
CD1  To calculate internodal areas using perpendicular bisectors.
CD1
C***********************************************************************
CD2
CD2  REVISION HISTORY
CD2
CD2 $Log:   /pvcs.config/fehm90/src/pebi.f_a  $
!D2
!D2    Rev 2.5   06 Jan 2004 10:43:36   pvcs
!D2 FEHM Version 2.21, STN 10086-2.21-00, Qualified October 2003
!D2
!D2    Rev 2.4   29 Jan 2003 09:11:48   pvcs
!D2 FEHM Version 2.20, STN 10086-2.20-00
!D2
!D2    Rev 2.3   14 Nov 2001 13:11:28   pvcs
!D2 FEHM Version 2.12, STN 10086-2.12-00
!D2
!D2    Rev 2.2   06 Jun 2001 13:25:44   pvcs
!D2 FEHM Version 2.11, STN 10086-2.11-00
!D2
!D2    Rev 2.1   30 Nov 2000 12:05:48   pvcs
!D2 FEHM Version 2.10, STN 10086-2.10-00
!D2
!D2    Rev 2.0   Fri May 07 14:44:12 1999   pvcs
!D2 FEHM Version 2.0, SC-194 (Fortran 90)
CD2
CD2    Rev 1.3   Thu Jan 18 10:06:00 1996   zvd
CD2 Modifications to prolog
CD2
CD2    Rev 1.2   Wed Jan 10 15:20:50 1996   hend
CD2 Added Prolog
CD2
CD2    Rev 1.1   03/18/94 15:57:50   gaz
CD2 Added solve_new and cleaned up memory management.
CD2
CD2    Rev 1.0   01/20/94 10:26:10   pvcs
CD2 original version in process of being certified
CD2
C***********************************************************************
CD3
CD3  REQUIREMENTS TRACEABILITY
CD3
CD3  2.2 Finite-Element Coefficient Generation
CD3
C***********************************************************************
CD4
CD4  SPECIAL COMMENTS AND REFERENCES
CD4
CD4 Requirements from SDN: 10086-RD-2.20-00
CD4   SOFTWARE REQUIREMENTS DOCUMENT (RD) for the
CD4   FEHM Application Version 2.20
CD4
C***********************************************************************
CD5
CD5 INTERFACES
CD5
CD5 Formal Calling Parameters
CD5
CD5   Identifier      Type     Use  Description
CD5
CD5   cpb             REAL*8   O    Perpindicular bisector
CD5   vol             REAL*8   O    Internodal area
CD5   x               REAL*8   I    X nodal coordinates
CD5   y               REAL*8   I    Y nodal coordinates
CD5
CD5 Interface Tables
CD5
CD5   None
CD5
CD5 Files
CD5
CD5   None
CD5
C***********************************************************************
CD6
CD6 GLOBAL OBJECTS
CD6
CD6 Global Constants
CD6
CD6   None
CD6
CD6 Global Types
CD6
CD6   None
CD6
CD6 Global Variables
CD6
CD6   None
CD6
CD6 Global Subprograms
CD6
CD6   None
CD6
C***********************************************************************
CD7
CD7 LOCAL IDENTIFIERS
CD7
CD7 Local Constants
CD7
CD7   None
CD7
CD7 Local Types
CD7
CD7   None
CD7
CD7 Local variables
CD7
CD7   Identifier      Type     Description
CD7
CD7   arx             REAL     Bisector calculation temporary variable
CD7   ary             REAL     Bisector calculation temporary variable
CD7   dsq             REAL     Bisector calculation temporary variable
CD7   i               INT      Loop index
CD7   i1              INT      Variable index
CD7   i2              INT      Variable index
CD7   i3              INT      Variable index
CD7   pneg            REAL     Area calculation temporary variable
CD7   pos             REAL     Area calculation temporary variable
CD7   xc              REAL     Bisector X coordinate
CD7   yc              REAL     Bisector Y coordinate
CD7
CD7 Local Subprograms
CD7
CD7   None
CD7
C***********************************************************************
CPS
CPS PSEUDOCODE
CPS
CPS BEGIN pebi
CPS
CPS   FOR each coordinate position
CPS
CPS       set indices
CPS
CPS       IF first and second ycoordinate are equal
CPS          compute bisector xcoordinate using first and second y
CPS       ELSE IF first and third ycoordinate are equal
CPS          compute bisector xcoordinate using first and third y
CPS       ELSE
CPS          compute bisector xcoordinate using all y's
CPS       END IF
CPS
CPS       IF first and second xcoordinate are equal
CPS          compute bisector ycoordinate using first and second x
CPS       ELSE IF first and third xcoordinate are equal
CPS          compute bisector ycoordinate using first and third x
CPS       ELSE
CPS          compute bisector ycoordinate using all x's
CPS       END IF
CPS
CPS       calculate internodal area using bisectors
CPS
CPS   END FOR
CPS
CPS END pebi
CPS
C***********************************************************************

      implicit none

      integer i, i1, i2, i3
      real*8 cpb(*), vol(*), x(*), y(*)
      real*8 arx, ary, dsq, pneg, pos, xc, yc

      do i = 0, 2
         i1 = i + 1
         i2 = mod(i+1,3) + 1
         i3 = mod(i+2,3) + 1
         if (y(i1).eq.y(i2)) then
	    xc = 0.5 * (x(i1)+x(i2))
         else if (y(i1).eq.y(i3)) then
	    xc = 0.5 * (x(i1)+x(i3))
         else
	    xc = 0.5 * ( (x(i3)*x(i3)-x(i1)*x(i1))/(y(i3)-y(i1)) -
     &           (x(i2)*x(i2)-x(i1)*x(i1))/(y(i2)-y(i1))+y(i3)-y(i2))/
     &           ( (x(i3)-x(i1))/(y(i3)-y(i1)) -
     &           (x(i2)-x(i1))/(y(i2)-y(i1)) )
         end if
         if (x(i1).eq.x(i2)) then
	    yc = 0.5 * (y(i1)+y(i2))
         else if (x(i1).eq.x(i3)) then
	    yc = 0.5 * (y(i1)+y(i3))
         else
	    yc = 0.5 * ( (y(i3)*y(i3)-y(i1)*y(i1))/(x(i3)-x(i1)) -
     &           (y(i2)*y(i2)-y(i1)*y(i1))/(x(i2)-x(i1))+x(i3)-x(i2))/
     &           ( (y(i3)-y(i1))/(x(i3)-x(i1)) -
     &           (y(i2)-y(i1))/(x(i2)-x(i1)) )
         end if

         arx = yc - 0.5*(y(i1)+y(i2))
         ary = xc - 0.5*(x(i1)+x(i2))
         dsq = (x(i2)-x(i1))**2 + (y(i2)-y(i1))**2
         cpb(i1) = sqrt( (arx*arx+ary*ary)/dsq )

         pos = x(i1)*0.5*(y(i1)+y(i2)) + 0.5*(x(i1)+x(i2))*yc +
     &        xc*0.5*(y(i1)+y(i3)) + 0.5*(x(i1)+x(i3))*y(i1)
         pneg = y(i1)*0.5*(x(i1)+x(i2)) + 0.5*(y(i1)+y(i2))*xc +
     &        yc*0.5*(x(i1)+x(i3)) + 0.5*(y(i1)+y(i3))*x(i1)
         vol(i1) = 0.5 * ( pos - pneg )
      enddo

      return
      end

max_stars_repo_name: satkarra/FEHM

max_issues_repo_name: satkarra/FEHM

max_forks_repo_name: satkarra/FEHM


module SourceM
  use UtilitiesM
  use FortranParser

  implicit none

  private
  public :: SourceDT, source

  type :: SourceDT
     integer(ikind)                                  :: nDim
     type(EquationParser), dimension(:), allocatable :: func
   contains
     procedure, public :: init
     procedure, public :: evaluateFirstFunc
     procedure, public :: evaluateAnyFunc
     generic           :: evaluate => evaluateFirstFunc, evaluateAnyFunc
  end type SourceDT

  interface source
     procedure :: constructor
  end interface source

contains

  type(SourceDT) function constructor(nVar, nDim, var, func)
    implicit none
    integer(ikind)                   , intent(in) :: nVar
    integer(ikind)                   , intent(in) :: nDim
    character(len=*), dimension(nVar), intent(in) :: var
    character(len=*), dimension(nDim), intent(in) :: func
    call constructor%init(nVar, nDim, var, func)
  end function constructor

  subroutine init(this, nVar, nDim, var, func)
    implicit none
    class(SourceDT)                  , intent(inout) :: this
    integer(ikind)                   , intent(in)    :: nVar
    integer(ikind)                   , intent(in)    :: nDim
    character(len=*), dimension(nVar), intent(in)    :: var
    character(len=*), dimension(nDim), intent(in)    :: func
    integer(ikind)                                   :: i
    allocate(this%func(nDim))
    do i = 1, nDim
       this%func(i) = EquationParser(func(i), var)
    end do
  end subroutine init

  real(rkind) function evaluateFirstFunc(this, val)
    implicit none
    class(SourceDT)              , intent(inout) :: this
    real(rkind)    , dimension(:), intent(in)    :: val
    evaluateFirstFunc = this%func(1)%evaluate(val)
  end function evaluateFirstFunc

  real(rkind) function evaluateAnyFunc(this, iFunc, val)
    implicit none
    class(SourceDT)              , intent(inout) :: this
    integer(ikind)               , intent(in)    :: iFunc
    real(rkind)    , dimension(:), intent(in)    :: val
    evaluateAnyFunc = this%func(iFunc)%evaluate(val)
  end function evaluateAnyFunc

end module SourceM


text:
Being here for 3 years, I really wish I had gotten involved in the Wiki a little earlier...

1 more Davis Timeline year to go.

Heh, you would have to bring a gun onto campus in order to store it at the station and if you get caught with it before you get to the station... lulz would ensue. Users/hankim

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


text: I grew up in Southern California and went to school at UC Davis. I have since moved to http://chicowiki.org Chico. I dont really know what else to put here, I have more on my profile there if you really care. This is my 2nd name on Daviswiki, I cant for the life of me remember my first.


text:
subroutine tfi2d(e0, e1, e2, e3, Nu, Nv, X)

  !***DESCRIPTION
  !
  !     Written by Gaetan Kenway
  !
  !     Abstract: Perform a simple 2D transfinite interpolation in 3
  !               spatial dimensions. This function is not used
  !               directly in pySpline, but from pyGeo.
  !
  !     Description of Arguments
  !     Input
  !     e0      - Real, Vector or coordinates along 0th edge. Size (3, Nu)
  !     e1      - Real, Vector or coordinates along 1st edge. Size (3, Nu)
  !     e2      - Real, Vector or coordinates along 2nd edge. Size (3, Nv)
  !     e3      - Real, Vector or coordinates along 3rd edge. Size (3, Nv)
  !
  !     Ouput
  !     X       - Real, Evaluated points, size 3 x Nu x Nv
  !

  use precision
  implicit none

  integer, intent(in) :: Nu,Nv
  real(kind=realType), intent(in) :: e0(3,Nu)
  real(kind=realType), intent(in) :: e1(3,Nu)

  real(kind=realType), intent(in) :: e2(3,Nv)
  real(kind=realType), intent(in) :: e3(3,Nv)

  real(kind=realType), intent(out) :: X(3,Nv,Nu)

  real(kind=realType) :: U(Nu),V(Nv)
  integer :: i,j

  do i=1,Nu
     U(i) = real(i-1)/(Nu-1)
  end do

  do j=1,Nv
     V(j) = real(j-1)/(Nv-1)
  end do

  do i=1,Nu
     do j=1,Nv
        X(:,j,i) = (1-V(j))*e0(:,i) + V(j)*e1(:,i) + &
                   (1-U(i))*e2(:,j) + U(i)*e3(:,j) - ( &
                   (  U(i))*(  V(j))*e1(:,Nu) + &
                   (  U(i))*(1-V(j))*e0(:,Nu) + &
                   (1-U(i))*(  V(j))*e1(:,1)  + &
                   (1-U(i))*(1-V(j))*e0(:,1) &
                   )
     end do
  end do

end subroutine tfi2d

max_stars_repo_name: eirikurj/pyspline

max_issues_repo_name: eirikurj/pyspline

max_forks_repo_name: eirikurj/pyspline


text:
! { dg-do compile }
! { dg-options "-std=f2003" }

! PR fortran/38936
! Check for errors with ASSOCIATE during parsing.

PROGRAM main
  IMPLICIT NONE

  ASSOCIATE ! { dg-error "Expected association list" }

  ASSOCIATE () ! { dg-error "Expected association" }

  ASSOCIATE (a => 1) 5 ! { dg-error "Junk after ASSOCIATE" }

  ASSOCIATE (x =>) ! { dg-error "Expected association" }

  ASSOCIATE (=> 5) ! { dg-error "Expected association" }

  ASSOCIATE (x => 5, ) ! { dg-error "Expected association" }

  myname: ASSOCIATE (a => 1)
  END ASSOCIATE ! { dg-error "Expected block name of 'myname'" }

  ASSOCIATE (b => 2)
  END ASSOCIATE myname ! { dg-error "Syntax error in END ASSOCIATE" }

  myname2: ASSOCIATE (c => 3)
  END ASSOCIATE myname3 ! { dg-error "Expected label 'myname2'" }

  ASSOCIATE (a => 1, b => 2, a => 3) ! { dg-error "Duplicate name 'a'" }

  ASSOCIATE (a => 5)
    INTEGER :: b ! { dg-error "Unexpected data declaration statement" }
  END ASSOCIATE
END PROGRAM main ! { dg-error "Expecting END ASSOCIATE" }
! { dg-excess-errors "Unexpected end of file" }


text:
C++
*    ssssssssssssssssssssssssssssssssssssssssssssssssssssssssssssss
*    s                                                            s
*    s                       UPPER_CASE.FOR                       s
*    s                                                            s
*    ssssssssssssssssssssssssssssssssssssssssssssssssssssssssssssss
*
*    Purpose :   Translate the program from lowercase to uppercase
*                (for the Fortran-77 compatibility)
*
*         in :   infile.for  ( the orignal program)
*        out :   outfile.for (= 'P'+INFILE)
*
*     Author :   Yanling FANG   ESO -- Garching
*    Language:   Fortran 77
*    Version :   1.1 10-Nov-1988
*
*       note :   Only translate the command part, not including the
*                comment and character to be written
*
C--

*
	IMPLICIT NONE
	INTEGER LENGTH,IPOSI

	CHARACTER*40 INFILE
	CHARACTER*80 LINE,SCHAR

*    Define read & write units, and file's channel to be opened
	INTEGER  IREAD,IWRIT,IIN,IOUT
        DATA     IREAD,IWRIT,IIN,IOUT,LINE / 5,6,10,20,' ' /

*    Open input& output  files
	WRITE(IWRIT,'(A)')

100	CONTINUE
	WRITE(IWRIT,'(A,$)') ' Name of the input file : '
	READ(IREAD,'(A)')INFILE

	OPEN ( UNIT=IIN,FILE=INFILE,ERR=999,STATUS='OLD')
	GOTO 150
999	WRITE(IWRIT,'(A)')' Wrong name, have another try !'
	GOTO 100

150	OPEN ( UNIT=IOUT,FILE='P'//INFILE,
     :          STATUS='NEW')

	WRITE(IWRIT,'(A)')' Working now ... '

*    Read the input file and transfer it
200	CONTINUE

	READ(IIN,400,ERR=350,END=300) LINE

        IF (.NOT.(LINE(1:1).EQ.'*'.OR.
     :            LINE(1:1).EQ.'C'.OR.
     :            LINE(1:1).EQ.'c'.OR.
     :            LINE(1:1).EQ.'!') ) CALL UP(LINE)

*    To remove the dummy blanks in the end of line
         LENGTH   = LEN(LINE)
         DO 250 I = LENGTH,1,-1
           SCHAR = LINE(I:I)
           IPOSI = I
           IF (SCHAR.NE.' ') GOTO 260
250      CONTINUE

*    To write the proper line in the output file
260     CONTINUE
	WRITE(IOUT,'(A)') LINE(:IPOSI)

	GOTO 200

300     CONTINUE
*    end transfer

        GOTO 360

350     WRITE(IWRIT,'(A)') ' Error while reading the program'

360     CONTINUE

400	FORMAT(A)

	WRITE(IWRIT,'(A)')

	CLOSE(UNIT=IIN)
	CLOSE(UNIT=IOUT)

	END

*----------------------------------------------------------------

      SUBROUTINE UP(C)

*+
*  purpose: change the lowercase into uppercase
*  note: 1. We don't change that the characters inside '...', this
*           is done by ICOUNT;
*        2. IRECORD is considered for the case '' when writing the
*           like :  WRITE(IWRIT,'('' simple way of writing '')')
*+

      CHARACTER*(*) C
      INTEGER J,ICOUNT,IRECORD

      IRECORD = 0
      ICOUNT = 1

      DO 40 J = 1,LEN(C)

         IF (C(J:J).EQ.'''') THEN
*        when  '  comes one by another, ICOUNT is not incremented
            IF ( J.NE.(IRECORD-1) ) ICOUNT = ICOUNT + 1
            IRECORD = J
         ENDIF

         GOTO (10,40,30), ICOUNT

*     Normal letter (not " ' ")
10       CONTINUE
         IF(C(J:J).GE.'a'.AND.C(J:J).LE.'z') THEN
             C(J:J) = CHAR(ICHAR(C(J:J)) + 'A' - 'a')
         ENDIF
         GOTO 40

*     Meet ' second time
30       CONTINUE
         ICOUNT = 1

40    CONTINUE

      END


!Definicion de la funcion a Evaluar
Real *8 Function f (x) !Declaracion
  Implicit None
  Real*8 :: x !Variables que se obtienen
  f =exp(-x*x/2) !funcion
end Function f
!
!Empieza Programa MAIN
Program ActI

  !Definicion de Variables
  Implicit None
  Real *8 :: f               !Funcion
  Real *8 :: x0              !Punto de interes
  Real *8 :: xi , xf, y, z   !Rango de valoración
  Real *4 :: k = 0.1         !Incrementos

  !Pedir Info a User
  Write(*,*) "==================================================="
  Write(*,*) "La función a Evaluar es exp(-x**2/2)"
  Write(*,*) "==================================================="
  Write(*,*) "Punto de Interes x0 a Evaluar"
  Read(*,*) x0

  !Evaluando Funcion en x0 e imprimeiendo en pantalla
  Write(*,*) "exp(-x**2/2) evaluado en", real(x0), "es" , f(x0)
 ! Write(*,*) "==================================================="

  !Calculo en un rango dado
  Write(*,*) "Rango de Medicion, Limite Inferior-Limite Superior"
  Read(*,*)  xi , xf

  if (xi > xf) then !Checar el orden de los limites
     Write(*,*) "Escribiste los limites en el orden incorrecto. Corrigiendo"
     z=xi
     xi=xf       !Reacomodo de los Limites
     xf=z
  end if

  !Abriendo documento para escribir salida
  Open(1, file="Dat.dat")

  do while (xi <= xf+k)
     y = f(xi)
     Write(1,*) xi, y
     xi =xi + k
  end do

  Close(1)
  Write(*,*) "Fin del Calculo"

End Program ActI


text:
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++!
!                          Futility Development Group                          !
!                             All rights reserved.                             !
!                                                                              !
! Futility is a jointly-maintained, open-source project between the University !
! of Michigan and Oak Ridge National Laboratory.  The copyright and license    !
! can be found in LICENSE.txt in the head directory of this repository.        !
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++!
!> @brief The global geometry module, collecting all public members of
!> other geometry modules. This is the module that should be used elsewhere
!> in the code.
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++!
MODULE Geom
USE IntrType
USE Constants_Conversion
USE Geom_Graph
USE Geom_Points
USE Geom_Line
USE Geom_Plane
USE Geom_CircCyl
USE Geom_Box
USE Geom_Poly
USE ParameterLists

IMPLICIT NONE
PRIVATE

PUBLIC :: newGeom

PUBLIC :: GraphType
PUBLIC :: MAX_COORD_STR_LEN
PUBLIC :: PointType
PUBLIC :: LinkedListPointType
PUBLIC :: LineType
PUBLIC :: PlaneType
PUBLIC :: CircleType
PUBLIC :: CylinderType
PUBLIC :: OBBoxType
PUBLIC :: ABBoxType
PUBLIC :: PolygonType
PUBLIC :: Distance
PUBLIC :: midPoint
PUBLIC :: ClearLinkedListPointType
PUBLIC :: Polygonize
PUBLIC :: OPERATOR(+)
PUBLIC :: OPERATOR(-)
PUBLIC :: OPERATOR(==)
PUBLIC :: OPERATOR(/=)
PUBLIC :: OPERATOR(.APPROXEQA.)
PUBLIC :: ASSIGNMENT(=)

INTERFACE newGeom
  !> @copybrief Geom::newGeom_line
  !> @copydetails Geom::newGeom_line
  MODULE PROCEDURE newGeom_line
  !> @copybrief Geom::newGeom_plane
  !> @copydetails Geom::newGeom_plane
  MODULE PROCEDURE newGeom_plane
  !> @copybrief Geom::newGeom_box
  !> @copydetails Geom::newGeom_box
  MODULE PROCEDURE newGeom_box
  !> @copybrief Geom::newGeom_circ
  !> @copydetails Geom::newGeom_circ
  MODULE PROCEDURE newGeom_circ
  !> @copybrief Geom::newGeom_cyl
  !> @copydetails Geom::newGeom_cyl
  MODULE PROCEDURE newGeom_cyl
  !> @copybrief Geom::newGeom_poly
  !> @copydetails Geom::newGeom_poly
  MODULE PROCEDURE newGeom_poly
ENDINTERFACE
!
!===============================================================================
CONTAINS
!
!-------------------------------------------------------------------------------
!> @brief
!> @param params
!> @param geom
!>
SUBROUTINE newGeom_line(params,geom)
  TYPE(ParamType),INTENT(IN) :: params
  TYPE(LineType),INTENT(INOUT) :: geom

  REAL(SRK),ALLOCATABLE :: c0(:),c1(:)
  TYPE(PointType) :: p0,p1

  CALL params%get('LineGeom -> StartPoint',c0)
  CALL params%get('LineGeom -> EndPoint',c1)
  IF(ALLOCATED(c0)) CALL p0%init(COORD=c0)
  IF(ALLOCATED(c1)) CALL p1%init(COORD=c1)
  CALL geom%set(p0,p1)
ENDSUBROUTINE newGeom_line
!
!-------------------------------------------------------------------------------
!> @brief
!> @param params
!> @param geom
!>
SUBROUTINE newGeom_plane(params,geom)
  TYPE(ParamType),INTENT(IN) :: params
  TYPE(PlaneType),INTENT(INOUT) :: geom

  REAL(SRK),ALLOCATABLE :: c0(:),v0(:)
  TYPE(PointType) :: p0

  CALL params%get('PlaneGeom -> Point',c0)
  CALL params%get('PlaneGeom -> NormalVector',v0)
  IF(ALLOCATED(c0)) CALL p0%init(COORD=c0)
  IF(ALLOCATED(v0)) CALL geom%set(v0,p0)
ENDSUBROUTINE newGeom_plane
!
!-------------------------------------------------------------------------------
!> @brief
!> @param params
!> @param geom
!>
SUBROUTINE newGeom_box(params,geom)
  TYPE(ParamType),INTENT(IN) :: params
  TYPE(OBBoxType),INTENT(INOUT) :: geom

  REAL(SRK),ALLOCATABLE :: c0(:),v1(:),v2(:),v3(:),extent(:)
  TYPE(PointType) :: p0

  CALL params%get('BoxGeom -> CornerPoint',c0)
  CALL params%get('BoxGeom -> Vector1',v1)
  CALL params%get('BoxGeom -> Vector2',v2)
  CALL params%get('BoxGeom -> ExtentVector',extent)
  IF(params%has('BoxGeom -> Vector3')) &
      CALL params%get('BoxGeom -> Vector3',v3)
  IF(ALLOCATED(c0)) CALL p0%init(COORD=c0)
  IF(ALLOCATED(v3)) THEN
    CALL geom%set(p0,extent,v1,v2,v3)
  ELSEIF(ALLOCATED(v1) .AND. ALLOCATED(v2)) THEN
    CALL geom%set(p0,extent,v1,v2)
  ENDIF
ENDSUBROUTINE newGeom_box
!
!-------------------------------------------------------------------------------
!> @brief
!> @param params
!> @param geom
!>
SUBROUTINE newGeom_circ(params,geom)
  TYPE(ParamType),INTENT(IN) :: params
  TYPE(CircleType),INTENT(INOUT) :: geom

  REAL(SRK) :: r,theta1,theta2
  REAL(SRK),ALLOCATABLE :: c0(:)
  TYPE(PointType) :: p0

  theta1=0.0_SRK
  theta2=TWOPI
  CALL params%get('CircleGeom -> Radius',r)
  CALL params%get('CircleGeom -> Centroid',c0)
  IF(params%has('CircleGeom -> StartingAngle')) &
      CALL params%get('CircleGeom -> StartingAngle',theta1)
  IF(params%has('CircleGeom -> StoppingAngle')) &
      CALL params%get('CircleGeom -> StoppingAngle',theta2)

  IF(ALLOCATED(c0)) CALL p0%init(COORD=c0)
  CALL geom%set(p0,r,theta1,theta2)
ENDSUBROUTINE newGeom_circ
!
!-------------------------------------------------------------------------------
!> @brief
!> @param params
!> @param geom
!>
SUBROUTINE newGeom_cyl(params,geom)
  TYPE(ParamType),INTENT(IN) :: params
  TYPE(CylinderType),INTENT(INOUT) :: geom

  REAL(SRK) :: r,theta1,theta2
  REAL(SRK),ALLOCATABLE :: c0(:),c1(:)
  TYPE(PointType) :: p0,p1

  theta1=0.0_SRK
  theta2=TWOPI
  CALL params%get('CylinderGeom -> Radius',r)
  CALL params%get('CylinderGeom -> BottomCentroid',c0)
  CALL params%get('CylinderGeom -> TopCentroid',c1)
  IF(params%has('CylinderGeom -> StartingAngle')) &
      CALL params%get('CylinderGeom -> StartingAngle',theta1)
  IF(params%has('CylinderGeom -> StoppingAngle')) &
      CALL params%get('CylinderGeom -> StoppingAngle',theta2)

  IF(ALLOCATED(c0)) CALL p0%init(COORD=c0)
  IF(ALLOCATED(c1)) CALL p1%init(COORD=c1)
  CALL geom%set(p0,p1,r,theta1,theta2)
ENDSUBROUTINE newGeom_cyl
!
!-------------------------------------------------------------------------------
!> @brief
!> @param params
!> @param geom
!>
SUBROUTINE newGeom_poly(params,geom)
  TYPE(ParamType),INTENT(IN) :: params
  TYPE(PolygonType),INTENT(INOUT) :: geom
  CHARACTER(LEN=8) :: ivchar
  INTEGER(SIK) :: i,n,m
  REAL(SRK),ALLOCATABLE :: vert(:),vPrev(:)
  TYPE(GraphType) :: tmpG

  n=0; m=0
  CALL params%get('PolygonGeom->nVert',n)
  CALL params%get('PolygonGeom->nQuad',m)
  IF(n > 2) THEN
    CALL params%get('PolygonGeom->vertex 1',vert)
    CALL tmpG%insertVertex(vert)
    DO i=2,n
      CALL MOVE_ALLOC(vert,vPrev)
      WRITE(ivchar,'(i8)') i; ivchar=ADJUSTL(ivchar)
      CALL params%get('PolygonGeom->vertex '//TRIM(ivchar),vert)
      CALL tmpG%insertVertex(vert)
      CALL tmpG%defineEdge(vPrev,vert)
    ENDDO
    CALL params%get('PolygonGeom->vertex 1',vPrev)
    CALL tmpG%defineEdge(vPrev,vert)
    DO i=1,m
      WRITE(ivchar,'(i8)') i; ivchar=ADJUSTL(ivchar)
      !CALL params%get('PolygonGeom->quad edge '//TRIM(ivchar))
      !Set the quadratic edge data
    ENDDO
    CALL geom%set(tmpG)
    CALL tmpG%clear()
  ENDIF
ENDSUBROUTINE newGeom_poly
!
ENDMODULE Geom

max_stars_repo_name: picmc/Futility

max_issues_repo_name: picmc/Futility

max_forks_repo_name: wgurecky/Futility


text:
	SUBROUTINE GDPGTS  ( iflno, time, ivcord, rgx, rgy,
     +			     vclsfc, havsfc, parm, havgfs, iret )
C************************************************************************
C* GDPGTS								*
C*									*
C* This subroutine gets the surface data for a cross section.		*
C*									*
C* GDPGTS  ( IFLNO, TIME, IVCORD, RGX, RGY,  VCLSFC, HAVSFC, 		*
C*	     PARM, HAVGFS, IRET )					*
C*									*
C* Input parameters:							*
C*	IFLNO		  INTEGER	Grid file number		*
C*	TIME  (2)	  CHAR*		Time to search for levels	*
C*	IVCORD		  INTEGER	Vertical coordinate for search	*
C*	RGX  		  REAL		X grid coordinate		*
C*	RGY  		  REAL		Y grid coordinate		*
C*									*
C* Output parameters:							*
C*      VCLSFC            REAL 		Vert coord location of sfc	*
C*      HAVSFC            LOGICAL       Flag for existence of sfc	*
C*	PARM		  CHAR*		Parameter name			*
C*      HAVGFS            LOGICAL       Flag for existence of sfc data	*
C*	IRET		  INTEGER	Return code			*
C*					  0 = normal return		*
C*					 -6 = GVCORD is invalid		*
C*					 +2 = no sfc value found	*
C**									*
C* Log:									*
C* D. McCann/NSSFC	12/94		Created from gdxgts.f		*
C* R. Tian/SAIC		 5/06		Changed to call new DG subs	*
C************************************************************************
	INCLUDE		'DGCMN.CMN'
C*
	CHARACTER*(*)	time ( 2 ), parm
        REAL		rgx , rgy , vclsfc
	LOGICAL		havsfc, havgfs, ermiss
C*
	REAL		grid ( LLMXGD )
	INTEGER		level ( 2 ), igrhdr ( LLGDHD )
C------------------------------------------------------------------------
	iret = 0
	havsfc = .false.
	havgfs = .false.
       ermiss = .false.
C
C*	Set the level and parameter for which to search.
C
	level (1) =  0
	level (2) = -1
	CALL LV_CCRD  ( ivcord, parm, ier )
C
C*	Try to read the surface data on IVCORD.
C
	CALL DG_NRDT  ( iflno, time, level, ivcord, parm,
     +			grid, igx, igy, igrhdr, ier )
	IF  ( ier .eq. 0 )  THEN
	    havsfc = .true.
	    havgfs = .true.
	ELSE
C
C*	    Try to read the surface data on JVCORD = 0.
C
	    jvcord = 0
	    CALL DG_NRDT  ( iflno, time, level, jvcord, parm,
     +			    grid, igx, igy, igrhdr, ier )
	    IF  ( ier .eq. 0 )  THEN
		havsfc = .true.
		ELSE
C
C*	    	Try to read the surface data on JVCORD = 4.
C
	    	jvcord = 4
	    	level (1) = 10000
	    	CALL DG_NRDT  ( iflno, time, level, jvcord, parm,
     +			    grid, igx, igy, igrhdr, ier )
		IF  ( ier .eq. 0 )  THEN
			havsfc = .true.
	    	END IF
	    END IF
	END IF
C
	if (havsfc) THEN
	    CALL GR_INTP ( 1, rgx, rgy, 1, igx, igy, grid, vclsfc, ierr )
	    IF  ( ierr .ne. 0 ) THEN
		vclsfc = RMISSD
		havsfc = .false.
	    END IF
	ELSE
            vclsfc = RMISSD
	END IF
C*
	RETURN
	END

max_stars_repo_name: oxelson/gempak

max_issues_repo_name: oxelson/gempak

max_forks_repo_name: oxelson/gempak


text:
C$Procedure      SETELM ( Calculate set of SPICE conic elements )

      SUBROUTINE SETELM ( CNGMVL, DVL, DSTCOF, ANGCOF,
     .                    PARAM,  NPARAM, ELTS, ERRFLG, ERRTXT )

C$ Abstract
C
C     Calculate set of standard conic elements used by SPICELIB
C     procedure CONIC from a set of different parameters allowed for
C     the input data type ELEMENTS.
C
C$ Disclaimer
C
C     THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE
C     CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S.
C     GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE
C     ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE
C     PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS"
C     TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY
C     WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A
C     PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC
C     SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE
C     SOFTWARE AND RELATED MATERIALS, HOWEVER USED.
C
C     IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA
C     BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT
C     LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND,
C     INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS,
C     REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE
C     REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY.
C
C     RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF
C     THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY
C     CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE
C     ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE.
C
C$ Required_Reading
C
C     MKSPK User's Guide
C
C$ Keywords
C
C     None.
C
C$ Declarations

      IMPLICIT   NONE

      DOUBLE PRECISION      CNGMVL
      DOUBLE PRECISION      DVL    ( * )
      DOUBLE PRECISION      ANGCOF
      DOUBLE PRECISION      DSTCOF
      INTEGER               PARAM  ( * )
      INTEGER               NPARAM
      DOUBLE PRECISION      ELTS   ( 8 )
      LOGICAL               ERRFLG
      CHARACTER*(*)         ERRTXT ( 2 )

C$ Brief_I/O
C
C     Variable  I/O  Description
C     --------  ---  --------------------------------------------------
C     CNGMVL     I   Gravitational parameter (GM) of primary body
C     DVL        I   Array of input parameters values
C     PARAM      I   Array of parameter indexes
C     NPARAM     I   Number elements in PARAM
C     ANGCOF     I   Angle unit coefficient
C     DSTCOF     I   Distance unit coefficient
C     ELTS       O   Set of conic elements
C     ERRFLG     O   Logical error flag
C     ERRTXT     O   Error explanation text
C
C$ Detailed_Input
C
C                                                       3    2
C     CNGMVL      is the gravitational parameter (GM, km /sec ) of
C                 the primary body.
C
C     DVL         is a set of parameter values obtained by parsing
C                 input text data of type ELEMENTS according
C                 to the specified DATA_ORDER. The parameters in the
C                 arrays are the following order (see array DOVAL in
C                 the main module):
C
C                    'EPOCH', 'X', 'Y', 'Z', 'VX' , 'VY', 'VZ',
C                    'A', 'E', 'RP', 'T', 'P', 'INC', 'PER', 'NOD',
C                    'MEAN', 'EXAN', 'TRAN', 'EPOCHP', 'TAU',
C                    'EQ_A', 'EQ_H', 'EQ_K', 'EQ_ML', 'EQ_P', 'EQ_Q',
C                    'DPER/DT', 'DMPN/DT', 'DNOD/DT' , 'SKIP'
C
C     PARAM       is an array of input parameter indexes in DVL
C                 arranged in accordance with the specified DATA_ORDER
C                 value. Index of each parameter is as ordinal number
C                 of the corresponding word in array DOVAL.
C
C     NPARAM      is the number of elements in array PARAM.
C
C     ANGCOF      is the coefficient that converts input angle units
C                 to radians.
C
C     DSTCOF      is the coefficient that converts input distance units
C                 to kilometers.
C
C$ Detailed_Output
C
C     ELTS        are conic elements describing the orbit of an object
C                 around a primary body. The elements are, in order:
C
C                            Perifocal distance.
C                            Eccentricity.
C                            Inclination.
C                            Longitude of the ascending node.
C                            Argument of periapsis.
C                            Mean anomaly at epoch.
C                            Epoch.
C                            Gravitational parameter.
C
C                 The epoch of the elements is the epoch of the input
C                 state. Units are km, rad, rad/sec. The same elements
C                 are used to describe all three types (elliptic,
C                 hyperbolic, and parabolic) of conic orbit.
C
C     ERRFLG      is the logical flag of the error. Set to .TRUE. value
C                 if calculation of the comic elements from a given
C                 input parameters wasn't possible.
C
C     ERRTXT      is the error explanation text.
C
C$ Parameters
C
C     None.
C
C$ Exceptions
C
C     If input elements are inconsistent or incomplete, this routine
C     doesn't signal an error but rather returns long anf short
C     error messages through its output parameters.
C
C$ Files
C
C     None.
C
C$ Particulars
C
C     None.
C
C$ Examples
C
C     None.
C
C$ Restrictions
C
C     None.
C
C$ Literature_References
C
C     None.
C
C$ Author_and_Institution
C
C     N.G. Khavenson (IKI RAS, Russia)
C     B.V. Semenov   (NAIF, JPL)
C
C$ Version
C
C-    Version 1.0.4, 18-JAN-2000 (BVS).
C
C        Added CHKOUT calls before each RETURN.
C
C-    Version 1.0.3, 19-MAR-1999 (BVS).
C
C        Corrected comments.
C
C-    Version 1.0.2, 13-JAN-1999 (BVS).
C
C        Modified error messages.
C
C-    Version 1.0.1, 22-NOV-1998 (NGK).
C
C
C-    Version 1.0.0, 8-SEP-1998 (NGK).
C
C-&

C$ Index_Entries
C
C     Calculate set of conic elements from MKSPK input parameters.
C
C-&

C
C     SPICELIB functions
C
      INTEGER               ISRCHI
      DOUBLE PRECISION      TWOPI
      DOUBLE PRECISION      PI
      DOUBLE PRECISION      DACOSH
      LOGICAL               RETURN

C
C     Local variables
C

      DOUBLE PRECISION      TWOTHR
      DOUBLE PRECISION      SENSE
      DOUBLE PRECISION      WRKEL
      DOUBLE PRECISION      MEANM

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

C
C     Initialize error flag.
C
      ERRFLG = .FALSE.

C
C     Assign constant value.
C
      TWOTHR = 2.D0 / 3.D0

C
C     All listed below DVL are present in the input data. We already
C     checked this when we checked specified DATA_ORDER against input
C     output type in CHCKDO.
C
      ELTS (2) = DVL (9)
      ELTS (3) = DVL (13) * ANGCOF
      ELTS (4) = DVL (15) * ANGCOF
      ELTS (5) = DVL (14) * ANGCOF
      ELTS (7) = DVL (1)
      ELTS (8) = CNGMVL

C
C     First we need to compute the distance of the periapsis.
C
      IF ( ISRCHI ( 10, NPARAM, PARAM ) .NE. 0 ) THEN

C
C        Distance of the peripsis is defined in input file.
C
         ELTS (1) = DVL (10) * DSTCOF

      ELSE IF ( ISRCHI ( 8, NPARAM, PARAM ) .NE. 0) THEN

C
C        Semi-major axis is defined in input file.
C
         IF ( DVL (9) .NE. 1.D0  ) THEN

C           Case of elliptic or hyperbolic orbits.
C
            ELTS (1) = DABS ( ( DVL (9) - 1.D0 ) * DVL (8) ) * DSTCOF

         ELSE

C
C           Case of parabolic orbit. Complain and exit.
C
            ERRTXT (1) = 'Semi-major axis can not be used as ' //
     .                   'input parameter in case of a parabolic ' //
     .                   'orbit (ecc = 1.0).'
            ERRTXT (2) = 'SPICE(INCONSISTELEMENTS)'

            ERRFLG = .TRUE.

            CALL CHKOUT ( 'SETELM' )
            RETURN

         END IF

      ELSE IF ( ISRCHI( 12, NPARAM, PARAM ) .NE. 0 ) THEN

C
C        Semi-latus rectum is defined in input file.
C
         ELTS (1) = DVL (12) / ( 1.D0 + DVL (9) ) * DSTCOF

      ELSE IF ( ISRCHI( 11, NPARAM, PARAM ) .NE. 0 ) THEN

C
C        Orbital period is defined in input file.
C
         IF ( DVL (9) .LT. 1.D0 )  THEN

C
C           Case of elliptic orbit.
C
            ELTS (1) = ( 1.D0 - DVL (9) ) *
     .            ( DVL (11) * DSQRT ( CNGMVL ) / TWOPI () ) ** TWOTHR
     .             * DSTCOF

         ELSE

C
C           Case of parabolic or hyperbolic orbit. Complain and exit.
C
            ERRTXT (1) = 'Orbital period can not be used as input ' //
     .                   'parameter in case of a parabolic ' //
     .                   '(ecc = 1.0) or hyperbolic (ecc > 1.0) ' //
     .                   'orbit.'
            ERRTXT (2) = 'SPICE(INCONSISTELEMENTS)'

            ERRFLG = .TRUE.

            CALL CHKOUT ( 'SETELM' )
            RETURN

         END IF

      ELSE

C
C        Not enough input data to calculate the periapsis distance.
C        Complain and exit.
C
         ERRTXT (1) = 'Perifocal distance can not be computed ' //
     .                'from given set of orbital elements.'
         ERRTXT (2) = 'SPICE(INCOMPLETEELEMENTS)'

         ERRFLG = .TRUE.

         CALL CHKOUT ( 'SETELM' )
         RETURN

      END IF

C
C     Second we need to compute the mean anomaly.
C
      IF ( ISRCHI ( 16, NPARAM, PARAM ) .NE. 0 ) THEN

C
C        Mean anomaly is defined in the input file.
C
         ELTS (6) = DVL (16) * ANGCOF

      ELSE IF ( ISRCHI( 17, NPARAM, PARAM ) .NE. 0 ) THEN

C
C        Eccentric anomaly is defined in the input file.
C
         WRKEL = DVL (17) * ANGCOF

         IF ( DVL (9) .LT. 1.D0  ) THEN

C
C           Case of elliptic orbit.
C
            ELTS (6) =  WRKEL - DVL (9) * DSIN ( WRKEL )

         ELSE IF ( DVL (9) .GT. 1.D0  ) THEN

C
C           Case of hyperbolic orbit.
C
            ELTS (6) = DVL (9) * DSINH ( WRKEL ) - WRKEL

         ELSE

C
C           Case of parabolic  orbit. Complain and exit.
C
            ERRTXT (1) = 'Eccentric anomaly can not be used as ' //
     .                   'input parameter in case of a parabolic ' //
     .                   'orbit (ecc = 1.0).'
            ERRTXT (2) = 'SPICE(INCONSISTELEMENTS)'
            ERRFLG = .TRUE.

            CALL CHKOUT ( 'SETELM' )
            RETURN

         END IF

      ELSE IF ( ISRCHI( 18, NPARAM, PARAM ) .NE. 0 ) THEN

C
C        True anomaly is defined in input file.
C
         WRKEL = DVL (18) * ANGCOF

         IF ( WRKEL .LT. PI () ) THEN

            SENSE = 1.D0

         ELSE

            SENSE = -1.D0

         END IF

         IF ( DVL (9) .LT. 1.D0  ) THEN

C
C           Case of elliptic orbit.
C
            WRKEL = DACOS (
     .           (   DVL (9) + DCOS ( WRKEL )      ) /
     .           ( 1.D0 + DVL (9) * DCOS ( WRKEL ) )   )
            ELTS (6) = SENSE * ( WRKEL - DVL (9) * DSIN ( WRKEL ) )

            IF ( ELTS (6) .LT. 0.D0 ) THEN

               ELTS (6) = ELTS (6) + TWOPI ()

            END IF

         ELSE  IF ( DVL (9) .GT. 1.D0  ) THEN

C
C           Case of hyperbolic orbit.
C
            WRKEL = DACOSH (
     .           (   DVL (9) + DCOS ( WRKEL )     ) /
     .           ( 1.D0 + DVL (9) * DCOS( WRKEL ) )    )
            ELTS (6) = SENSE*( ( DVL (9) * DSINH ( WRKEL ) - WRKEL ) )

         ELSE

C
C           Case of parabolic orbit.
C
            WRKEL = DTAN ( WRKEL / 2.D0 )
            ELTS (6) = SENSE * ( WRKEL + WRKEL ** 3 / 3.D0 )

         END IF

      ELSE IF ( ISRCHI( 19, NPARAM, PARAM ) .NE. 0 ) THEN

C
C        Epoch of periapsis is defined in input file.
C
         IF ( DVL (9) .LT. 1.D0 ) THEN

C
C           Case of elliptic orbit.
C
            WRKEL = ( 1.D0 - DVL (9) ) / ELTS (1)
            MEANM = DSQRT ( ELTS (8) * WRKEL ) *  WRKEL
            ELTS (6) = MEANM * ( DVL (1) - DVL (19) )

         ELSE  IF ( DVL (9) .GT. 1.D0  ) THEN

C
C           Case of hyperbolic orbit.
C
            WRKEL = ( DVL (9) - 1.D0 ) / ELTS (1)
            MEANM = DSQRT ( ELTS (8) * WRKEL ) * WRKEL
            ELTS (6) = MEANM * ( DVL (1) - DVL (19) )

         ELSE

C
C           Case of parabolic orbit.
C
            MEANM = DSQRT ( ELTS (8) / 2.D0 / ELTS (1) ) / ELTS (1)
            ELTS (6) = MEANM * ( DVL (1) - DVL (19) )

         END IF

      ELSE IF ( ISRCHI( 20, NPARAM, PARAM ) .NE. 0 ) THEN

C
C        Time interval between current epoch and periapsis epoch
C        is defined in input file.
C
         IF ( DVL (9) .LT. 1.D0 ) THEN

C
C           Case of elliptic orbit.
C
            WRKEL = ( 1.D0 - DVL (9) ) / ELTS (1)
            MEANM = DSQRT ( ELTS (8) * WRKEL ) * WRKEL
            ELTS (6) = MEANM * DVL (20)

         ELSE  IF ( DVL (9) .GT. 1.D0  ) THEN

C
C           Case of hyperbolic orbit.
C
            WRKEL = ( DVL (9) - 1.D0 ) / ELTS (1)
            MEANM = DSQRT ( ELTS (8) * WRKEL ) * WRKEL
            ELTS ( 6 ) = MEANM * DVL (20)

         ELSE

C
C           Case of parabolic orbit.
C
            MEANM = DSQRT ( ELTS (8) / 2.D0 / ELTS (1) ) / ELTS (1)
            ELTS (6) = MEANM * DVL (20)

         END IF

      ELSE

C
C        Not enough input data to calculate mean anomaly.
C        Complain and exit.
C
         ERRTXT (1) = 'Mean anomaly can not be computed ' //
     .                'from given set of orbital elements.'
         ERRTXT (2) = 'SPICE(INCOMPLETEELEMENTS)'
         ERRFLG = .TRUE.

         CALL CHKOUT ( 'SETELM' )
         RETURN

      END IF

      CALL CHKOUT ( 'SETELM' )
      RETURN

      END

max_stars_repo_name: agforero/FTFramework

max_issues_repo_name: agforero/fortran-testing-framework

max_forks_repo_name: agforero/fortran-testing-framework


text:
C
C
*+ mapcat_remove

       subroutine mapcat_remove(imap,check_status,status)
C      --------------------------------------------------
C
C Remove an entry from the map catalogue
C Input:
C    Map entry
       integer      imap
C    Check status of catalogue entry before removing it
       logical      check_status
C Returned:
C    Status
       integer      status
C
C Remove an entry from the map catalogue, the disc file will remain.
C If check_status is set to true then the entry must be a permanent
C map file. If set to false this check is not performed.
*-
       include 'mapcat_pars.inc'

       integer   map_status(10)

       if (status.ne.0) return
       call mapcat_end(imap,status)
       call mapcat_enqst(imap,map_status,status)
       if (check_status) then
         if (map_status(ip_page).eq.true) then
           call io_wrout(
     *     '***(CATALOGUE) Catalogue entry is not a permanent file')
           call io_wrout(
     *     '***(CATALOGUE) Use Export-Permanent-Map or Delete-Map-File')
           return
         end if
       end if
       map_status(ip_access) = access_clear
       map_status(ip_data)   = false
       map_status(ip_open)   = false
       map_status(ip_unit)   = 0
       map_status(ip_page)   = false
       call mapcat_setst(imap,map_status,status)
       call stack_remove(imap,status)
       end

max_stars_repo_name: CavendishAstrophysics/anmap

max_issues_repo_name: CavendishAstrophysics/anmap

max_forks_repo_name: CavendishAstrophysics/anmap


 !***********************************************************************************
    SUBROUTINE upscale(ix,iy,iz)
    !c----------------------------------------------------------------------
    !c
    !c----------------------------------------------------------------------
    use varmod
    implicit none
    integer :: upix,upiy,upiz,Fineix,Fineiy,Fineiz
    INTEGER, INTENT(IN) :: ix,iy,iz

    upix = (ix-1)/ncoarse+1
    upiy = (iy-1)/ncoarse+1
    upiz = (iz-1)/ncoarse+1

    do Fineix = max((upix-1)*ncoarse+1,1),min(upix*ncoarse,nx)
        do Fineiy = max((upiy-1)*ncoarse+1,1),min(upiy*ncoarse,ny)
            do Fineiz = max((upiz-1)*ncoarse+1,1),min(upiz*ncoarse,nz)

                if(simim(Fineix,Fineiy,Fineiz)>0) then
                    simim(ix,iy,iz) = simim(Fineix,Fineiy,Fineiz)
                    DenSim(ix,iy,iz) = DenSim(Fineix,Fineiy,Fineiz)
                    VpSim(ix,iy,iz) = VpSim(Fineix,Fineiy,Fineiz)
                    SnySim(ix,iy,iz) = SnySim(Fineix,Fineiy,Fineiz)
                    go to 56
                end if
            end do
        end do
    end do

56  continue

    END SUBROUTINE upscale

    !***********************************************************************************


text:
*
*     ------------------------------------------------------------------
*       O R T H O G G
*     ------------------------------------------------------------------
*
      SUBROUTINE ORTHOGG(LET)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC(4),IALL,JSC(3),ISCW
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/MEDEFN/IHSH,NJ(16),LJ(16),NOSH1(16),NOSH2(16),J1QN1(31,3),
     :     J1QN2(31,3),IJFUL(16)
      COMMON/SIGNF /SIGNFA
*
*     THIS SUBROUTINE CHECKS FOR POSSIBLE ORTHOGONALITY DUE TO
*     COUPLING DIFFERENCES OR UNEVEN PARITY
*
  102 FORMAT(52H ORTHOGONALITY IN COUPLING SCHEMES OF CONFIGURATIONS)
  103 FORMAT(59H THE TWO CONFIGURATIONS HAVE DIFFERING NUMBERS OF ELECTR
     :ONS)
  104 FORMAT(51H THE TWO CONFIGURATIONS HAVE DIFFERING TOTAL PARITY)
*
* --- DO PSI AND PSIP CONTAIN THE SAME NUMBERS OF ELECTRONS
*     DO PSI AND PSIP HAVE THE SAME TOTAL PARITY
*
      N5=0
      N6=0
      N7=0
      DO 1 I=1,IHSH
        L1=LJ(I)
        L2=NOSH1(I)
        L3=NOSH2(I)
        N5=N5+L2
        N6=N6+L3
        N7=N7+L1*(L2-L3)
    1 CONTINUE
*
*     CHECK ON NUMBER OF ELECTRONS
*
      IF ((N5-N6).NE.0) THEN
        WRITE(IWRITE,103)
        LET=0
*
*     CHECK ON PARITY
*
      ELSE
        IF((N7-N7/2*2).NE.0) THEN
          WRITE(IWRITE,104)
          LET=0
        ELSE
          N72=N7/2
          SIGNFA=1.D0
          IF( (N72-(N72/2)*2).NE.0 ) SIGNFA=-SIGNFA
*
* --- COUPLING ORTHOGONALITY TEST FOR FIRST TWO SHELLS
*
          LET=1
        ENDIF
      ENDIF
      RETURN
      END
*
*     ------------------------------------------------------------------
*       O R T H O G G B
*     ------------------------------------------------------------------
*
      SUBROUTINE ORTHOGGB(LET,INCL)
*
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
*
      COMMON/INFORM/IREAD,IWRITE,IOUT,ISC(8),ISCW
      COMMON/DEBUG/IBUG1,IBUG2,IBUG3,NBUG6,NBUG7,IFULL
      COMMON/IMAGNT/CONST,CONSOO,CONSS,ISPORB,ISOORB,ISPSPN,
     :     IREL,ISTRICT,IZOUT,IELST,ITENPR
      COMMON/MEDEFN/IHSH,NJ(16),LJ(16),NOSH1(16),NOSH2(16),J1QN1(31,3),
     :     J1QN2(31,3),IJFUL(16)
      LOGICAL INCL
*
*     THIS SUBROUTINE CHECKS FOR POSSIBLE ORTHOGONALITY DUE TO
*     COUPLING DIFFERENCES OR UNEVEN PARITY
*
  101 FORMAT(37H DIFFERING RESULTANT ANGULAR MOMENTUM)
  102 FORMAT(52H ORTHOGONALITY IN COUPLING SCHEMES OF CONFIGURATIONS)
  103 FORMAT(59H THE TWO CONFIGURATIONS HAVE DIFFERING NUMBERS OF ELECTR
     :ONS)
  104 FORMAT(51H THE TWO CONFIGURATIONS HAVE DIFFERING TOTAL PARITY)
*
* --- DO PSI AND PSIP CONTAIN THE SAME NUMBERS OF ELECTRONS
*     DO PSI AND PSIP HAVE THE SAME TOTAL PARITY
*
      N5=0
      N6=0
      N7=0
      IELST=1
      DO 20 I=1,IHSH
      L1=LJ(I)
      L2=NOSH1(I)
      L3=NOSH2(I)
      N5=N5+L2
      N6=N6+L3
      N7=N7+L1*(L2-L3)
   20 CONTINUE
*
*     CHECK ON NUMBER OF ELECTRONS
*
      IF (N5-N6) 21,22,21
   21 IF(IBUG2-1) 11,28,28
   28 WRITE(IWRITE,103)
      GO TO 11
*
*     CHECK ON PARITY
*
   22 IF(N7-N7/2*2) 23,24,23
   23 IF(IBUG2-1) 11,25,25
   25 WRITE(IWRITE,104)
      GO TO 11
   24 N1=2*IHSH-1
      N2=IHSH+1
      N3=IHSH-1
      N4=IHSH-2
      GO TO 3
*
* --- THE TWO CONFIGURATIONS WILL HAVE ZERO HAMILTONIAN MATRIX ELEMENT
*
   11 LET=0
      RETURN
    3 CONTINUE
*
* --- NO OBVIOUS ANGULAR MOMENTUM ORTHOGONALITY
*
   12 LET=1
      IF (IELST.EQ.0.AND. .NOT.INCL) LET = 0
      RETURN
      END
*
*     ------------------------------------------------------------------
*       S H E L L S
*     ------------------------------------------------------------------
*
      SUBROUTINE SHELLS(JA,JB,LET)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
*
      POINTER(QNOC,NOCCSH(1)),(QNELCSH,NELCSH(8,1)),
     :       (QNOCORB,NOCORB(8,1)),(QJ1,J1QNRD(15,1))
      POINTER(QIAJCMP,IAJCMP(1)),(QLJCOMP,LJCOMP(1)),
     :       (QNJCOMP,NJCOMP(1)),(QIAJCLD,IAJCLD(1)),
     :       (QLJCLSD,LJCLSD(1))
      COMMON /NDIMS/ QNOC,QNELCSH,QNOCORB,QJ1,NCFG
      COMMON /NON30/ QIAJCMP,QNJCOMP,QLJCOMP,QIAJCLD,QLJCLSD,MAXORB
      COMMON/MEDEFN/IHSH,NJ(16),LJ(16),NOSH(16,2),J1QN(31,3,2),
     :IJFUL(16)
      COMMON/OCCUPATION/NCG(16,2),ICG(16,2)
      COMMON/JB/IHSHJB,NJB(16),LJB(16),NOSHJB(16),IJFULJB(16)
C      WRITE(6,'(/A)') '  **************  iena  ******** '
C      WRITE(6,'(/A)') '  JA JB '
C      WRITE(6,'(2I5)') JA,JB
      LET=1
      IF(JA.EQ.JB) THEN
      CALL SHELLSAME(JA)
      RETURN
      ENDIF
      IA=NOCCSH(JA)
      IB=NOCCSH(JB)
      IF(IABS(IA-IB).GT.4) THEN
      LET=0
      RETURN
      ENDIF
      IAM=0
      IBM=0
      IGA=0
      IGGA=0
      IHSH=0
      JG1=1
      JG2=1
    1 J1=NOCORB(JG1,JA)
      J2=IJFULJB(JG2)
      IHSH=IHSH+1
C      WRITE(6,'(/A)') '  IHSH '
C      WRITE(6,'(I3)')  IHSH
C      WRITE(6,'(A)') '  JG1, JG2 '
C      WRITE(6,'(2I5)') JG1,JG2
C      WRITE(6,'(A)') '        J1, J2 '
C      WRITE(6,'(5X,2I5)') J1,J2
      IF(J1.EQ.J2) THEN
        NA=2
        NB=2
        JJ=J1
      ELSEIF(J1.GT.J2) THEN
        IF(IBM.EQ.0) THEN
          NA=1
          NB=2
          JJ=J2
	ELSE
          NA=2
          NB=1
          JJ=J1
	ENDIF
      ELSEIF(J1.LT.J2) THEN
        IF(IAM.EQ.0) THEN
          NA=2
          NB=1
          JJ=J1
	ELSE
          NA=1
          NB=2
          JJ=J2
	ENDIF
      ENDIF
      NCG(IHSH,1)=NA
      NCG(IHSH,2)=NB
      ICG(IHSH,1)=JG1
      ICG(IHSH,2)=JG2
      IF(NA.EQ.1) THEN
        NOSH(IHSH,1)=0
      ELSE
        NOSH(IHSH,1)=NELCSH(JG1,JA)
      ENDIF
      IF(NB.EQ.1) THEN
        NOSH(IHSH,2)=0
      ELSE
        NOSH(IHSH,2)=NOSHJB(JG2)
      ENDIF
C      WRITE(6,'(A)') '  nosh(1),nosh(2) '
C      WRITE(6,'(2I5)') NOSH(IHSH,1),NOSH(IHSH,2)
      NJ(IHSH)=NJCOMP(JJ)
C      WRITE(6,'(A)') '  NJ '
C      WRITE(6,'(I5)') NJ(IHSH)
      LJ(IHSH)=LJCOMP(JJ)
C      WRITE(6,'(A)') '  LJ '
C      WRITE(6,'(I5)') LJ(IHSH)
      IJFUL(IHSH)=JJ
      IF(JG1.LT.IA.OR.JG2.LT.IB) THEN
        IF(JG1.LT.IA) THEN
	  IAM=0
          IF(NA.EQ.2) JG1=JG1+1
        ELSE
          IF(NA.EQ.2)  THEN
	  IF(IGGA.EQ.0) THEN
	  IAM=1
	  IGG1=IGG1+1
        ENDIF
        ENDIF
        ENDIF
        IF(JG2.LT.IB) THEN
	  IBM=0
          IF(NB.EQ.2) JG2=JG2+1
        ELSE
          IF(NB.EQ.2)  THEN
	  IF(IGGA.EQ.0) THEN
	  IBM=1
	  IGG1=IGG1+1
        ENDIF
        ENDIF
        ENDIF
        GO TO 1
      ENDIF
      IF(NA.EQ.1) THEN
	 IF(J1.LT.J2)RETURN
	 IF(IAM.EQ.1) RETURN
         IAM=1
	 IBM=1
         IGA=IGA+1
         IF(IGA.EQ.1)GO TO 1
      ENDIF
      IF(NB.EQ.1) THEN
	 IF(J2.LT.J1)RETURN
	 IF(IBM.EQ.1) RETURN
         IAM=1
	 IBM=1
         IGA=IGA+1
         IF(IGA.EQ.1)GO TO 1
      ENDIF
      RETURN
      END
*
*     ------------------------------------------------------------------
*       S H E L L S J B
*     ------------------------------------------------------------------
*
      SUBROUTINE SHELLSJB(JB)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
*
      POINTER(QNOC,NOCCSH(1)),(QNELCSH,NELCSH(8,1)),
     :       (QNOCORB,NOCORB(8,1)),(QJ1,J1QNRD(15,1))
      POINTER(QIAJCMP,IAJCMP(1)),(QLJCOMP,LJCOMP(1)),
     :       (QNJCOMP,NJCOMP(1)),(QIAJCLD,IAJCLD(1)),
     :       (QLJCLSD,LJCLSD(1))
      COMMON /NDIMS/ QNOC,QNELCSH,QNOCORB,QJ1,NCFG
      COMMON /NON30/ QIAJCMP,QNJCOMP,QLJCOMP,QIAJCLD,QLJCLSD,MAXORB
      COMMON/OCCUPATION/NCG(16,2),ICG(16,2)
      COMMON/JB/IHSHJB,NJB(16),LJB(16),NOSHJB(16),IJFULJB(16)
      IB=NOCCSH(JB)
      IHSHJB=IB
      DO 1 J=1,IB
        I1=NOCORB(J,JB)
        NOSHJB(J)=NELCSH(J,JB)
        NJB(J)=NJCOMP(I1)
        LJB(J)=LJCOMP(I1)
        IJFULJB(J)=I1
    1 CONTINUE
      RETURN
      END
*
*     ------------------------------------------------------------------
*       S H E L L S A M E
*     ------------------------------------------------------------------
*
      SUBROUTINE SHELLSAME(JA)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      COMMON/JB/IHSHJB,NJB(16),LJB(16),NOSHJB(16),IJFULJB(16)
      COMMON/MEDEFN/IHSH,NJ(16),LJ(16),NOSH(16,2),J1QN(31,3,2),
     :IJFUL(16)
      COMMON/OCCUPATION/NCG(16,2),ICG(16,2)
      DO 2 J=1,IHSHJB
        NCG(J,1)=2
        NCG(J,2)=2
        ICG(J,1)=J
        ICG(J,2)=J
        NOSH(J,1)=NOSHJB(J)
        NOSH(J,2)=NOSHJB(J)
        NJ(J)=NJB(J)
        LJ(J)=LJB(J)
        IJFUL(J)=IJFULJB(J)
    2 CONTINUE
      IHSH=IHSHJB
      RETURN
      END
*
*     ------------------------------------------------------------------
*      C O U P L I N G
*     ------------------------------------------------------------------
*
      SUBROUTINE COUPLING(JA,JB)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
*
      POINTER(QNOC,NOCCSH(1)),(QNELCSH,NELCSH(8,1)),
     :       (QNOCORB,NOCORB(8,1)),(QJ1,J1QNRD(15,1))
      COMMON /NDIMS/ QNOC,QNELCSH,QNOCORB,QJ1,NCFG
      COMMON/MEDEFN/IHSH,NJ(16),LJ(16),NOSH(16,2),J1QN(31,3,2),
     :IJFUL(16)
      COMMON/OCCUPATION/NCG(16,2),ICG(16,2)
      DO 1 IH=1,IHSH
        JC=JA
        DO 2 I=1,2
          IC=ICG(IH,I)
          NC=NCG(IH,I)
          I2H=IHSH+IH-1
*
* --- FIRST CONSIDER THE L.H.S. (I=1) OF THE MATRIX ELEMENT. NC=1 MEANS
*     UNOCCUPIED, REPRESENTED BY A DUMMY SINGLET S SHELL, AND THE
*    ADDITIONAL SET OF COUPLING QUANTUM NUMBERS WILL BE THE SAME AS THE
*     LAST SET OF COUPLING QUANTUM NUMBERS ALREADY OBTAINED.
*     NC=2 MEANS OCCUPIED.  THEN ALL THE NEW QUANTUM NUMBERS (BOTH FOR
*     THE SHELL AND FOR THE COUPLING OF THIS SHELL TO THE RESULTANT OF
*     THE PREVIOUS ONES) ARE DEFINED IN THE CORRESPONDING J1QNRD ARRAY.
*     NOSH - THE NUMBER OF ELECTRONS IN THIS SHELL, IS DEFINED BY THE
*     APPROPRIATE ENTRY IN NELCSH .  THE R.H.S. IS THEN CONSIDERED
*     SIMILARLY (I=2)
*
          IF(NC.EQ.1) THEN
            J1QN(IH,1,I)=0
            J1QN(IH,2,I)=1
            J1QN(IH,3,I)=1
            IF(IH.EQ.2) THEN
              J1QN(I2H,1,I)=0
              J1QN(I2H,2,I)=J1QN(1,2,I)
              J1QN(I2H,3,I)=J1QN(1,3,I)
            ELSEIF(IH.GT.2) THEN
	      I2H1=I2H-1
              J1QN(I2H,1,I)=J1QN(I2H1,1,I)
              J1QN(I2H,2,I)=J1QN(I2H1,2,I)
              J1QN(I2H,3,I)=J1QN(I2H1,3,I)
            END IF
          ELSE
          JD = J1QNRD(IC,JC)
          J1QN(IH,1,I)=MOD(JD,64)
          JD = JD/64
          J1QN(IH,2,I) = MOD(JD,64)
          J1QN(IH,3,I) = JD/64
*
*     IS THIS THE FIRST OCCUPIED SHELL OF EITHER CONFIGURATION. IF SO,
*    THEN THERE ARE NO INTERMEDIATE COUPLINGS TO CONSIDER AT THIS STAGE
*
          IF(IH .GT. 1) THEN
*
*    IS THIS THE FIRST OCCUPIED SHELL OF THIS CONFIGURATION, THOUGH NOT
*     THE FIRST OF THE OTHER CONFIGURATION.  IF SO, THE INTERMEDIATE
*     COUPLING FORMED HAS THE SAME  L,S  VALUES AS THIS OCCUPIED SHELL,
*     SINCE WE COUPLE THE SHELL TO A DUMMY SINGLET S.
*
            IF(IC .LE.1) THEN
              I2 = 1
            ELSE
              I2 = NOCCSH(JC)+IC-1
            END IF
            JD = J1QNRD(I2,JC)
            IF (IC .LE. 1) THEN
              J1QN(I2H,1,I) = 0
            ELSE
              J1QN(I2H,1,I) = MOD(JD,64)
            END IF
            JD = JD/64
            J1QN(I2H,2,I) = MOD(JD, 64)
            J1QN(I2H,3,I) = JD/64
          END IF
        END IF
        JC=JB
    2   CONTINUE
    1 CONTINUE
      END
*
*     -------------------------------------------------------------
*      S A V E N O N
*     -------------------------------------------------------------
*                                                                  *
*     THIS SUBROUTINE FOR         G N O N H                        *
*                                                                  *
*     Written by G. Gaigalas,                                      *
*     Vanderbilt University,  Nashville            February 1994   *
*
      SUBROUTINE SAVENON(I,A,KL,LA,LB,LC,LD,JA,JB,IPTR)
      IMPLICIT DOUBLEPRECISION (A-H,O-Z)
      COMMON/SIGNF/SIGNFA
      A=A*SIGNFA
      CALL SAVELS(I,A,KL,LA,LB,LC,LD,JA,JB,IPTR)
      RETURN
      END


text:
      subroutine harris(jant,fileant,mode,idx,antfile,fs,fe,
     +                  beam_main,offazim)
c*******************************************************************
c     11/98 Greg Hand modifications for Harris
c
c     Given an external antenna type xx (90-99),
c     1. create a data file ANTTYPxx.DAT
c     2. execute the program ..\BIN_WIN\ANTTYPxx.EXE rundirectory mode
c        where:
c              rundirectory = C:\ITSHFBC\RUN (default installation)
c              mode = (blank) = point-to-point
c                   = a = area coverage
c*******************************************************************
      use voacapl_defs
      use crun_directory
      CHARACTER fileant*10,antfile*21,mode*1
c jw      common /crun_directory/ run_directory
c jw         character run_directory*50
      dimension areagain(91)
      character program*300,alf*1,exe*12,dat*12
      logical*1 iexe_exist
      nch_run=lcount(run_directory,50)

#ifdef _WIN
      write(exe,'(6hANTTYP,i2.2,4h.EXE)') jant
      write(dat,'(6hANTTYP,i2.2,4h.DAT)') jant
#else
      write(exe,'(6hanttyp,i2.2)') jant
      write(dat,'(6hanttyp,i2.2,4h.dat)') jant
#endif

C
C**********Harris modification M Packer for Area Coverage 11 May 1998 *******
C
c        This handles the Harris Type 99 antenna files.  Call the external
c        routines (written in C) to prepare the
c        GAINxx.DAT for area coverage.  Then read in the file.
c        We're gonna use numbers 250-299 if necessary.

c           Follow the same routine necessary for AREA type 13 antenna files.

c           The difficulty here is that antcalc is writing the GAINxx.DAT
c           file on lu 22 and we want to call a routine to generate the
c           GAINxx.DAT file.

c           Prepare ANTTYPxx.DAT to communicate with ANTTYPxx.EXE:
      open( 99, file=run_directory(1:nch_run)//PATH_SEPARATOR//dat,err=900)
      write( 99, 289 ) idx, antfile, fs,fe,beam_main,offazim
289   format( i2, 28x,   '=antenna number (idx)' /
     2       a21,  9x,  '=antenna file filename' /
     3      f7.3, 23x, '=antenna starting frequency (MHz)' /
     4      f7.3, 23x, '=antenna ending frequency (MHz)' /
     5      f7.2, 23x, '=antenna main beam (deg)' /
     6      f7.2, 23x, '=antenna off azimuth angle (deg)' )
      close( 99 )

c**********************************************************************
c          Make sure ANTTYPxx.EXE exists before executing it.
c**********************************************************************
c jw      iexe_exist=it_exist(run_directory(1:nch_run-3)//'bin_win\'//exe)
#ifdef _WIN
      inquire(file=run_directory(1:nch_run-3)//'bin_win'//PATH_SEPARATOR//exe,exist=iexe_exist)
c jw      if(iexe_exist.ne.1) then    !  antenna program does not exist
      if(.NOT.iexe_exist) then    !  antenna program does not exist
         write(*,1) jant,run_directory(1:nch_run-3)//'bin_win\'//exe
1        format('You are attempting to use an antenna of type=',i4,'.',/
     +   'That requires an external antenna calculation program.',/
     +   'The required program=',a,'  does not exist.'/
     +   'Please correct the problem and try again.'//)
         stop 'Required external antenna program does not exist.'
      end if
c**********************************************************************
C        Call ANTTYPxx.EXE to create GAINxx.DAT file:
C
      PROGRAM=run_directory(1:nch_run-3)//'bin_win'//PATH_SEPARATOR//exe//' '//
     +             run_directory(1:nch_run)//' '//mode
#else
      PROGRAM=exe//' '//trim(run_directory)//' '//trim(root_directory)//' '//mode
#endif
      nch=lcount(PROGRAM,300)
c jw      call gh_exec(PROGRAM,nch,1)   !  execute and wait for ANTTYPxx.EXE
      call system(PROGRAM)
      if(mode.eq.' ') return     !  point-to-point is done

C           Open GAINxx.DAT file
      open(22,file=run_directory(1:nch_run)//PATH_SEPARATOR//fileant,
     1            status='old', err=910)
      rewind(22)
      read(22,'(a)') alf
      read(22,'(a)') alf
      read(22, 292 ) freq, eff
292   format( 10x, f7.3, 10x, f10.3 )

c           Read the rest of the file which contains antenna gain info
      do iazim=0,359
c jw         call yieldit
         read( 22, 293 ) ( areagain(ielev), ielev=1,91 )
293      format( 9x, 10f7.3 )
         call antsave( idx, iazim, freq, eff, areagain )
      end do

c           We are done.  Close up, take your ball, and go home.
      close( 22 )
      return
c*****End Harris mod      *************************************************
900   write(*,901) run_directory(1:nch_run)//dat
901   format('Could not create file:',a)
      go to 920
910   write(*,911) run_directory(1:nch_run)//PATH_SEPARATOR//fileant
911   format('Could not OPEN file:',a)
920   write(*,921)
921   format('Failure in subroutine harris')
      stop 'Fatal error in subroutine harris'
      end
c-------------------------------------------------------

max_stars_repo_name: jawatson/voacapl

max_issues_repo_name: jawatson/voacapl

max_forks_repo_name: jawatson/voacapl


text:
Image(briggs 4.jpg, right, thumbnail, Delta, 2007.) Image(briggs 2.jpg, thumbnail, left, At Delta after a solo performance.)Image(briggs 3.jpg, thumbnail, right, Performing at Delta of Venus in 2007.) Image(briggs.JPG, left, thumbnail, 200, Teddy performing at Delta of Venus with Users/JoeFinkel and Users/AndyPastalaniec on July 14, 2007 (not under his Chief moniker.)) Chief Briggum is the stage name of Users/TeddyBriggs when he performs solo. The Chief has performed at Delta of Venus, Operation Restore Maximum Freedom, and at other venues in Davis and Sacramento. The Chief plays drums and sings in Boss the Big Bit when not venturing out on his own. He is/was also a KDVS DJ.

More information on Chief Briggum may be found http://www.myspace.com/chiefbriggum here.