PROGRAM PLOTFILE
C=======================================================================
C FINITE ELEMENT GRAPHICS FOR TWO-DIMENSIONAL CFD SOLUTIONS
C ELEMENT TYPE: FOUR-NODED ISOPARAMETRIC ELEMENT
C EIJI FUKUMORI, JULY 1985, JAN. 2013
C----------------- DIVV(I)=MISES-HENCKY CRITERION ----------------------
C=======================================================================
INCLUDE 'PARAM.DAT'
CCCCC PARAMETER ( MXE=24000,MXN=25000,MXB=6000, ND=4 )
C=======================================================================
IMPLICIT REAL*8 ( A-H , O-Z )
C ARRAYS
CHARACTER*12 INPFILE
DIMENSION NODEX(MXE,ND)
DIMENSION IBNDFX(MXB),IBNDFY(MXB),BVX(MXB),BVY(MXB)
DIMENSION XCOORD(MXN),YCOORD(MXN),U(MXN),V(MXN),
* DIVV(MXN),TAUXX(MXN),TAUYY(MXN),TAUXY(MXN),TAU1(MXN),TAU2(MXN)
C=======================================================================
DATA INPFILE /'STATIC04.DAT'/
C=======================================================================
WRITE(*,*)' STATIC GRAPHICS PROGRAM'
WRITE (*,*)' READING IN DATA FILES'
CALL INPUT ( MXE,MXN,MXB,ND,NE,NNODE,NBFX,NBFY,YOUNG,POISSON,
* NODEX,XCOORD,YCOORD,U,V,IBNDFX,IBNDFY,BVX,BVY,NF,INPFILE,
* DIVV, TAUXX,TAUYY,TAUXY,TAU1, TAU2 )
C=======================================================================
CALL PLTEL4 ( ND,MXE,MXN,NODEX,NE,XCOORD,YCOORD, NNODE )
C
CALL PLTUV ( MXN,NNODE,XCOORD,YCOORD,U,V,NE,MXE,ND,NODEX )
C
CALL PLTLGO ("MISES000.OUT")
CALL CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,YCOORD,DIVV)
CALL PLTEXT
C
CALL PLTLGO ("STRESSXX.OUT")
CALL CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,YCOORD,TAUXX)
C
CALL PLTLGO ("STRESSYY.OUT")
CALL CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,YCOORD,TAUYY)
CALL PLTEXT
C
CALL PLTLGO ("STRESSXY.OUT")
CALL CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,YCOORD,TAUXY)
CALL PLTEXT
C
CALL PLTLGO ("STRESSXY.OUT")
CALL CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,YCOORD,TAUXY)
CALL PLTEXT
C
CALL PLTLGO ("STRESS01.OUT")
CALL CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,YCOORD,TAU1)
C
CALL PLTLGO ("STRESS02.OUT")
CALL CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,YCOORD,TAU2)
CALL PLTEXT
C
STOP 'TERMINATION'
END
C
C
SUBROUTINE INPUT ( MXE,MXN,MXB,ND,NE,NNODE,NBFX,NBFY,YOUNG,
* POISSON,NODEX,XCOORD,YCOORD,U,V,IBNDFX,IBNDFY,BVX,BVY,
* NF,INPFILE,DIVV, TAUXX,TAUYY,TAUXY,TAU1,TAU2 )
IMPLICIT REAL*8 ( A-H , O-Z )
DIMENSION NODEX(MXE,ND),XCOORD(MXN),YCOORD(MXN),U(MXN),
* V(MXN),IBNDFX(MXB),IBNDFY(MXB),BVX(MXB),BVY(MXB)
DIMENSION DIVV(MXN), TAUXX(MXN), TAUYY(MXN), TAUXY(MXN),
* TAU1(MXN), TAU2(MXN)
CHARACTER*12 INPFILE, STRESSFL, DSPFILE
CHARACTER EXFILE*4
LOGICAL YES
WRITE (*,*)' READING IN STATIC.DAT DATA FILES'
C========> FILENAME INPFILE SEE MAIN PROGRAM
INQUIRE ( FILE=INPFILE, EXIST=YES )
IF ( YES ) THEN
OPEN ( 1, FILE=INPFILE, STATUS='OLD' )
ELSE
WRITE (*,*)' INPUT FILE DOES NOT EXIST'
STOP
ENDIF
C========> PARAMETERS
READ (1,*) YOUNG, POISSON
C========> ELEMENTS
READ (1,*) NE
IF ( NE .GT. MXE) STOP'ERROR #1'
IF ( NE .LE. 0 ) STOP'NE =0'
DO I = 1 , NE
READ (1,*) IEL, ( NODEX(IEL,J), J = 1 , ND )
END DO
C========> FILENAME COORDINATES OF NODAL POINTS
READ (1,*) NNODE
IF ( NNODE .GT. MXN) STOP'ERROR #2'
IF ( NNODE .LE. 0 ) STOP'NNODE =0'
DO I = 1 , NNODE
READ (1,*) NODE, XCOORD(NODE) , YCOORD(NODE)
END DO
C========> DIRICHLET TYPE BOUNDARY CONDITIONS
READ (1,*) NBFX
IF ( NBFX .GT. MXB) STOP'NBFX>MXB'
IF ( NBFX .LT. 1 ) STOP'NBFX<1'
DO I = 1 , NBFX
READ (1,*) IBNDFX(I) , BVX(I)
END DO
READ (1,*) NBFY
IF ( NBFY .GT. MXB) STOP'NBFY .GT. MXB'
IF ( NBFY .LT. 1 ) STOP'NBFY<1'
DO I = 1 , NBFY
READ (1,*) IBNDFY(I) , BVY(I)
END DO
CLOSE (1)
C========> FILENAME STRESS DATA: TAUXX,TAUYY,TAUXY,DIVV,TAI1,TAU2
C------------------ DIVV=MISES-HENCKY CRITERION -------------------
WRITE (*,*)' READING IN STRESS.DAT DATA FILES'
STRESSFL = 'STRESS04.DAT'
INQUIRE ( FILE=STRESSFL, EXIST=YES )
IF ( YES ) THEN
OPEN ( 1, FILE=STRESSFL, STATUS='OLD' )
ELSE
WRITE (*,*)' STRESS FILE DOES NOT EXIST'
STOP
ENDIF
READ(1,*)
READ(1,*)
READ(1,*)
DO NODE = 1 , NNODE
READ (1,*)I,TAUXX(I),TAUYY(I),TAUXY(I),DIVV(I),TAU1(I),TAU2(I)
END DO
CLOSE (1)
C========> FILENAME DISPLACE.MNT
DSPFILE = 'DISPLACE.MNT'
EXFILE = 'NEW'
INQUIRE ( FILE=DSPFILE, EXIST=YES )
IF ( YES ) THEN
OPEN ( 1, FILE=DSPFILE, STATUS='OLD',FORM='UNFORMATTED' )
ELSE
WRITE (*,*)' DISPLACE.MNT FILE DOES NOT EXIST'
STOP
ENDIF
READ (1) ( U(I) , I = 1 , NNODE )
READ (1) ( V(I) , I = 1 , NNODE )
CLOSE (1)
RETURN
END
C
C
SUBROUTINE PLTUV (MXN,NNODE,XCOORD,YCOORD,U,V,NE,MXE,ND,NODEX )
IMPLICIT REAL*8 ( A-H , O-Z )
DIMENSION XCOORD(MXN),YCOORD(MXN),U(MXN),V(MXN),NODEX(MXE,ND)
C
CALL MINMAX ( MXN, NNODE, XCOORD, XMIN , XMAX )
CALL MINMAX ( MXN, NNODE, YCOORD, YMIN , YMAX )
DATA AL / 0.9D0 /
RATIO = 0.05
RMAX = 0.
DO I = 1 , NNODE
RMAX = DMAX1 ( RMAX, U(I)*U(I)+V(I)*V(I) )
END DO
IF ( RMAX .EQ. 0. ) RETURN
RMAX = DSQRT ( RMAX )
CALL MINMAX ( MXN, NNODE, U, UMIN, UMAX )
CALL MINMAX ( MXN, NNODE, V, VMIN, VMAX )
FACT = RATIO * DMAX1 ( (YMAX-YMIN) , (XMAX-XMIN) ) / RMAX
C-------- START PLOT
CALL PLTLGO ("DISPLACE.OUT")
DO IEL = 1 , NE
DO I = 1 , ND
J = I+1
IF ( I .EQ. ND ) J = 1
CALL XMOVE ( XCOORD(NODEX(IEL,I)),YCOORD(NODEX(IEL,I)) )
CALL XDRAW ( XCOORD(NODEX(IEL,J)),YCOORD(NODEX(IEL,J)) )
END DO
END DO
C
DO IEL = 1 , NE
DO I = 1 , ND
J = I+1
IF ( I .EQ. ND ) J = 1
UUI = XCOORD(NODEX(IEL,I))+U(NODEX(IEL,I))*FACT*AL
VVI = YCOORD(NODEX(IEL,I))+V(NODEX(IEL,I))*FACT*AL
UUJ = XCOORD(NODEX(IEL,J))+U(NODEX(IEL,J))*FACT*AL
VVJ = YCOORD(NODEX(IEL,J))+V(NODEX(IEL,J))*FACT*AL
CALL XMOVE ( UUI,VVI )
CALL XDRAW ( UUJ,VVJ )
END DO
END DO
CALL PLTEXT
C
RETURN
END
C
C
SUBROUTINE BOUND (MXE,MXN,ND,NE,NODEX,XCOORD,YCOORD)
IMPLICIT REAL*8 ( A-H , O-Z )
DIMENSION XCOORD(MXN),YCOORD(MXN),NODEX(MXE,ND)
LOGICAL LINE
C
DO IEL = 1 , NE
DO I = 1 , ND
LINE = .TRUE.
J = I + 1
IF ( I .EQ. ND ) J = 1
C
DO JEL = 1 , NE
IF ( IEL .NE. JEL ) THEN
DO II = 1 , ND
JJ = II + 1
IF ( II .EQ. ND ) JJ = 1
IF ( NODEX(IEL,I).EQ.NODEX(JEL,JJ) .AND.
* NODEX(IEL,J).EQ.NODEX(JEL,II) ) THEN
LINE = .FALSE.
EXIT
END IF
END DO
END IF
IF ( .NOT. LINE ) EXIT
END DO
IF ( LINE ) THEN
CALL XMOVE ( XCOORD(NODEX(IEL,I)) , YCOORD(NODEX(IEL,I)) )
CALL XDRAW ( XCOORD(NODEX(IEL,J)) , YCOORD(NODEX(IEL,J)) )
END IF
C
END DO
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
C
C
SUBROUTINE PLTEL4 (ND,MXE,MXN,NODEX,NE,XCOORD,YCOORD,NNODE)
IMPLICIT REAL*8 ( A-H , O-Z )
DIMENSION NODEX(MXE,ND),XCOORD(MXN),YCOORD(MXN)
CALL PLTLGO ("ELEMENT4.OUT")
DO I = 1 , NE
DO J = 1 , ND-1
CALL XMOVE ( XCOORD(NODEX(I,J )),YCOORD(NODEX(I,J )) )
CALL XDRAW ( XCOORD(NODEX(I,J+1)),YCOORD(NODEX(I,J+1)) )
END DO
CALL XMOVE ( XCOORD(NODEX(I,ND )),YCOORD(NODEX(I,ND )) )
CALL XDRAW ( XCOORD(NODEX(I, 1 )),YCOORD(NODEX(I, 1 )) )
END DO
CALL PLTEXT
RETURN
END
C
C
SUBROUTINE CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,YCOORD,P)
IMPLICIT REAL*8 ( A-H , O-Z )
DIMENSION NODEX(MXE,ND),XCOORD(MXN),YCOORD(MXN),P(MXN)
CALL MINMAX ( MXN, NNODE, P, PPMIN, PPMAX )
IF ( PPMAX .EQ. PPMIN ) STOP 'PPMAX=PPMIN AT CONTOUR'
C-----------------------------------------------------------------------
NSTEP = 22
C-----------------------------------------------------------------------
NSTEP = NSTEP/2*2 + 1
DS = ( PPMAX - PPMIN ) / NSTEP
PPMIN = INT ( PPMIN/DS ) * DS
CALL BOUND ( MXE,MXN,ND,NE,NODEX,XCOORD,YCOORD )
C
DO IEL = 1 , NE
C
DO LEVEL = 1 , NSTEP
SXY = PPMIN + (LEVEL-1) * DS
K = 1
DO I = 1 , 4
J = I + 1
IF ( I .EQ. 4 ) J = 1
S1 = P(NODEX(IEL,I))
S2 = P(NODEX(IEL,J))
IF ( (S1-SXY)*(S2-SXY) .LT. 0 ) THEN
T=(SXY-S1)/(S2-S1)
X1 = XCOORD(NODEX(IEL,I))
X2 = XCOORD(NODEX(IEL,J))
Y1 = YCOORD(NODEX(IEL,I))
Y2 = YCOORD(NODEX(IEL,J))
X0=X2*T+(1.D0-T)*X1
Y0=Y2*T+(1.D0-T)*Y1
IF ( K .GT. 0 ) THEN
CALL XMOVE ( X0, Y0 )
ELSE
CALL XDRAW ( X0, Y0 )
END IF
K = -K
END IF
END DO
END DO
C
END DO
RETURN
END
C
C
SUBROUTINE PLTLGO (FILENAME)
IMPLICIT REAL*8 ( A-H , O-Z )
CHARACTER*12 FILENAME
OPEN ( 1, FILE=FILENAME, STATUS='UNKNOWN' )
RETURN
END
C
C
SUBROUTINE PLTEXT
CLOSE (1)
RETURN
END
C
C
SUBROUTINE XMOVE ( X , Y )
IMPLICIT REAL*8 ( A-H , O-Z )
WRITE(1,*) X , Y
RETURN
END
C
C
SUBROUTINE XDRAW ( X , Y )
IMPLICIT REAL*8 ( A-H , O-Z )
WRITE(1,*) X , Y
WRITE(1,*)
RETURN
END