PROGRAM GRA9FILE
C=======================================================================
C   FINITE ELEMENT GRAPHICS FOR TWO-DIMENSIONAL TORSION SOLUTIONS
C            ELEMENT TYPE: 9-NODED ISOPARAMETRIC ELEMENT
C             CREATE FILES CONTAINING XY-COORDINTE VALUES
C         ORIGINAL CODE: NSEQ8GRA.FOR, EIJI FUKUMORI, JULY 1985
C 15 FEB. 2013, 03 NOV. 2024
C=======================================================================
      PARAMETER ( ND=9, MXE=3130, MXN=3441, MXB=2100 )
CCCCCCCCCCCC      INCLUDE 'PARAM.DAT'
      IMPLICIT REAL*8 ( A-H , O-Z )
C                                ARRAYS
      CHARACTER*14 INPFILE, FN
      DIMENSION NODEX(MXE,ND), XCOORD(2,MXN), U(MXN)
      DIMENSION TAUZX(MXN),TAUZY(MXN)
C=======================================================================
      INPFILE = 'TRS9DATA.DAT'
      WRITE(*,*)'TORSION9 GRAPHICS PROGRAM'
      WRITE (*,*)' READING IN DATA FILES'
      CALL INPUT ( INPFILE,ND,MXE,MXN,NE,NNODE,GM, THETA,NODEX,
     * XCOORD,U, TAUZX,TAUZY  )
      WRITE (*,*)'PROJECT NAME =======>', INPFILE
C=======================================================================
      NSTEP = 20
      FN = "TORSION9.OUT"
      CALL TEMPA (FN,MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,U,NSTEP)
      FN = "STRESSZX9.OUT"
      CALL TEMPA (FN,MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,TAUZX,NSTEP)
      FN = "STRESSZY9.OUT"
      CALL TEMPA (FN,MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,TAUZY,NSTEP)
C=======================================================================
      STOP 'NORMAL TERMINATION'
      END
C
C
      SUBROUTINE INPUT ( INPFILE,ND,MXE,MXN,NE,NNODE,GM, THETA,NODEX,
     * XCOORD,U, TAUZX,TAUZY  )
      IMPLICIT REAL*8 ( A-H , O-Z )
      DIMENSION NODEX(MXE,ND),XCOORD(2,MXN),U(MXN),TAUZX(MXN),TAUZY(MXN)
      CHARACTER INPFILE*14
      OPEN ( 1, FILE = INPFILE, STATUS = 'UNKNOWN' )
      READ (1,*) GM, THETA
      READ (1,*) NE
      IF ( NE .GT. MXE ) STOP 'NE > MXE'
      READ (1,*) (IEL,(NODEX(IEL,J),J=1,ND), I=1,NE)
      READ (1,*) NNODE
      IF ( NNODE .GT. MXN ) STOP 'NNODE > MXN'
      READ (1,*) (NODE,XCOORD(1,NODE),XCOORD(2,NODE),J=1,NNODE)
      CLOSE (1)
C=======================================================================
C========> FILENAME SOLUTION9.BIN
      OPEN (1,FILE="SOLUTION9.BIN",STATUS='UNKNOWN',FORM='UNFORMATTED')
      READ (1) ( U(I) , I = 1 , NNODE )
      READ (1) ( TAUZX(I) , I = 1 , NNODE )
      READ (1) ( TAUZY(I) , I = 1 , NNODE )
      CLOSE (1)
      RETURN
      END
C
C
      SUBROUTINE BOUND (MXE,MXN,ND,NE,NODEX,XCOORD)
      IMPLICIT REAL*8 ( A-H , O-Z )
      DIMENSION XCOORD(2,MXN),NODEX(MXE,ND)
      LOGICAL LINE
C--------- FOR USE OF ND=8 OR 9 ELEMENTS ---------------
      NDD = 7
      DO IEL = 1 , NE
C
      DO I = 1 , NDD, 2
      LINE = .TRUE.
      J = I + 1
      K = J + 1
      IF ( I .EQ. 7 ) K = 1
C
      DO JEL = 1 , NE
      IF ( IEL .NE. JEL ) THEN
      DO II = 1 , NDD, 2
      JJ = II + 1
      KK = JJ + 1
      IF ( II .EQ. 7 ) KK = 1
      IF ( NODEX(IEL,I).EQ.NODEX(JEL,KK) .AND.
     *     NODEX(IEL,K).EQ.NODEX(JEL,II)        ) THEN
      LINE = .FALSE.
      EXIT
      END IF
      END DO
      END IF
      IF ( .NOT. LINE ) EXIT
      END DO
C
      IF ( LINE ) THEN
      WRITE (1,*) XCOORD(1,NODEX(IEL,I)) , XCOORD(2,NODEX(IEL,I))
      WRITE (1,*) XCOORD(1,NODEX(IEL,J)) , XCOORD(2,NODEX(IEL,J))
      WRITE (1,*) XCOORD(1,NODEX(IEL,K)) , XCOORD(2,NODEX(IEL,K))
      WRITE (1,*)
      END IF
      END DO
C
      END DO
      RETURN
      END
C
C
      SUBROUTINE TEMPA (FN,MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,P,NSTEP )
      IMPLICIT REAL*8 ( A-H , O-Z )
      DIMENSION NODEX(MXE,ND),XCOORD(2,MXN),P(MXN)
      CHARACTER FN*14
      OPEN ( 1, FILE=FN, STATUS='UNKNOWN' )
      CALL CONTOUR ( MXE, MXN,ND,NE,NNODE,NODEX,XCOORD,P,NSTEP)
      CLOSE (1)
      RETURN
      END
C
C
      SUBROUTINE CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,P,NSTEP)
      IMPLICIT REAL*8 ( A-H , O-Z )
      DIMENSION NODEX(MXE,ND),XCOORD(2,MXN),P(MXN),S(4),CRD(2,4),
     *   IFOUR(4,4)
      DATA (IFOUR(1,J),J=1,4) /1, 2, 9, 8 /
      DATA (IFOUR(2,J),J=1,4) /2, 3, 4, 9 /
      DATA (IFOUR(3,J),J=1,4) /8, 9, 6, 7 /
      DATA (IFOUR(4,J),J=1,4) /9, 4, 5, 6 /
C
      CALL MINMAX ( MXN, NNODE, P, PPMIN, PPMAX )
      IF ( PPMAX .EQ. PPMIN ) RETURN
      NSTEP = NSTEP/2*2 + 1
      DS = ( PPMAX - PPMIN ) / NSTEP
      I = PPMIN/DS
      START= I*DS
C
      CALL BOUND ( MXE,MXN,ND,NE,NODEX,XCOORD )
      DO IEL = 1 , NE
C********* 9-NODED ELEMENT DIVIDED INTO FOUR SUB-ELEMENTS ************
      DO K = 1 , 4
      DO L = 1 , 4
      I = IFOUR(K,L)
      CRD(1,L) = XCOORD(1,NODEX(IEL,I))
      CRD(2,L) = XCOORD(2,NODEX(IEL,I))
      S(L)   =        P(NODEX(IEL,I))
      END DO
      CALL PLTSAI ( DS, NSTEP, START, CRD, S )
      END DO
C********************************************************************
      END DO
      RETURN
      END
C
C
      SUBROUTINE PLTSAI ( DS, NSTEP, START, CRD, S )
      IMPLICIT REAL*8 ( A-H , O-Z )
      DIMENSION CRD(2,4), S(4)
      IF ( NSTEP .EQ. 0 ) RETURN
      DO LEVEL = 1 , NSTEP
        SXY = START + (LEVEL-1) * DS
        K = 1
         DO I = 1 , 4
          J = I + 1
          IF ( I .EQ. 4 ) J = 1
           IF ( (S(I)-SXY)*(S(J)-SXY) .LT. 0.D0 ) THEN
           T=(SXY-S(I))/(S(J)-S(I))
           X0=CRD(1,J)*T+(1.D0-T)*CRD(1,I)
           Y0=CRD(2,J)*T+(1.D0-T)*CRD(2,I)
           WRITE (1,*) X0, Y0
           K = 0
          END IF
         END DO
       IF ( K .EQ. 0 ) WRITE (1,*)
      END DO
      RETURN
      END
C
C
      SUBROUTINE MINMAX ( MXN, NNODE, Q, QMIN, QMAX )
      IMPLICIT REAL*8 ( A-H , O-Z )
      DIMENSION Q(MXN)
      QMIN = Q(1)
      QMAX = Q(1)
      DO I = 1 , NNODE
      QMIN = DMIN1 ( QMIN , Q(I) )
      QMAX = DMAX1 ( QMAX , Q(I) )
      END DO
      RETURN
      END