⭐ 欢迎来到虫虫下载站! | 📦 资源下载 📁 资源专辑 ℹ️ 关于我们
⭐ 虫虫下载站

📄 coneq.f

📁 油田化学驱模拟的经典 fortran源代码
💻 F
📖 第 1 页 / 共 5 页
字号:
                        COE1(I) = 0.5*(VX(LZ(I,L),L)+VX(LZ(I,L)-1,L))
                        CG1(I) = 2.*
     *                        (C(LZ(I,L)+1,KC,L)-C(LZ(I,L)-1,KC,L))/
     *                  (DDX(LZ(I,L))*(ONE+ALX1(LZ(I,L)))+DDX(LZ(I,L)
     *                  -1)* ( ONE + ALX1(LZ(I,L)-1) ) )
                        IF (S(LZ(I,L)+1,L).LE.ONEM8.AND.S(LZ(I,L)-1,L)
     *                     .GT.ONEM8) THEN 
                           CG1(I) = 2.*(C(LZ(I,L),KC,L)-C(LZ(I,L)-1,KC
     *                      ,L))/(DDX(LZ(I,L)-1)*(ONE+ALX1(LZ(I,L)-1)))
                           LWKSP1(I) = .FALSE.
                        ENDIF
  420                CONTINUE
                     DO 425 I = 2,NBLM3
                        IF (LWKSP1(I).AND.S(LZ(I,L)-1,L).LE.ONEM8
     *                     .AND.S(LZ(I,L)+1,L).GT.ONEM8) THEN 
                           CG1(I) = 2.*(C(LZ(I,L)+1,KC,L)-C(LZ(I,L),KC
     *                         ,L))/(DDX(LZ(I,L))*(ONE+ALX1(LZ(I,L))))
                        ENDIF
  425                CONTINUE
                     DO 430 I = 1,NBLM3,NX
                        COE1(I) = 0.5*VX(LZ(I,L),L)
                        CG1(I) = 2.*(C(LZ(I,L)+1,KC,L)-C(LZ(I,L),KC,L))
     *                           /(DDX(LZ(I,L))*(ONE+ALX1(LZ(I,L))))
                        IF (S(LZ(I,L),L).LE.ONEM8.OR.S(LZ(I,L)+1,L).LE.
     *                         ONEM8) CG1(I)=0.0
  430                CONTINUE
                     DO 435 I = NX,NBLM3,NX
                        COE1(I) = 0.5*VX(LZ(I,L)-1,L)
                        CG1(I) = 2.*(C(LZ(I,L),KC,L)-C(LZ(I,L)-1,KC,L))
     *                          /(DDX(LZ(I,L)-1)*(ONE+ALX1(LZ(I,L)-1)))
                        IF (S(LZ(I,L),L).LE.ONEM8.OR.S(LZ(I,L)-1,L).LE.
     *                         ONEM8) CG1(I)=0.0
  435                CONTINUE
                  ELSE
                     DO 423 I = 2,NBLM3
                        LWKSP1(I) = .TRUE.
                        COE1(I) = 0.5*(VX(LZ(I,L),L)+VX(LZ(I,L)-1,L))
                        CG1(I) = (C(LZ(I,L)+1,KC,L)-C(LZ(I,L)-1,KC,L))/
     *                (.5*(DDX(LZ(I,L)+1)+DDX(LZ(I,L)-1))+DDX(LZ(I,L)))
                        IF (S(LZ(I,L)+1,L).LE.ONEM8.AND.S(LZ(I,L)-1,L) 
     *                     .GT.ONEM8) THEN 
                           CG1(I) = 2.*(C(LZ(I,L),KC,L)-C(LZ(I,L)-1,KC
     *                         ,L))/(DDX(LZ(I,L))+DDX(LZ(I,L)-1)) 
                           LWKSP1(I) = .FALSE.
                        ENDIF
  423                CONTINUE
                     DO 428 I = 2,NBLM3
                        IF (LWKSP1(I).AND.S(LZ(I,L)-1,L).LE.ONEM8
     *                    .AND.S(LZ(I,L)+1,L).GT.ONEM8) THEN 
                           CG1(I) = 2.*(C(LZ(I,L)+1,KC,L)-C(LZ(I,L),KC
     *                         ,L))/(DDX(LZ(I,L))+DDX(LZ(I,L)+1))
                        ENDIF
  428                CONTINUE
                     DO 433 I = 1,NBLM3,NX
                        COE1(I) = 0.5*VX(LZ(I,L),L)
                        CG1(I) = 2.*(C(LZ(I,L)+1,KC,L)-C(LZ(I,L),KC,L))
     *                         /(DDX(LZ(I,L))+DDX(LZ(I,L)+1))
                        IF (S(LZ(I,L),L).LE.ONEM8.OR.S(LZ(I,L)+1,L).LE.
     *                   ONEM8) CG1(I)=0.0
  433                CONTINUE
                     DO 437 I = NX,NBLM3,NX
                        COE1(I) = 0.5*VX(LZ(I,L)-1,L)
                        CG1(I) = 2.*(C(LZ(I,L),KC,L)-C(LZ(I,L)-1,KC,L))
     *                         /(DDX(LZ(I,L))+DDX(LZ(I,L)-1))
                        IF (S(LZ(I,L),L).LE.ONEM8.OR.S(LZ(I,L)-1,L).LE.
     *                  ONEM8) CG1(I)=0.0  
  437                CONTINUE
                  ENDIF
               ENDIF
C
C     FLUX AND CONCENTRATION GRADIENT FOR THE CALCULATION
C     OF Z-Y CROSS DISPERSION 
C
               IF (NY.GT.1) THEN
               DO 440 I = NX+1,NBLM3
                  LWKSP1(I) = .TRUE.
                  COE2(I) = 0.5*(VY(LZ(I,L),L)+VY(LZ(I,L)-NX,L))
C
                  CG2(I)=(C(LZ(I,L)+NX,KC,L)-C(LZ(I,L)-NX,KC,L))/
     *              (.5*(DDY(LZ(I,L)+NX)+DDY(LZ(I,L)-NX))+DDY(LZ(I,L)))
                  IF (S(LZ(I,L)+NX,L).LE.ONEM8.AND.S(LZ(I,L)-NX,L).GT. 
     *              ONEM8) THEN 
                    CG2(I) = 2.*(C(LZ(I,L),KC,L)-C(LZ(I,L)-NX,KC,L))/
     *                     (DDY(LZ(I,L))+DDY(LZ(I,L)-NX))
                    LWKSP1(I) = .FALSE.
                  ENDIF
  440          CONTINUE
               DO 450 I = NX+1,NBLM3
                  IF (LWKSP1(I).AND.S(LZ(I,L)-NX,L).LE.ONEM8
     *              .AND.S(LZ(I,L)+NX,L).GT.ONEM8) THEN 
                    CG2(I) = 2.*(C(LZ(I,L)+NX,KC,L)-C(LZ(I,L),KC,L))/
     *                     (DDY(LZ(I,L))+DDY(LZ(I,L)+NX))
                  ENDIF
  450          CONTINUE
                  DO 480 K = 1,NZ-1
                     IBGN = (K-1)*NXNY+1
                     IEND = IBGN+NX-1 
                     DO 460 I = IBGN,IEND
                        COE2(I) = 0.5*VY(LZ(I,L),L)
                        CG2(I) = 2.*(C(LZ(I,L)+NX,KC,L)-C(LZ(I,L),KC,L)
     *                      )/(DDY(LZ(I,L))+DDY(LZ(I,L)+NX))
                        IF (S(LZ(I,L),L).LE.ONEM8.OR.S(LZ(I,L)+NX,L).LE.
     *                     ONEM8) CG2(I)=0.0 
  460                CONTINUE                            
                     IBGN = K*NXNY-NX+1
                     IEND = IBGN+NX-1
                     DO 470 I = IBGN,IEND
                        COE2(I) = 0.5*VY(LZ(I,L)-NX,L)
                        CG2(I)=2.*(C(LZ(I,L),KC,L)-C(LZ(I,L)-NX,KC,L))/
     *                         (DDY(LZ(I,L))+DDY(LZ(I,L)-NX))
                        IF (S(LZ(I,L),L).LE.ONEM8.OR.S(LZ(I,L)-NX,L)
     *                     .LE.ONEM8) CG2(I) = 0.0 
  470                CONTINUE 
  480             CONTINUE
               ELSE 
                  DO 490 I = 1,NBLM3
                     COE2(I) = 0.0
  490             CONTINUE
               ENDIF
C
C     CALCULATE DISPERSION DIAGONAL TERMS
C 
               DO 500 I = 1,NBLM3
                  COE3(I) = VZ(I,L)**2
                  IF (ABS(VZ(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)                                
  500          CONTINUE
               DO 510 I = 1,NBLM3
                  IF (KC.GE.MTWM1+1.AND.KC.LE.MTWM1+NTW) THEN  
                     COE7(I) = SF(LZ(I,L),L)*PORC(LZ(I,L))
                  ELSE
                     COE7(I) = S(LZ(I,L),L)*PORC(LZ(I,L))
                  ENDIF
  510          CONTINUE
C
C     ADD THE MOLECULAR DIFFUSION
C                                  
               DO 520 I = 1,NBLM3
                  COE7(I) = COE7(I)*D(KC,L)+(ALPHAL(L)*COE3(I)+ALPHAT(L)
     *                      *(COE4(I)+COE5(I)))*COE6(I)
                  COE1(I) = ALPHAX*(COE1(I)*VZ(I,L))*COE6(I)
                  COE2(I) = ALPHAX*(COE2(I)*VZ(I,L))*COE6(I)
  520          CONTINUE
C
C     CONC. GRADIENT IN Z-DIRECTION
C
C
               IF (ICOORD.EQ.4) THEN
C
C  ---FLEXI GRID 
C
                  DO 530 I = 1,NBLM3
                     CG3(I) = 2.*(C(I+NXNY,KC,L)-C(I,KC,L))
     *                      /(DDZ(I)*(ONE+ALZ1(I)))
                     IF (S(I+NXNY,L).LE.ONEM8.OR.S(I,L).LE.ONEM8)
     *                  CG3(I)=0.0
  530             CONTINUE          
               ELSE
                  DO 535 I = 1,NBLM3
                     CG3(I) = 2.*(C(I+NXNY,KC,L)-C(I,KC,L))/
     *                        (DDZ(I+NXNY)+DDZ(I))
                     IF (S(I+NXNY,L).LE.ONEM8.OR.S(I,L).LE.ONEM8)
     *                  CG3(I)=0.0
  535             CONTINUE          
               ENDIF    
C
C     CALCULATE NUMERICAL DISPERSION COEFFICIENT TERM
C
               DO 537 I = NXNY+1,NBLM3
                  COE6(I) = 0.0
                  IF (ITC.GE.1) COE6(I)=-.5*DT*VZ(I,L)*VZ(I,L)/PORC(I)
  537          CONTINUE
               DO 538 I = 1,NXNY
                  COE6(I) = 0.0
                  IF (ITC.GE.1) COE6(I)=-.5*DT*VZ(I,L)*VZ(I,L)/PORC(I)
  538          CONTINUE
C  
C     CALCULATE NUMERICAL DISPERSION CONTROL AND
C     DIAGONAL DISPERSION TERMS IN Z-DIRECTION
C
               DO 540 I = 1,NBLM3
                  COE3(I) = (COE6(I)-COE7(I))*CG3(I) 
  540          CONTINUE
C
C     CALCULATE Z-X CROSS DISPERSION TERM
C
               IF (NX.GT.1) THEN 
                  DO 550 I = 1,NBLM3 
                     COE4(I) = -COE1(I)*CG1(I)
                     IF (ITC.GE.1) COE4(I) = COE4(I)-.5*DT*VZ(I,L)*
     *                  VX(I,L)/PORC(I)*CG1(I)
  550             CONTINUE
               ELSE
                  DO 560 I = 1,NBLM3
                     COE4(I) = 0.0
  560             CONTINUE
               ENDIF
C  
C     CALCULATE Z-Y CROSS DISPERSION TERM
C
               DO 570 I = 1,NBLM3
                  COE5(I) = 0.0
  570          CONTINUE
               IF (NY.GT.1) THEN 
                  DO 580 I = 1,NBLM3 
                     COE5(I) = -COE2(I)*CG2(I)
                     IF (ITC.GE.1) COE5(I) = COE5(I)-.5*DT*VY(I,L)*
     *                 VZ(I,L)/PORC(I)*CG2(I)
  580             CONTINUE
               ENDIF
C                             
C     Z-DISPERSION FOR OIL COMPONENT 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
                  DO 585 I = 1,NBLM3
                     ODT(I+NXNY,3) = COE3(I)+COE4(I)+COE5(I)
  585             CONTINUE
               ENDIF
C                             
C     Z-DISPERSION ACCUMULATE OVER PHASES
C
               DO 590 I = 1,NBLM3
                  COEZT(I+NXNY) = COEZT(I+NXNY)+COE3(I)+COE4(I)+COE5(I)
  590          CONTINUE
            ENDIF  
  600    CONTINUE
C
         DO 610 I = 1,NBL,NX
            COEX(I) = 0.0
            ODT(I,1) = 0.0
  610    CONTINUE
         DO 620 K = 1,NZ
            IBGN = (K-1)*NXNY+1
            IEND = IBGN+NXM1
            DO 620 I = IBGN,IEND
               COEYU(I) = 0.0
               ODT(I,2) = 0.0 
  620    CONTINUE
         DO 630 I = 1,NXNY
            COEZT(I) = 0.0
            ODT(I,3) = 0.0
  630    CONTINUE
C
C     CONVECTION TERM FOR OIL COMPONENT IN WATER, OR SINGLE M.E.
C     OR ME (TYPE II)
         IF (IMASS.EQ.2.AND.(KC.EQ.2.OR.
     &           (KC.GT.MMOM1.AND.KC.LE.MMOM1+NO))) THEN
            DO 631 I = 1,NBL
               CONVO(I,1) = CONVX(I,KC,1)
               CONVO(I,2) = CONVY(I,KC,1)
               CONVO(I,3) = CONVZ(I,KC,1)
 631        CONTINUE
         ENDIF
C
C     CONVECTION ACCUMULATE OVER PHASES
C
         DO 632 I = 1,NBL
            CONVX(I,KC,1) = CONVX(I,KC,1)+CONVX(I,KC,2)+CONVX(I,KC,3) 
            CONVY(I,KC,1) = CONVY(I,KC,1)+CONVY(I,KC,2)+CONVY(I,KC,3) 
            CONVZ(I,KC,1) = CONVZ(I,KC,1)+CONVZ(I,KC,2)+CONVZ(I,KC,3) 
  632    CONTINUE
C
         IF (IGAS.GE.1) THEN
            DO 633 I = 1,NBL
               CONVX(I,KC,1) = CONVX(I,KC,1)+CONVX(I,KC,NPHAS)
               CONVY(I,KC,1) = CONVY(I,KC,1)+CONVY(I,KC,NPHAS)
               CONVZ(I,KC,1) = CONVZ(I,KC,1)+CONVZ(I,KC,NPHAS)
  633       CONTINUE
         ENDIF
C
         DO 634 I = 1,NBL
            COE4(I) = 0.0
            COE5(I) = 0.0
  634    CONTINUE
C
C     OIL FLUX DIFFERENCE
         IF (IMASS.EQ.2.AND.(KC.EQ.2.OR.
     &           (KC.GT.MMOM1.AND.KC.LE.MMOM1+NO))) THEN
            IF (ICOORD.EQ.4) THEN
               DO 635 I = 1,NBL-1
                  ODT(I+1,1) = ODT(I+1,1)*DUM5(I)
  635          CONTINUE
               DO 636 I = 1,NBL-NX
                  ODT(I+NX,2) = ODT(I+NX,2)*DUM6(I)
  636          CONTINUE
               DO 637 I = 1,NBL-NXNY
                  ODT(I+NXNY,3) = ODT(I+NXNY,3)*DUM7(I)
  637          CONTINUE
               DO 638 I = 1,NBL
                  CONVO(I,1) = CONVO(I,1)*DUM5(I)
                  CONVO(I,2) = CONVO(I,2)*DUM6(I)
                  CONVO(I,3) = CONVO(I,3)*DUM7(I)
  638          CONTINUE
C
C     CALC' TRUE SURFACE FLUX DIFFERENCES
C     WHERE COEX = DISPERSION :   CONVX = CONVECTION
C
               COE3(1) = ODT(2,1)+CONVO(1,1)
               DO 639 I = 2,NBL-1
                  COE3(I) = ODT(I+1,1)-ODT(I,1)+CONVO(I,1)-CONVO(I-1,1)
  639          CONTINUE
               COE3(NBL) = -ODT(NBL,1)+CONVO(NBL,1)-CONVO(NBL-1,1)
               IF (NY.GT.1) THEN
                  DO 640 I = 1,NX
                     COE4(I) = ODT(I+NX,2)+CONVO(I,2)
  640             CONTINUE
                  DO 641 I = NX+1,NBL-NX
                     COE4(I) = ODT(I+NX,2)-ODT(I,2)+CONVO(I,2)-
     *                        CONVO(I-NX,2)
  641             CONTINUE
                  DO 642 I = NBL-NX+1,NBL
                     COE4(I) = -ODT(I,2)+CONVO(I,2)-CONVO(I-NX,2)
  642             CONTINUE
               ENDIF
               IF (NZ.GT.1) THEN
                  DO 643 I = 1,NXNY
                     COE5(I) = ODT(I+NXNY,3)+CONVO(I,3)
  643             CONTINUE
                  IF (NZ.EQ.2) THEN
                     DO 644 I = NXNY+1,NBL
                        COE5(I) = -ODT(I,3)+CONVO(I,3)-CONVO(I-NXNY,3)  
  644                CONTINUE
                  ELSE
                     DO 645 I = NXNY+1,NBLM3
                        COE5(I) = ODT(I+NXNY,3)-ODT(I,3)+CONVO(I,3)-
     *                          CONVO(I-NXNY,3)
  645                CONTINUE
                     DO 646 I = NBLM3+1,NBL
                        COE5(I) = -ODT(I,3)+CONVO(I,3)-CONVO(I-NXNY,3)
  646                CONTINUE
                  ENDIF
               ENDIF
               DO 647 I = 1,NBL
                  DCO(I) = -DT*(COE3(I)+COE4(I)+COE5(I))/VB(I)
  647          CONTINUE
            ELSE
               IF (ICOORD.NE.2) THEN

⌨️ 快捷键说明

复制代码 Ctrl + C
搜索代码 Ctrl + F
全屏模式 F11
切换主题 Ctrl + Shift + D
显示快捷键 ?
增大字号 Ctrl + =
减小字号 Ctrl + -