PROGRAM MAKING_TEST_DATA27 C======================================================================= C 27-NODED HEXA PARABOLIC ELEMENT C THREE DIMENSIONAL DATA GENERATING PROGRAM FOR LANCZOS-SOUNDX.FOR C PROJECT NAME: EIGENVALUE C BOUNDARY CONDITIONS C ALL BOUNDARY: NATURAL C DOMAIN: RECTANGULAR C ELEMENT: 27-NODED ISOPARAMETRIC HEXA 3D PARABOLIC ELEMENT C EIJI FUKUMORI DECEMBER 9, 2010 C======================================================================= IMPLICIT REAL*8 ( A-H , O-Z ) PARAMETER ( ND=27, INTEPT=3 ) PARAMETER ( MXE=220000, MXN=4000000, MXB=200000 ) C======================================================================= DIMENSION NODEX(MXE,ND), XCOORD(MXN), YCOORD(MXN) ,ZCOORD(MXN) CHARACTER INPFILE*13,EXFILE*3 CHARACTER INFNAME*12 LOGICAL YES, RSTRAINT C======================================================================= DATA INFNAME / 'EIGN3D27.DAT' / C======================================================================= C------------------------ SPEED OF SOUND ------------------------------- WSPD = 340.D0 C----------------------------------------------------------------------- C TLX=LENGTH, TLY=HEIGHT, TLZ=DEPTH TLX = 4.D0*DATAN(1.D0) TLY = TLX/2.D0 TLZ = TLX/3.D0 C----------------------------------------------------------------------- NELEMENT = 10 ALMIN = DMIN1 ( TLX, TLY, TLZ ) NEX = NELEMENT * (TLX/ALMIN + 0.001D0) NEY = NELEMENT * (TLY/ALMIN + 0.001D0) NEZ = NELEMENT * (TLZ/ALMIN + 0.001D0) WRITE (*,*) NEX, NEY, NEZ C======================================================================= DX=TLX/(2*NEX) DY=TLY/(2*NEY) DZ=TLZ/(2*NEZ) NDX = 2*NEX+1 NDY = 2*NEY+1 NDZ = 2*NEZ+1 C======================================================================= C NEZ: NEMBER OF VERTICAL ELEMENTS C NEY: NEMBER OF DEPTH ELEMENTS C NEX: NEMBER OF HORIZONTAL ELEMENTS C======================================================================= C--------- EIGRNVALUE PARAMETER MXNEIGEN=MIN0 ( 20, NNODE/10 ) IF ( MXNEIGEN .LT. 9 ) MXNEIGEN = 9 C----------------------------------------------------------------------- C--------- DELTA = SHIFT PARAMETER THAT MUST BE > 0. AL = DMAX1 ( TLX, TLY, TLZ ) DELTA = (4.D0*DATAN(1.D0)/AL)**2 * 0.9 C======================================================================= C ELEMENT CREATION NE = 0 DO I = 1 , NEX DO J = 1 , NEZ DO K = 1 , NEY NE = NE + 1 IF ( NE .GT. MXE ) STOP 'NE > MXE' C------- BOTTOM FLOOR NODEX(NE,1) = (I-1)*NDY*NDZ*2 + (J-1)*NDY*2 + (K-1)*2 + 1 NODEX(NE,2) = NODEX(NE,1) + NDY*NDZ NODEX(NE,3) = NODEX(NE,2) + NDY*NDZ NODEX(NE,4) = NODEX(NE,3) + 1 NODEX(NE,5) = NODEX(NE,4) + 1 NODEX(NE,8) = NODEX(NE,1) + 1 NODEX(NE,7) = NODEX(NE,8) + 1 NODEX(NE,9) = NODEX(NE,2) + 1 NODEX(NE,6) = NODEX(NE,9) + 1 C------- MIDIUM FLOOR NODEX(NE,10) = NODEX(NE,1) + NDY NODEX(NE,17) = NODEX(NE,10) + 1 NODEX(NE,16) = NODEX(NE,17) + 1 NODEX(NE,11) = NODEX(NE,10) + NDY*NDZ NODEX(NE,18) = NODEX(NE,11) + 1 NODEX(NE,15) = NODEX(NE,18) + 1 NODEX(NE,12) = NODEX(NE,11) + NDY*NDZ NODEX(NE,13) = NODEX(NE,12) + 1 NODEX(NE,14) = NODEX(NE,13) + 1 C------- TOP FLOOR NODEX(NE,19) = NODEX(NE,10) + NDY NODEX(NE,26) = NODEX(NE,19) + 1 NODEX(NE,25) = NODEX(NE,26) + 1 NODEX(NE,20) = NODEX(NE,19) + NDY*NDZ NODEX(NE,27) = NODEX(NE,20) + 1 NODEX(NE,24) = NODEX(NE,27) + 1 NODEX(NE,21) = NODEX(NE,20) + NDY*NDZ NODEX(NE,22) = NODEX(NE,21) + 1 NODEX(NE,23) = NODEX(NE,22) + 1 END DO END DO END DO WRITE (*,*) 'NUMBER OF ELEMENTS =', NE C======================================================================= C--------0---------0---------0---------0---------0---------0---------0-- C NODAL COORDINATE CREATION NNODE = 0 DO I = 1 , NDX DO K = 1 , NDZ DO J = 1 , NDY NNODE = NNODE + 1 IF ( NNODE .GT. MXN ) STOP 'NNODE > MXN' XCOORD(NNODE) = (I-1)*DX YCOORD(NNODE) = (J-1)*DY ZCOORD(NNODE) = (K-1)*DZ END DO END DO END DO WRITE (*,*) 'NUMBER OF NODAL POINTS =', NNODE C======================================================================= IWW = 2 OPEN (IWW, FILE='DOMAIN.SUM', STATUS='UNKNOWN') WRITE (IWW,*) 'SPEED OF SOUND =', WSPD WRITE (IWW,*) 'TLX, TLY, TLZ =', TLX, TLY, TLZ WRITE (IWW,*) 'TOTAL VOLUME =', TOTALVOL WRITE (IWW,*) 'NEX, NEY, NEZ =', NEX, NEY, NEZ WRITE (IWW,*) 'MAXIMUM SURFACE LENGTH =', AL WRITE (IWW,*) 'SHIFT PARAMETER (DELTA) =', DELTA WRITE (IWW,*) 'MAX EIGEN (ALPHA SQ) =', EIGENMX WRITE (IWW,*) 'NE =', NE WRITE (IWW,*) 'NNODE =', NNODE WRITE (IWW,*) 'MXNEIGEN =', MXNEIGEN CLOSE (IWW) C======================================================================= C CREATING DATA FILES IW = 1 OPEN (IW, FILE=INFNAME, STATUS='UNKNOWN') WRITE (IW,*) DELTA, MXNEIGEN, WSPD C---------- ELEMENTS WRITE(IW,*) NE DO I = 1 , NE WRITE (IW,*) I, (NODEX(I,J), J = 1 , ND ) END DO C---------- COORDINATES OF NONAL POINTS WRITE(IW,*) NNODE DO I = 1 , NNODE WRITE(IW,*) I,XCOORD(I), YCOORD(I), ZCOORD(I) END DO C======================================================================= C---------- FINAL CLOSE (IW) C======================================================================= C------ CREATION OF PARAMETER FILE TO BE USED IN INCLUDE STATEMENT CALL BANDWID ( ND, MXE, NE, NODEX , NBW ) OPEN ( 1, FILE='PARAM.DAT', STATUS='UNKNOWN' ) WRITE (1,*) ' PARAMETER ( ND=',ND,' )' WRITE (1,*) ' PARAMETER ( INTEPT=',INTEPT,' )' WRITE (1,*) ' PARAMETER ( MXE=',NE,' )' WRITE (1,*) ' PARAMETER ( MXN=',NNODE,' )' WRITE (1,*) ' PARAMETER ( MXW=',NBW,' )' WRITE (1,*) ' PARAMETER ( MXENGN=',MXNEIGEN,')' CLOSE (1) C======================================================================= STOP 'NORMAL TERMINATION' END C C SUBROUTINE BANDWID ( ND, MXE, NE, NODEX , NBW ) DIMENSION NODEX(MXE,ND) 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