📄 coneq.f
字号:
RRSQ(IK) = R(I)**2
RRP(IK) = RP(I)
15 CONTINUE
ENDIF
C
IF (IMES.NE.1) THEN
DCMAX = 0.0
DCMIN = ONEP50
ENDIF
C
IF (NG.NE.0) THEN
XOH = ONEM7
DO 16 I = 1,NBL
IF (X16.LE.ZERO) THEN
IF (KGOPT.NE.3) THEN
LWKSP2(I) = C(I,NGP3,1).LT.EPS.OR.C(I,4,1).LT.EPS1
ELSE
LWKSP2(I) = C(I,NGP3,1).LT.XOH.OR.C(I,4,1).LT.EPS1
ENDIF
ELSE
LWKSP2(I) = C(I,NGP3,1).LT.EPS.OR.C(I,4,1).LT.EPS1.OR.
* HION(I).LT.EPS2
ENDIF
C MALONATE CONC. IS BELOW EPS1 OR NOT
LWKSP3(I) = C(I,NGP1,1).LT.EPS1
LWKSP4(I) = C(I,NGP2,1).LT.EPS1
16 CONTINUE
ENDIF
C
C CALCULATE THE POROSITY FOR THE PREVIOUS TIME STEP AND CORRECTION
C FACTORS DUE TO PORE COMPRESSIBILITIES
C
DO 17 I = 1,NBL
PORC(I) = POR(I)*(1.+COMPR*(POLD(I)-PR(I)))
17 CONTINUE
IF (IMASS.EQ.2) THEN
DO 19 J = 1,3
DO 19 I = 1,NBL
ODT(I,J) = 0.0
19 CONTINUE
ENDIF
C
C SET ICF(2) = 0 TO SKIP CALCULATION
C
ICFSV2 = ICF(2)
IF (LMO) THEN
ICF(2) = 0
ENDIF
C
DO 900 KCR = 1,N
KC = KCR
IF (ICF(KC).EQ.0) GO TO 890
EPHI1 = 1.0
IF (KC.EQ.3) EPHI1 = EPHI3
IF (KC.EQ.4) EPHI1 = EPHI4
DO 20 I = 1,NBL
COEX(I) = 0.0
COEYU(I) = 0.0
COEZT(I) = 0.0
ODT(I,1) = 0.0
ODT(I,2) = 0.0
ODT(I,3) = 0.0
20 CONTINUE
C
C CALCULATE THE CORRECTION FACTOR DUE FLUID COMPRESSIBILITIES
C
DO 26 I = 1,NBL
CORFF(I) = 1.+COMPC(KC)*(P(I,1)-PSTAND)
26 CONTINUE
C
C CALCULATE TOTAL FLUXES IN AND OUT OF THE GRID-BLOCKS
C
DO 600 L = 1,NPHAS
ALPHAX = ALPHAL(L)-ALPHAT(L)
DO 25 I = 1,NBL
LWKSP1(I) = .TRUE.
25 CONTINUE
C
C X-DIRECTION DERIVATIVES
C
IF (NY.GT.1) THEN
C
C FLUX AND CONCENTRATION GRADIENT FOR THE CALCULATION
C OF X-Y CROSS DISPERSION
C
DO 30 I = NX+1,NBLM2
COE1(I) = 0.5*(VY(LX(I,L),L)+VY(LX(I,L)-NX,L))
CG1(I) = (C(LX(I,L)+NX,KC,L)-C(LX(I,L)-NX,KC,L))/
* (.5*(DDY(LX(I,L)+NX)+DDY(LX(I,L)-NX))+DDY(LX(I,L)))
IF (S(LX(I,L)+NX,L).LE.ONEM8.AND.S(LX(I,L)-NX,L).GT.
* ONEM8) THEN
CG1(I) = 2.*(C(LX(I,L),KC,L)-C(LX(I,L)-NX,KC,L))/
* (DDY(LX(I,L))+DDY(LX(I,L)-NX))
LWKSP1(I) = .FALSE.
ENDIF
30 CONTINUE
DO 35 I = NX+1,NBLM2
IF (LWKSP1(I).AND.S(LX(I,L)-NX,L).LE.ONEM8
* .AND.S(LX(I,L)+NX,L).GT.ONEM8) THEN
CG1(I) = 2.*(C(LX(I,L)+NX,KC,L)-C(LX(I,L),KC,L))/
* (DDY(LX(I,L))+DDY(LX(I,L)+NX))
ENDIF
35 CONTINUE
DO 60 K = 1,NZ
IBGN = (K-1)*NXNY+1
IEND = IBGN+NXM1
DO 40 I = IBGN,IEND
COE1(I) = 0.5*VY(LX(I,L),L)
CG1(I) = 2.*(C(LX(I,L)+NX,KC,L)-C(LX(I,L),KC,L))/
* (DDY(LX(I,L))+DDY(LX(I,L)+NX))
IF (S(LX(I,L),L).LE.ONEM8.OR.S(LX(I,L)+NX,L).LE.
* ONEM8) CG1(I)=0.0
40 CONTINUE
IBGN = K*NXNY-NXM1
IEND = IBGN+NXM1
DO 50 I = IBGN,IEND
COE1(I) = 0.5*VY(LX(I,L)-NX,L)
CG1(I) = 2.*(C(LX(I,L),KC,L)-C(LX(I,L)-NX,KC,L))/
* (DDY(LX(I,L))+DDY(LX(I,L)-NX))
IF (S(LX(I,L),L).LE.ONEM8.OR.S(LX(I,L)-NX,L).LE.
* ONEM8) CG1(I)=0.0
50 CONTINUE
60 CONTINUE
ELSE
DO 70 I = 1,NBL
COE1(I) = 0.0
70 CONTINUE
ENDIF
C
C FLUX AND CONCENTRATION GRADIENT FOR THE CALCULATION
C OF X-Z CROSS DISPERSION
C
IF (NZ.GT.1) THEN
IF (ICOORD.EQ.4) THEN
C
C FLEXI GRID: NOTE INTODUCTION OF RATIO ALZ1
C BOUNDARY CONDITIONS: DCDZ = 0 | Z = TOP
C
DO 80 I = 1,NXNY
COE2(I) = 0.5*VZ(LX(I,L),L)
CG2(I) = 2.*(C(LX(I,L)+NXNY,KC,L)-C(LX(I,L),KC,L))
* /(DDZ(LX(I,L))*(ONE+ALZ1(LX(I,L))))
IF (S(LX(I,L),L).LE.ONEM8.OR.S(LX(I,L)+NXNY,L).LE.
* ONEM8) CG2(I)=0.0
80 CONTINUE
C
C BOUNDARY CONDITIONS: DCDZ = 0 | Z = BOTTOM
C
IBGN = (NZ-1)*NXNY+1
DO 90 I = IBGN,NBLM1
COE2(I) = 0.5*VZ(LX(I,L)-NXNY,L)
CG2(I) = 2.*(C(LX(I,L),KC,L)-C(LX(I,L)-NXNY,KC,L))
* /(DDZ(LX(I,L)-NXNY)*(ONE+ALZ1(LX(I,L)-NXNY)))
IF (S(LX(I,L),L).LE.ONEM8.OR.S(LX(I,L)-NXNY,L).LE.
* ONEM8) CG2(I)=0.0
90 CONTINUE
C
C DCDZ IN FIELD
C
IF (NZ.GT.2) THEN
DO 100 I = NXNY+1,NBLM3
LWKSP1(I) = .TRUE.
COE2(I) = 0.5*(VZ(LX(I,L),L)+VZ(LX(I,L)-NXNY,L))
CG2(I)=2.*
* (C(LX(I,L)+NXNY,KC,L)-C(LX(I,L)-NXNY,KC,L))
* /(DDZ(LX(I,L))*(ONE+ALZ1(LX(I,L)))+
* DDZ(LX(I,L)-NXNY)*(ONE+ALZ1(LX(I,L)-NXNY)))
IF (S(LX(I,L)+NXNY,L).LE.ONEM8.AND.
* S(LX(I,L)-NXNY,L).GT.ONEM8) THEN
CG2(I) = 2.*(C(LX(I,L),KC,L)-
* C(LX(I,L)-NXNY,KC,L))/(DDZ(LX(I,L)-NXNY)*
* (ONE+ALZ1(LX(I,L)-NXNY)))
LWKSP1(I) = .FALSE.
ENDIF
100 CONTINUE
DO 105 I = NXNY+1,NBLM3
IF (LWKSP1(I).AND.S(LX(I,L)-NXNY,L).LE.ONEM8
* .AND.S(LX(I,L)+NXNY,L).GT.ONEM8) THEN
CG2(I) = 2.*(C(LX(I,L)+NXNY,KC,L)-C(LX(I,L)
* ,KC,L))/(DDZ(LX(I,L))*(ONE+ALZ1(LX(I,L))))
ENDIF
105 CONTINUE
ENDIF
ELSE
C
C FLUX AND CONCENTRATION GRADIENT FOR THE CALCULATION
C OF X-Z CROSS DISPERSION
C
DO 85 I = 1,NXNY
COE2(I) = 0.5*VZ(LX(I,L),L)
CG2(I) = 2.*(C(LX(I,L)+NXNY,KC,L)-C(LX(I,L),KC,L
* ))/(DDZ(LX(I,L))+DDZ(LX(I,L)+NXNY))
IF (S(LX(I,L),L).LE.ONEM8.OR.S(LX(I,L)+NXNY,L).LE.
* ONEM8) CG2(I)=0.0
85 CONTINUE
IBGN = (NZ-1)*NXNY+1
DO 95 I = IBGN,NBLM1
COE2(I) = 0.5*VZ(LX(I,L)-NXNY,L)
CG2(I) = 2.*(C(LX(I,L),KC,L)-C(LX(I,L)-NXNY,KC,L))/
* (DDZ(LX(I,L))+DDZ(LX(I,L)-NXNY))
IF (S(LX(I,L),L).LE.ONEM8.OR.S(LX(I,L)-NXNY,L).LE.
* ONEM8) CG2(I)=0.0
95 CONTINUE
IF (NZ.GT.2) THEN
DO 103 I = NXNY+1,NBLM3
LWKSP1(I) = .TRUE.
COE2(I) = 0.5*(VZ(LX(I,L),L)+VZ(LX(I,L)-NXNY,L))
CG2(I) = (C(LX(I,L)+NXNY,KC,L)-C(LX(I,L)-NXNY,
* KC,L))/(.5*(DDZ(LX(I,L)+NXNY)+
* DDZ(LX(I,L)-NXNY))+DDZ(LX(I,L)))
IF (S(LX(I,L)+NXNY,L).LE.ONEM8.AND.
* S(LX(I,L)-NXNY,L).GT.ONEM8) THEN
CG2(I) = 2.*(C(LX(I,L),KC,L)-C(LX(I,L)-
* NXNY,KC,L))/(DDZ(LX(I,L))+DDZ(LX(I,L)-NXNY))
LWKSP1(I) = .FALSE.
ENDIF
103 CONTINUE
DO 106 I = NXNY+1,NBLM3
IF (LWKSP1(I).AND.S(LX(I,L)-NXNY,L).LE.ONEM8
* .AND.S(LX(I,L)+NXNY,L).GT.ONEM8) THEN
CG2(I) = 2.*(C(LX(I,L)+NXNY,KC,L)-C(LX(I,L)
* ,KC,L))/(DDZ(LX(I,L))+DDZ(LX(I,L)+NXNY))
ENDIF
106 CONTINUE
ENDIF
ENDIF
ELSE
DO 110 I = 1,NBLM1
COE2(I) = 0.0
110 CONTINUE
ENDIF
C
C CALCULATE DISPERSION DIAGONAL TERMS
C
DO 120 I = 1,NBLM1
COE3(I) = VX(I,L)**2
IF (ABS(VX(I,L)).LE.ONEM10) THEN
COE1(I) = 0.0
COE2(I) = 0.0
ENDIF
COE4(I) = COE1(I)**2
COE5(I) = COE2(I)**2
COE6(I) = SQRT(COE3(I)+COE4(I)+COE5(I))
COE6(I) = CVMGT(ONE199,COE6(I),COE6(I).LE.ONEM10)
COE6(I) = ONE/COE6(I)
120 CONTINUE
DO 130 I = 1,NBLM1
IF (KC.GE.MTWM1+1.AND.KC.LE.MTWM1+NTW) THEN
COE7(I) = SF(LX(I,L),L)*PORC(LX(I,L))
ELSE
COE7(I) = S(LX(I,L),L)*PORC(LX(I,L))
ENDIF
130 CONTINUE
C
C ADD THE MOLECULAR DIFFUSION
C
DO 140 I = 1,NBLM1
COE7(I) = COE7(I)*D(KC,L)+(ALPHAL(L)*COE3(I)+ALPHAT(L)*
* (COE4(I)+COE5(I)))*COE6(I)
COE1(I) = ALPHAX*(COE1(I)*VX(I,L))*COE6(I)
COE2(I) = ALPHAX*(COE2(I)*VX(I,L))*COE6(I)
140 CONTINUE
C
C CONC. GRADIENT IN X-DIRECTION
C
IF (ICOORD.NE.2) THEN
C
IF (ICOORD.EQ.4) THEN
C
C ---FLEXI GRID
C
DO 143 I = 1,NBLM1
CG3(I) = 2.*(C(I+1,KC,L)-C(I,KC,L))
* /(DDX(I)*(ONE+ALX1(I)))
IF (S(I+1,L).LE.ONEM8.OR.S(I,L).LE.ONEM8)
* CG3(I)=0.0
143 CONTINUE
ELSE
DO 145 I = 1,NBLM1
CG3(I) = 2.*(C(I+1,KC,L)-C(I,KC,L))/(DDX(I+1)
* +DDX(I))
IF (S(I+1,L).LE.ONEM8.OR.S(I,L).LE.ONEM8)
* CG3(I)=0.0
145 CONTINUE
ENDIF
ELSE
DO 147 I = 1,NBL-1
CG3(I) = 2.*RRP(I)*(C(I+1,KC,L)-C(I,KC,L))/
* (RRSQ(I+1)-RRSQ(I))
IF (S(I+1,L).LE.ONEM8.OR.S(I,L).LE.ONEM8)CG3(I)=0.0
147 CONTINUE
ENDIF
C
C CALCULATE NUMERICAL DISPERSION COEFFICIENT TERM
C
DO 148 I = 2,NBLM1
COE6(I) = 0.0
C
C SECOND ORDER TIME CORRECTION
C
IF (ITC.GE.1) COE6(I) = -.5*DT*VX(I,L)*VX(I,L)/
* PORC(I)
148 CONTINUE
COE6(1) = 0.0
IF (ITC.GE.1) COE6(1) = -.5*DT*VX(1,L)*VX(1,L)/
* PORC(1)
C
C CALCULATE NUMERICAL DISPERSION CONTROL AND
C DIAGONAL DISPERSION TERMS IN X-DIRECTION
C
DO 150 I = 1,NBLM1
COE3(I) = (COE6(I)-COE7(I))*CG3(I)
150 CONTINUE
C
C CALCULATE X-Y CROSS DISPERSION TERM
C
IF (NY.GT.1) THEN
DO 160 I = 1,NBLM1
COE4(I) = -COE1(I)*CG1(I)
IF (ITC.GE.1) COE4(I) = COE4(I)-.5*DT*VX(I,L)*VY(I,L)
* /PORC(I)*CG1(I)
160 CONTINUE
ELSE
DO 200 I = 1,NBLM1
COE4(I) = 0.0
200 CONTINUE
ENDIF
C
C CALCLATE X-Z CROSS DISPERSION TERM
C
DO 210 I = 1,NBLM1
COE5(I) = 0.0
210 CONTINUE
IF (NZ.GT.1) THEN
DO 220 I = 1,NBLM1
COE5(I) = -COE2(I)*CG2(I)
IF (ITC.GE.1) COE5(I) = COE5(I)-.5*DT*VX(I,L)*VZ(I,L)/
* PORC(I)*CG2(I)
220 CONTINUE
ENDIF
C
C X-DISPERSION FOR COMPONENT 2 IN WATER PHASE
IF (IMASS.EQ.2.AND.(KC.EQ.2.OR.(KC.GT.MMOM1.AND.
* KC.LE.MMOM1+NO)).AND.L.EQ.1) THEN
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -