PROGRAM GRA8FILE
C=======================================================================
C FINITE ELEMENT GRAPHICS FOR TWO-DIMENSIONAL FEM9 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
C=======================================================================
INCLUDE 'PARAM.DAT'
IMPLICIT REAL*8 ( A-H , O-Z )
CHARACTER*14 INPFILE
DIMENSION NODEX(MXE,ND), XCOORD(2,MXN), U(MXN), F(ND)
C=======================================================================
WRITE(*,*)' FEM8 GRAPHICS PROGRAM'
WRITE (*,*)' READING IN DATA FILES'
CALL INPUT ( INPFILE,ND,MXE,MXN,NE,NNODE,NODEX,XCOORD,U )
WRITE (*,*)'PROJECT NAME =======>', INPFILE
C=======================================================================
NSTEP = 20
CALL TEMPA (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,U,NSTEP,F )
C=======================================================================
STOP 'TERMINATION'
END
C
C
SUBROUTINE INPUT ( INPFILE,ND,MXE,MXN,NE,NNODE,NODEX,XCOORD,U )
IMPLICIT REAL*8 ( A-H , O-Z )
DIMENSION NODEX(MXE,ND),XCOORD(2,MXN), U(MXN)
CHARACTER INPFILE*14
IF ( ND .LE. 2 ) STOP 'NOT ACCEPTABLE ND'
IF ( ND .EQ. 3 ) INPFILE = 'FEM03INPUT.DAT'
IF ( ND .EQ. 4 ) INPFILE = 'FEM04INPUT.DAT'
IF ( ND .EQ. 8 ) INPFILE = 'FEM08INPUT.DAT'
IF ( ND .EQ. 9 ) INPFILE = 'FEM09INPUT.DAT'
IF ( ND .EQ. 12 ) INPFILE = 'FEM12INPUT.DAT'
OPEN ( 1, FILE = INPFILE, STATUS = 'UNKNOWN' )
READ (1,*) EXX, EXY, EYY
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 XXXXXXX.BIN
OPEN (1,FILE="SOLUTION.BIN",STATUS='UNKNOWN',FORM='UNFORMATTED')
READ (1) ( U(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 (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,P,NSTEP,F )
IMPLICIT REAL*8 ( A-H , O-Z )
DIMENSION NODEX(MXE,ND),XCOORD(2,MXN),P(MXN),F(ND)
CHARACTER FILENAME*12
FILENAME = "TEMPEAT.OUT"
CALL PLTLGO ( FILENAME )
CALL CONTOUR ( MXE, MXN,ND,NE,NNODE,NODEX,XCOORD,P,NSTEP,F)
CALL PLTEXT
RETURN
END
C
C
SUBROUTINE CONTOUR (MXE,MXN,ND,NE,NNODE,NODEX,XCOORD,P,NSTEP,F)
IMPLICIT REAL*8 ( A-H , O-Z )
DIMENSION NODEX(MXE,ND),XCOORD(2,MXN),P(MXN),S(4),B(2,4),
* IFOUR(4,4),F(ND)
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
PPMIN= I*DS
C
CALL BOUND ( MXE,MXN,ND,NE,NODEX,XCOORD )
CALL ISOPARA ( ND, 0.D0 , 0.D0 , F )
DO IEL = 1 , NE
X0 = 0.D0
Y0 = 0.D0
P0 = 0.D0
DO J = 1 , ND
X0 = X0 + F(J)*XCOORD(1,NODEX(IEL,J))
Y0 = Y0 + F(J)*XCOORD(2,NODEX(IEL,J))
P0 = P0 + F(J)*P(NODEX(IEL,J))
END DO
C********* 8-NODED ELEMENT DIVIDED INTO FOUR SUB-ELEMENTS ************
DO K = 1 , 4
DO L = 1 , 4
I = IFOUR(K,L)
IF ( I .EQ. 9 ) THEN
B(1,L) = X0
B(2,L) = Y0
S(L) = P0
ELSE
B(1,L) = XCOORD(1,NODEX(IEL,I))
B(2,L) = XCOORD(2,NODEX(IEL,I))
S(L) = P(NODEX(IEL,I))
END IF
END DO
CALL PLTSAI ( DS, NSTEP, PPMIN, B, 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 ) 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)
IF ( K .GT. 0 ) THEN
CALL XMOVE ( X0, Y0 )
ELSE
CALL XDRAW ( X0, Y0 )
END IF
K = -K
END IF
END DO
END DO
RETURN
END
C
C
SUBROUTINE ISOPARA ( ND, E1 , E2 , F )
IMPLICIT REAL*8 ( A-H , O-Z )
DIMENSION F(ND)
F(1) = -0.25D0*(1.D0- E1 )*(1.D0- E2 )*(E1+E2+1.D0)
F(2) = 0.50D0*(1.D0- E1*E1)*(1.D0- E2 )
F(3) = 0.25D0*(1.D0+ E1 )*(1.D0- E2 )*(E1-E2-1.D0)
F(4) = 0.50D0*(1.D0+ E1 )*(1.D0- E2*E2)
F(5) = 0.25D0*(1.D0+ E1 )*(1.D0+ E2 )*(E1+E2-1.D0)
F(6) = 0.50D0*(1.D0- E1*E1)*(1.D0+ E2 )
F(7) = -0.25D0*(1.D0- E1 )*(1.D0+ E2 )*(E1-E2+1.D0)
F(8) = 0.50D0*(1.D0- E1 )*(1.D0- E2*E2)
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======================== GRAPHICS RUTINES =============================
SUBROUTINE PLTLGO ( FILENAME )
IMPLICIT REAL*8 ( A-H , O-Z )
CHARACTER FILENAME*12
OPEN ( 1, FILE=FILENAME, STATUS='UNKNOWN' )
RETURN
END
C
C
SUBROUTINE PLTEXT
IMPLICIT REAL*8 ( A-H , O-Z )
CLOSE (1)
RETURN
END
C
C
SUBROUTINE JCOLOR ( I )
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