PROGRAM SETR4NSEQ1 C======================================================================= C DATA GENERATING PROGRAM FOR NSEQ8XX.FOR C PROJECT: DRIVEN CAVITY FLOW PROJECT NAME: SEE BELOW DATA PROJECT C DOMAIN: SQUARE ( 1 X 1 ) C ELEMENT: FOUR-NODED ISOPARAMETRIC ELEMENT C DOMAIN DISCRETIZATION: UNEVEN ELEMENTS WITH VERTICAL SCAN C REYNOLDS NUMBER: SPECIFIED IN PARAMETER AS REYNLD C VISCOSITY RATIO [VRATIO]: 2ND VISCOSITY / 1ST VISCOSITY C EIJI FUKUMORI DECEMBER 29, 1993 C JAN. 29, 2013 C======================================================================= IMPLICIT REAL*8 ( A-H , O-Z ) PARAMETER ( ND=4, MXE=12000, MXN=13000, MXB=1610, NF=3 ) PARAMETER ( REYNLD=400.D0, VRATIO=1.D7, HEIGHT=1.D0, UAVE=1.D0, * VISCO=HEIGHT*UAVE/REYNLD, FLMDA=VISCO*VRATIO, * DT = 0.01D0, TMAX=200.D0, MAXITE=10, ERMAX=0.001D0, C1=0.D0, * TLX = HEIGHT, TLY = HEIGHT, NEY = 20, NEX = 20, * DX = TLX / NEX, DY = TLY / NEY, * NDX=NEX+1, NDY=NEY+1, TIME = 0.D0, DTMAX=0.05D0, ITEFIX=3 ) PARAMETER ( INTEPT=2 ) C======================================================================= DIMENSION NODEX(MXE,ND), XCOORD(MXN), YCOORD(MXN), U(MXN), V(MXN), * IBNDFX(MXB),IBNDFY(MXB),BVX(MXB),BVY(MXB),IBNDS(MXB),BVS(MXB) CHARACTER FNAME(NF)*19,EXTEN(NF)*4,PROJECT*15,EXFILE*3 LOGICAL YES C======================================================================= DATA PROJECT / 'DCF4REYNLD-0400' / DATA EXTEN / '.DAT', '.BIN', '.STM' / C======================================================================= F(X) = X * X * ( 3.- 2.* X ) G(X) = X * ( (0.1-1.)* X + (2.- 0.1) ) C======================================================================= FNAME(1)=PROJECT//EXTEN(1) FNAME(2)=PROJECT//EXTEN(2) FNAME(3)=PROJECT//EXTEN(3) OPEN ( UNIT=2, FILE='NSDATA4.FLN', STATUS='UNKNOWN') WRITE (2,'(A15)') PROJECT WRITE (2,'(A19)') FNAME(1) WRITE (2,'(A19)') FNAME(2) WRITE (2,'(A19)') FNAME(3) CLOSE (2) WRITE (*,*) PROJECT WRITE (*,*) FNAME(1) WRITE (*,*) FNAME(2) WRITE (*,*) FNAME(3) C======================================================================= C C1: ADDED VISCOSITY COMPUTATION PARAMETER C C1=0. FOR NO ADDED VISCOSITY, C1=1. FOR UPWIND (OPTIMUM) ADDED C C1 OPTION IS AVALABLE IN NSEQ8DD.FOR SOLVER C NEY: NEMBER OF VERTICAL ELEMENTS (NUMBER OF NODES: NEY+1) C NEX: NEMBER OF HORIZONTAL ELEMENTS (NUMBER OF NODES: NEX+1) C MAXITE: MAXIMUM NUMBER OF ITERATIONS AT EACH TIME INCREMENT C ITERATION PROCEDURE IS REQUIRED BECAUSE OF NON-LINEAR CONVECTIVE TERM. C ERMAX: MAXIMUM DIFFERENCE BETWEEN GUESSED {U&V} AND COMPUTED {U&V} C HEIGHT: HEIGHT OF DRIVEN CAVITY FLOW DOMAIN C======================================================================= WRITE (*,*)' DRIVEN CAVITY FLOW: REYNOLDS NUMBER OF ',REYNLD WRITE (*,*)' VISCOSITY RATIO OF ',VRATIO C======================================================================= C ELEMENT CREATION NE = 0 DO J = 1 , NEX DO I = 1 , NEY NE = NE + 1 IF ( NE .GT. MXE ) STOP 'NE > MXE' NODEX(NE,1) = NDY*(J-1) + I NODEX(NE,2) = NODEX(NE,1) + NDY NODEX(NE,3) = NODEX(NE,2) + 1 NODEX(NE,4) = NODEX(NE,1) + 1 END DO END DO WRITE (*,*) 'NE=', NE C======================================================================= C NODAL COORDINATE CREATION NNODE = 0 DO J = 1 , NDX DO I = 1 , NDY NNODE = NNODE + 1 IF ( NNODE .GT. MXN ) STOP 'NNODE > MXN' XCOORD(NNODE) = TLX * F( DX*(J-1)/TLX ) YCOORD(NNODE) = TLY * G( DY*(I-1)/TLY ) END DO END DO WRITE (*,*) 'NNODE=', NNODE C======================================================================= C BOUNDARY CONDITIONS C--------- MOMENTUM EQUATIONS C--------- FACE OF -X NBFX = 0 NBFY = 0 U0 = 0. V0 = 0. DO J = 2 , NDY-1 NBFX = NBFX + 1 NBFY = NBFY + 1 IBNDFX(NBFX) = J IBNDFY(NBFY) = J BVX(NBFX) = U0 BVY(NBFY) = V0 END DO C--------- FACE OF -Y U0 = 0. V0 = 0. DO I = 1 , NDX NBFX = NBFX + 1 IBNDFX(NBFX) = 1 + NDY*(I-1) NBFY = NBFY + 1 IBNDFY(NBFY) = 1 + NDY*(I-1) BVX(NBFX) = U0 BVY(NBFY) = V0 END DO C--------- FACE OF +Y U0 = UAVE V0 = 0. DO I = 1 , NDX NBFX = NBFX + 1 IBNDFX(NBFX) = I*NDY NBFY = NBFY + 1 IBNDFY(NBFY) = I*NDY BVX(NBFX) = U0 BVY(NBFY) = V0 END DO C--------- FACE OF +X U0 = 0. V0 = 0. DO I = 2 , NEY NBFX = NBFX + 1 IBNDFX(NBFX) = NEX*NDY + I NBFY = NBFY + 1 IBNDFY(NBFY) = NEX*NDY + I BVX(NBFX) = U0 BVY(NBFY) = V0 END DO WRITE (*,*) 'NBFX=', NBFX WRITE (*,*) 'NBFY=', NBFY C======================================================================= C INITIAL CONDITIONS U0 = 0. V0 = 0. DO I = 1 , NNODE U(I) = U0 V(I) = V0 END DO C======================================================================= C B. C. FOR STREAM FUNCTION DDY = YCOORD(NDY) - YCOORD(NDY-1) Q = DDY * UAVE / 2. NBS = NBFX DO I = 1 , NBS IBNDS(I) = IBNDFX(I) BVS(I) = BVX(I) IF ( BVX(I) .EQ. UAVE ) BVS(I) = Q END DO C======================================================================= C MAKING DATA FILES C---------- HEADER OPEN ( 1, FILE=FNAME(1), STATUS = 'UNKNOWN' ) WRITE(1,*) VISCO, FLMDA WRITE(1,*) TMAX , DTMAX WRITE(1,*) MAXITE, ERMAX, ITEFIX C---------- ELEMENT WRITE(1,*) NE DO I = 1 , NE WRITE (1,*) I, (NODEX(I,J), J = 1 , ND ) END DO C---------- NODE WRITE(1,*) NNODE DO I = 1 , NNODE WRITE(1,*) I, XCOORD(I), YCOORD(I) END DO C---------- BOUNDARY 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 CLOSE (1) C======================================================================= C---------- 'PROJECT'.BIN OPEN ( 1, FILE=FNAME(2),STATUS ='UNKNOWN',FORM='UNFORMATTED' ) WRITE(1) TIME, DT WRITE(1) ( U(I) , I = 1 , NNODE ) WRITE(1) ( V(I) , I = 1 , NNODE ) CLOSE (1) C======================================================================= C---------- 'PROJECT'.STM OPEN ( 1, FILE=FNAME(3), STATUS = 'UNKNOWN' ) WRITE(1,*) NBS WRITE (1,*) ( IBNDS(I) , BVS(I) , I = 1 , NBS ) CLOSE (1) C======================================================================= CALL BANDWID ( ND, MXE, NE, NODEX , NBW ) C------ CREATION OF PARAMETER FILE TO BE USED IN INCLUDE STATEMENT OPEN ( 1, FILE='PARAM.DAT', STATUS='UNKNOWN' ) WRITE (1,*) ' PARAMETER ( NF=',NF,' )' 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,')' NB = MAX0 (NBFX, NBFY, NBS ) WRITE (1,*) ' PARAMETER ( MXB=',NB,' )' CLOSE (1) C======================================================================= STOP "NORMAL TERMINATION" END C C SUBROUTINE BANDWID ( ND, MXE, NE, NODEX , NBW ) IMPLICIT REAL*8 ( A-H , O-Z ) 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(*,*) 'HALF BANDWIDTH + 1 =', NBW RETURN END