PROGRAM SETAXISYMM4 C======================================================================= C ********** 4-NODED ELEMENT ********** C------------------- ELEMENT NUMBERING: HORIZONTAL SCAN ---------------- C DATA GENERATING PROGRAM FOR STATIC8QFXAXISYMM-NEW.FOR C DOMAIN SIZE: TLX BY TLY C BOUNDARY CONDITIONS C------------ NORMAL DIRECTION C ON FACE OF -X: P=P0 C ON FACE OF +X: P=0 C ON FACE OF -Y: V=0 C ON FACE OF +Y: V=0 C------------ TANGENTIAL DIRECTION C ON FACE OF -X: NATURAL BOUNDARY C ON FACE OF +X: NATURAL BOUNDARY C ON FACE OF -Y: NATURAL BOUNDARY C ON FACE OF +Y: NATURAL BOUNDARY C NOV. 28, 2012 C EIJI FUKUMORI C======================================================================= IMPLICIT REAL*8 ( A-H , O-Z ) PARAMETER ( ND=4, INTEPT=2, MXE=31300, MXN=34410, MXB=21000 ) PARAMETER ( TLX0=0.1D0,TLX1=0.15D0,TLY=0.2D0 ) PARAMETER ( YOUNG=210.D9, POISSON=0.3D0 ) PARAMETER ( NEX=11,NEY=22 ) PARAMETER ( P0=10.D6, ALPHA = 1.D0 ) C======================================================================= DIMENSION NODEX(MXE,ND), XCOORD(MXN), YCOORD(MXN), * IBNDFX(MXB),IBNDFY(MXB),BVX(MXB),BVY(MXB), * IBFORCEX(MXB), IBFORCEY(MXB), BVFORCEX(MXB), BVFORCEY(MXB) DIMENSION NEUTRAL(MXN) CHARACTER PROJECT*12,EXFILE*3 LOGICAL YES C======================================================================= DATA PROJECT / 'STATIC4.DAT' / C======================================================================= C NEX: NEMBER OF ELEMENTS IN RADIAL DIRECTION C NEY: NEMBER OF ELEMENTS IN VERTICAL DIRECTION C P0, P1, P2: APPLIED HYDRAULIC PRESSURE C======================================================================= THICKNSS = TLX1-TLX0 DX = THICKNSS / NEX DY = TLY / NEY NDX = NEX + 1 NDY = NEY + 1 C======================================================================= WRITE (*,*)' YOUNG MODULUS = ',YOUNG WRITE (*,*)' POISSON RATIO = ',POISSON C======================================================================= C ELEMENT CREATION NE = 0 DO I = 1 , NEY DO J = 1 , NEX NE = NE + 1 IF ( NE .GT. MXE ) STOP 'NE > MXE' NODEX(NE,1) = (I-1)*NDX + J NODEX(NE,2) = NODEX(NE,1) + 1 NODEX(NE,3) = NODEX(NE,2) + NDX NODEX(NE,4) = NODEX(NE,1) + NDX END DO END DO C======================================================================= C NODAL COORDINATE CREATION NNODE = 0 DO I = 1 , NDY DO J = 1 , NDX NNODE = NNODE + 1 IF ( NNODE .GT. MXN ) STOP 'NNODE > MXN' XCOORD(NNODE) = THICKNSS * F(DX*(J-1)/THICKNSS,ALPHA) + TLX0 YCOORD(NNODE) = DY*(I-1) END DO END DO C======================================================================= PI = 4.D0* ATAN (1.D0) C BOUNDARY CONDITIONS C==== FIRST KIND C--------- NAVIER EQUATIONS NBFX = 0 NBFY = 0 NFORCEX = 0 NFORCEY = 0 C C------- BOTTOM DO I = 1 , NDX NBFY = NBFY + 1 IBNDFY(NBFY) = I BVY(NBFY) = 0.D0 END DO C------- TOP DO I = 1 , NDX NBFY = NBFY + 1 IBNDFY(NBFY) = I + NEY*NDX BVY(NBFY) = 0.D0 END DO C C==== SECOND KIND C--------- FACE OF INNER SIDE DF0 = 2.D0*PI*TLX0*DY*P0 C------ RESET MEMORY FOR NEUMANN TYPE BOUNDARY VALUE DO I = 1 , NDY BVFORCEX(I) = 0.D0 END DO C--------- LOAD DISTRIBUTION FOR PARABOLIC ELEMENT UNDER UNIFORM LOAD WEIGHT = 1.D0/2.D0 DO I = 1 , NEY NODE = (I-1)*NDX + 1 NFORCEX = (I-1) + 1 IBFORCEX(NFORCEX) = NODE BVFORCEX(NFORCEX) = BVFORCEX(NFORCEX) + WEIGHT*DF0 NODE = NODE + NDX NFORCEX = I + 1 IBFORCEX(NFORCEX) = NODE BVFORCEX(NFORCEX) = BVFORCEX(NFORCEX) +WEIGHT*DF0 END DO C======================================================================= WRITE (*,*) ' NUMBER OF ELEMENTS (NE) = ',NE WRITE (*,*) ' NUMBER OF NODAL POINTS (NNODE) = ',NNODE WRITE (*,*) ' NUMBER OF DIRICHLET R (NBFX) = ',NBFX WRITE (*,*) ' NUMBER OF DIRICHLET Z (NBFY) = ',NBFY WRITE (*,*) ' NUMBER OF NEUMANN R (NFORCEX) = ',NFORCEX WRITE (*,*) ' NUMBER OF NEUMANN Z (NFORCEY) = ',NFORCEY C======================================================================= C DATA FILE INQUIRY EXFILE = 'NEW' INQUIRE ( FILE = PROJECT, EXIST = YES ) IF ( YES ) EXFILE='OLD' C======================================================================= C MAKING DATA FILES C---------- 'PROJECT'.JNK IR = 1 C---------- PARAMETERS OPEN ( IR, FILE=PROJECT, STATUS = EXFILE ) WRITE(1,*) YOUNG , POISSON C---------- ELEMENTS C---------- ELEMENTS WRITE(1,*) NE DO I = 1 , NE WRITE (1,*) I, (NODEX(I,J), J = 1 , ND ) END DO C---------- COORDINATES OF NONAL POINTS WRITE(1,*) NNODE DO I = 1 , NNODE WRITE(1,*) I,XCOORD(I), YCOORD(I) END DO C---------- DIRICHLET TYPE BOUNDARY CONDITIONS WRITE(1,*) NBFX DO I = 1 , NBFX WRITE (1,*) IBNDFX(I), BVX(I) END DO WRITE(1,*) NBFY DO I = 1 , NBFY WRITE (1,*) IBNDFY(I), BVY(I) END DO C---------- NUEMANN TYPE BOUNDARY CONDITIONS WRITE(1,*) NFORCEX IF ( NFORCEX .GT. 0 ) THEN DO I = 1 , NFORCEX WRITE(1,*) IBFORCEX(I), BVFORCEX(I) END DO ENDIF WRITE(1,*) NFORCEY IF ( NFORCEY .GT. 0 ) THEN DO I = 1 , NFORCEY WRITE(1,*) IBFORCEY(I), BVFORCEY(I) END DO ENDIF C---------- FINAL CLOSE (1) C---------- ELEMENT DRAWING OPEN ( 1, FILE='ELEMENT4.DAT', STATUS = 'UNKNOWN') DO I = 1, NE DO J = 1, ND WRITE (1,*) XCOORD(NODEX(I,J)), YCOORD(NODEX(I,J)) END DO WRITE (1,*) XCOORD(NODEX(I,1)), YCOORD(NODEX(I,1)) WRITE (1,*) END DO CLOSE (1) C======================================================================= C------ CREATION OF PARAMETER FILE TO BE USED IN INCLUDE STATEMENT CALL BANDWID ( MXE, ND, NE, NODEX, NBW ) OPEN ( 1, FILE='PARAM.DAT', STATUS='UNKNOWN' ) WRITE (1,*) ' PARAMETER ( ND=',ND,', INTEPT=',INTEPT,' )' WRITE (1,*) ' PARAMETER ( MXE=',NE,', MXN=',NNODE, * ', MXB=',MAX0(NBFX,NBFY,NFORCEX,NFORCEY), ', MXW=',NBW, ' )' CLOSE (1) STOP 'NORMAL TERMINATION' END C C FUNCTION F (X, A) IMPLICIT REAL*8 ( A-H , O-Z ) F = X**A RETURN END C C SUBROUTINE BANDWID ( MXE , ND , NE , NODEX , NBW ) DIMENSION NODEX(MXE,ND) C------- RETURN VALUE: NBW NBW = 0 DO I = 1 , NE DO J = 1 , ND-1 DO K = J+1 , ND NBW = MAX0 ( NBW , IABS(NODEX(I,J)-NODEX(I,K)) ) END DO END DO END DO NBW = NBW + 1 WRITE (*,*) ' HALH BANDWIDTH =', NBW RETURN END