Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- SUBROUTINE VORTX1
- c Elaborado por Serguei A. Lonin, 1999-2000 Version 2.1, 2000
- INCLUDE 'COMML.FOR'
- epsp=0.0001
- MOL=0
- DO 19 I=1,M1
- DO 19 J=1,N1
- IF(KT.EQ.1)PSI(J,I)=PSI0(J,I)
- IF(H(J,I).EQ.PF)GOTO 19
- MOL=MOL+1
- 19 CONTINUE
- TETA=1.
- L=0
- 35 L=L+1
- IF(L.GE.2000000)GOTO 325
- CALL BOUND
- DO 778 I=2,M
- DO 778 J=2,N
- IF(H(J,I).EQ.PF)GOTO 778
- C IF(H(J+1,I).EQ.PF.OR.H(J,I+1).EQ.PF.OR.H(J,I-1).EQ.PF
- C &.OR.H(J-1,I).EQ.PF)GOTO 778
- C IF(H(J+1,I+1).EQ.PF.OR.H(J-1,I+1).EQ.PF.OR.H(J+1,I-1).EQ.PF
- C &.OR.H(J-1,I-1).EQ.PF)GOTO 778
- DDT=(AFV(J,I)+AFV(J-1,I)-AFV(J,I-1)-AFV(J-1,I-1)
- &-AFU(J,I)-AFU(J,I-1)+AFU(J-1,I)+AFU(J-1,I-1))
- A1=RK(J,I)*HL1(J,I)/HL2(J,I)/H1(J,I)
- A4=RK(J,I-1)*HL1(J,I-1)/HL2(J,I-1)/H1(J,I-1)
- A8=RK(J-1,I)*HL1(J-1,I)/HL2(J-1,I)/H1(J-1,I)
- A9=RK(J-1,I-1)*HL1(J-1,I-1)/HL2(J-1,I-1)/H1(J-1,I-1)
- B4=YL*HH1(J+1,I)/HH2(J+1,I)/H(J+1,I)
- B8=YL*HH1(J-1,I)/HH2(J-1,I)/H(J-1,I)
- G1=RK(J,I)*HL2(J,I)/HL1(J,I)/H1(J,I)
- G2=RK(J-1,I)*HL2(J-1,I)/HL1(J-1,I)/H1(J-1,I)
- G6=RK(J,I-1)*HL2(J,I-1)/HL1(J,I-1)/H1(J,I-1)
- G7=RK(J-1,I-1)*HL2(J-1,I-1)/HL1(J-1,I-1)/H1(J-1,I-1)
- O2=YL*HH2(J,I+1)/HH1(J,I+1)/H(J,I+1)
- O6=YL*HH2(J,I-1)/HH1(J,I-1)/H(J,I-1)
- PSI1=((G2+G1)*PSI0(J,I+1)+(G7+G6)*PSI(J,I-1)+
- & (A4+A1)*PSI0(J+1,I)+(A9+A8)*PSI(J-1,I)+
- & 0.5*(B4*(PSI0(J+1,I+1)-PSI0(J+1,I-1))-B8*(PSI0(J-1,I+1)-
- & PSI(J-1,I-1))-O2*(PSI0(J+1,I+1)-PSI0(J-1,I+1))+
- & O6*(PSI0(J+1,I-1)-PSI(J-1,I-1)))-HAG*DDT)/ ! DOESN'T DIVIDE BY 2 AND ONE HAG!
- & (G1+G2+G6+G7+A1+A4+A8+A9)
- PSI(J,I)=TETA*PSI1+(1.-TETA)*PSI0(J,I)
- 778 CONTINUE
- CALL VELCTY
- CALL FRIC
- CALL RHS
- if(L.eq.1)go to 35
- KAP=0
- DO 97 I=1,M1
- DO 97 J=1,N1
- IF(H(J,I).EQ.PF)GOTO 97
- IF(PSI0(J,I).EQ.0..AND.PSI(J,I).EQ.0.)GOTO 322
- IF(ABS((PSI(J,I)-PSI0(J,I))/(PSI0(J,I)+0.1E-09)).GT.EPSP)GOTO 97
- 322 KAP=KAP+1
- 97 CONTINUE
- DO 2 I=1,M1
- DO 2 J=1,N1
- IF(H(J,I).EQ.PF)GOTO 2
- PSI0(J,I)=PSI(J,I)
- 2 CONTINUE
- IF(KAP.EQ.MOL)GOTO 444
- IF(L/ITER1*ITER1.EQ.L)WRITE(*,*)
- &'Acercamiento esta en ',KAP/FLOAT(MOL)*100,' %'
- IF(KAP/FLOAT(MOL)*100..GE.PERC)GO TO 444
- GOTO 35
- 444 CONTINUE
- CCC write(*,*)l,kap,mol
- RETURN
- 325 WRITE(*,764)
- RETURN
- END
Advertisement
Add Comment
Please, Sign In to add comment