home *** CD-ROM | disk | FTP | other *** search
/ Hall of Fame / HallofFameCDROM.cdr / oilfield / spe-46-2.lzh / BLOCK7.FOR < prev    next >
Text File  |  1988-06-21  |  18KB  |  489 lines
  1. $DO66
  2. C................................................................SOLONE
  3.       SUBROUTINE SOLONE(II,JJ,KK,DIV1,D288,
  4.      &KSM,KSM1,N,NN,KCOFF)
  5. C      MACHINE DEPENDENT INCLUDE STATEMENT
  6. $INCLUDE:'PARAMS.FOR'
  7. C      FLOW EQ COEF WITH SINGLE POINT UPSTREAM REL PERMS
  8.       REAL KROT,KROGT,KRWT,KRGT,MUOT,MUWT,MUGT,KX,KY,KZ
  9.       REAL MUO4,MUW4,MUG4,MUO5,MUW5,MUG5,MUO6,MUW6,MUG6
  10.      &,KRO1,KRW1,KRG1,KRO2,KRW2,KRG2,KRO3,KRW3,KRG3
  11.      &,MUO1,MUW1,MUG1,MUO2,MUW2,MUG2,MUO3,MUW3,MUG3
  12.      &,KRO4,KRW4,KRG4,KRO5,KRW5,KRG5,KRO6,KRW6,KRG6
  13.      &,MO1,MW1,MG1,MO2,MW2,MG2,MO3,MW3,MG3
  14.      &,MO4,MW4,MG4,MO5,MW5,MG5,MO6,MW6,MG6
  15.      &,MUO,MUW,MUG
  16.       COMMON /BUBBLE/ PBO,VSLOPE(LP8),BSLOPE(LP8),RSLOPE(LP8),PMAXT,
  17.      & IREPRS,MPGT(LP8),
  18.      & RHOSCO(LP8),RHOSCW(LP8),RHOSCG(LP8),MSAT(LP7),MPOT(LP8),
  19.      & MPWT(LP8),PBOT(LP1,LP2,LP3),PBOTN(LP1,LP2,LP3)
  20.       COMMON /COEF/ AW(LP1,LP2,LP3),AE(LP1,LP2,LP3),AN(LP1,LP2,LP3),
  21.      & AS(LP1,LP2,LP3),AB(LP1,LP2,LP3),AT(LP1,LP2,LP3),E(LP1,LP2,LP3),
  22.      & B(LP1,LP2,LP3)
  23.       COMMON /SARRAY/ PN(LP1,LP2,LP3),IOCODE,IDMAX,
  24.      & SON(LP1,LP2,LP3),SWN(LP1,LP2,LP3),SGN(LP1,LP2,LP3),
  25.      & A1(LP1,LP2,LP3),A2(LP1,LP2,LP3),A3(LP1,LP2,LP3),
  26.      & SUM(LP1,LP2,LP3),GAM(LP1,LP2,LP3),QS(LP1,LP2,LP3)
  27.       COMMON /SPARM/ KX(LP1,LP2,LP3),KY(LP1,LP2,LP3),KZ(LP1,LP2,LP3),
  28.      & EL(LP1,LP2,LP3),TX(LP4,LP2,LP3),TY(LP1,LP5,LP3),TZ(LP1,LP2,LP6),
  29.      & PDAT(LP1,LP2,LP3),PDATUM,GRAD
  30.       COMMON /SPRTPS/ P(LP1,LP2,LP3),SO(LP1,LP2,LP3),SW(LP1,LP2,LP3),
  31.      & SG(LP1,LP2,LP3)
  32.       COMMON /SPVT/ SAT(LP7,LP9),KROT(LP7,LP9),KRWT(LP7,LP9),
  33.      & BGT(LP7,LP9),
  34.      & KRGT(LP7,LP9),ITHREE(LP7),RSOT(LP7,LP9),BWPT(LP7,LP9),
  35.      & PCOWT(LP7,LP9),PCGOT(LP7,LP9),KROGT(LP7,LP9),SWR(LP7),
  36.      & POT(LP7,LP9),MUOT(LP7,LP9),BOT(LP7,LP9),BOPT(LP7,LP9),
  37.      & RSOPT(LP7,LP9),PWT(LP7,LP9),MUWT(LP7,LP9),BWT(LP7,LP9),
  38.      & RSWT(LP7,LP9),RSWPT(LP7,LP9),PGT(LP7,LP9),MUGT(LP7,LP9),
  39.      & BGPT(LP7,LP9),CRT(LP7,LP9),IPVT(LP1,LP2,LP3),IROCK(LP1,LP2,LP3),
  40.      & NROCK,NPVT,PSIT(LP7,LP9),PRT(LP7,LP9),WOROCK(LP7),GOROCK(LP7)
  41.       COMMON /SSOLN/ BO(LP1,LP2,LP3),BW(LP1,LP2,LP3),BG(LP1,LP2,LP3),
  42.      & QO(LP1,LP2,LP3),QW(LP1,LP2,LP3),QG(LP1,LP2,LP3),
  43.      & GOWT(LP1,LP2,LP3),GWWT(LP1,LP2,LP3),GGWT(LP1,LP2,LP3),
  44.      & OW(LP4,LP2,LP3),OE(LP4,LP2,LP3),WW(LP4,LP2,LP3),WE(LP4,LP2,LP3),
  45.      & OS(LP1,LP5,LP3),ON(LP1,LP5,LP3),WS(LP1,LP5,LP3),WN(LP1,LP5,LP3),
  46.      & OT(LP1,LP2,LP6),OB(LP1,LP2,LP6),WT(LP1,LP2,LP6),WB(LP1,LP2,LP6),
  47.      & QOWG(LP1,LP2,LP3),VP(LP1,LP2,LP3),CT(LP1,LP2,LP3)
  48.       COMMON /VECTOR/ DX(LP1,LP2,LP3),DY(LP1,LP2,LP3),DZ(LP1,LP2,LP3),
  49.      & DZNET(LP1,LP2,LP3),IQN1(LP11),IQN2(LP11),IQN3(LP11),IHEDIN(80)
  50.       DIMENSION RPW(LP1,LP2,LP3),RPG(LP1,LP2,LP3),RPOW(LP1,LP2,LP3),
  51.      & RPO3(LP1,LP2,LP3),CAPOW(LP1,LP2,LP3),CAPGO(LP1,LP2,LP3)
  52.       EQUIVALENCE (RPW,DX), (RPG,DY), (RPOW,DZ), (RPO3,KX),
  53.      & (CAPOW,KY), (CAPGO,KZ)
  54.       DATA RSO1,RSO2,RSO3,RSO4,RSO5,RSO6/6*0.0/
  55.       DATA RSW1,RSW2,RSW3,RSW4,RSW5,RSW6/6*0.0/
  56.       DO 2000 K=1,KK
  57.       DO 2000 J=1,JJ
  58.       DO 2000 I=1,II
  59.       RPW(I,J,K)=0.
  60.       RPG(I,J,K)=0.
  61.       RPOW(I,J,K)=0.
  62.       RPO3(I,J,K)=0.
  63.       CAPOW(I,J,K)=0.
  64.       CAPGO(I,J,K)=0.
  65.       IF(VP(I,J,K).LE.0.0) GO TO 2000
  66.       SSO=SO(I,J,K)
  67.       SSW=SW(I,J,K)
  68.       SSG=SG(I,J,K)
  69.       IROCKB=IROCK(I,J,K)
  70.       CALL INTERP(IROCKB,SAT,KROT,MSAT(IROCKB),SSO,KRO1)
  71.       RPOW(I,J,K)=KRO1
  72.       CALL INTERP(IROCKB,SAT,KRWT,MSAT(IROCKB),SSW,KRW1)
  73.       RPW(I,J,K)=KRW1
  74.       CALL INTERP(IROCKB,SAT,KRGT,MSAT(IROCKB),SSG,KRG1)
  75.       RPG(I,J,K)=KRG1
  76.       CALL TRIKRO(IROCKB,SSO,SSW,KRO3)
  77.       RPO3(I,J,K)=KRO3
  78.       CALL INTERP(IROCKB,SAT,PCOWT,MSAT(IROCKB),SSW,PCOW)
  79.       CAPOW(I,J,K)=PCOW
  80.       CALL INTERP(IROCKB,SAT,PCGOT,MSAT(IROCKB),SSG,PCGO)
  81.       CAPGO(I,J,K)=PCGO
  82.  2000 CONTINUE
  83.       DO 200 K=1,KK
  84.       DO 200 J=1,JJ
  85.       DO 200 I=1,II
  86.       IF(VP(I,J,K).LE.0.0) GO TO 200
  87. C      I,J,K BLOCK
  88.       PP=P(I,J,K)
  89.       BPT=PBOT(I,J,K)
  90.       IPVTR=IPVT(I,J,K)
  91.       CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),PP,RSO)
  92.       CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),PP,MUO)
  93.       CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),PP,RSW)
  94.       CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),PP,MUW)
  95.       CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),PP,MUG)
  96.       SSO=SO(I,J,K)
  97.       SSW=SW(I,J,K)
  98.       SSG=SG(I,J,K)
  99.       IROCKB=IROCK(I,J,K)
  100.       PCOW=CAPOW(I,J,K)
  101.       PCGO=CAPGO(I,J,K)
  102.       RO=(RHOSCO(IPVTR) + RSO*RHOSCG(IPVTR))/BO(I,J,K)
  103.       RW=(RHOSCW(IPVTR) + RSW*RHOSCG(IPVTR))/BW(I,J,K)
  104.       RG=RHOSCG(IPVTR)/BG(I,J,K)
  105.       IF(I.EQ.1)GO TO 115
  106.       IF(VP(I-1,J,K).LE.0.0) GO TO 115
  107. C      I-1,J,K BLOCK
  108.       P1=P(I-1,J,K)
  109.       BPT=PBOT(I-1,J,K)
  110.       IPVTR=IPVT(I-1,J,K)
  111.       CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P1,RSO1)
  112.       CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P1,MUO1)
  113.       CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P1,RSW1)
  114.       CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P1,MUW1)
  115.       CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P1,MUG1)
  116.       SO1S=SO(I-1,J,K)
  117.       SW1S=SW(I-1,J,K)
  118.       SG1S=SG(I-1,J,K)
  119.       IROCKR=IROCK(I-1,J,K)
  120.       PCOW1=CAPOW(I-1,J,K)
  121.       PCGO1=CAPGO(I-1,J,K)
  122.       RO1=(RHOSCO(IPVTR) + RSO1*RHOSCG(IPVTR))/BO(I-1,J,K)
  123.       RW1=(RHOSCW(IPVTR) + RSW1*RHOSCG(IPVTR))/BW(I-1,J,K)
  124.       RG1=RHOSCG(IPVTR)/BG(I-1,J,K)
  125.       FACT=-D288*(EL(I-1,J,K)-EL(I,J,K))
  126.       GOW1=(RO1+RO)*FACT
  127.       GWW1=(RW1+RW)*FACT + PCOW-PCOW1
  128.       GGW1=(RG1+RG)*FACT + PCGO1-PCGO
  129.       P11=P1-PP
  130.       HO1=P11+GOW1
  131.       HW1=P11+GWW1
  132.       HG1=P11+GGW1
  133.       IF(ITHREE(IROCKB).EQ.1) GO TO 10
  134.       KRO1=RPOW(I,J,K)
  135.       IF(HO1.GE.0.) KRO1=RPOW(I-1,J,K)
  136.       GO TO 15
  137.    10 KRO1=RPO3(I,J,K)
  138.       IF(HO1.GE.0.) KRO1=RPO3(I-1,J,K)
  139.    15 CONTINUE
  140.       KRW1=RPW(I,J,K)
  141.       IF(HW1.GE.0.) KRW1=RPW(I-1,J,K)
  142.       KRG1=RPG(I,J,K)
  143.       IF(HG1.GE.0.) KRG1=RPG(I-1,J,K)
  144.       MO1=4.0*KRO1/((BO(I-1,J,K)+BO(I,J,K)) * (MUO1+MUO))
  145.       MW1=4.0*KRW1/((BW(I-1,J,K)+BW(I,J,K)) * (MUW1+MUW))
  146.       MG1=4.0*KRG1/((BG(I-1,J,K)+BG(I,J,K)) * (MUG1+MUG))
  147. 115   AOW=TX(I,J,K)*MO1
  148.       AWW=TX(I,J,K)*MW1
  149.       AGW=TX(I,J,K)*MG1
  150.       IF(I.EQ.II)GO TO 125
  151.       IF(VP(I+1,J,K).LE.0.0) GO TO 125
  152. C      I+1,J,K BLOCK
  153.       P2=P(I+1,J,K)
  154.       BPT=PBOT(I+1,J,K)
  155.       IPVTR=IPVT(I+1,J,K)
  156.       CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P2,RSO2)
  157.       CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P2,MUO2)
  158.       CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P2,RSW2)
  159.       CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P2,MUW2)
  160.       CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P2,MUG2)
  161.       SO2=SO(I+1,J,K)
  162.       SW2=SW(I+1,J,K)
  163.       SG2=SG(I+1,J,K)
  164.       IROCKR=IROCK(I+1,J,K)
  165.       PCOW2=CAPOW(I+1,J,K)
  166.       PCGO2=CAPGO(I+1,J,K)
  167.       RO2=(RHOSCO(IPVTR) + RSO2*RHOSCG(IPVTR))/BO(I+1,J,K)
  168.       RW2=(RHOSCW(IPVTR) + RSW2*RHOSCG(IPVTR))/BW(I+1,J,K)
  169.       RG2=RHOSCG(IPVTR)/BG(I+1,J,K)
  170.       FACT=-D288*(EL(I+1,J,K)-EL(I,J,K))
  171.       GOW2=(RO2+RO)*FACT
  172.       GWW2=(RW2+RW)*FACT + PCOW-PCOW2
  173.       GGW2=(RG2+RG)*FACT + PCGO2-PCGO
  174.       P22=P2-PP
  175.       HO2=P22+GOW2
  176.       HW2=P22+GWW2
  177.       HG2=P22+GGW2
  178.       IF(ITHREE(IROCKB).EQ.1) GO TO 20
  179.       KRO2=RPOW(I,J,K)
  180.       IF(HO2.GE.0.) KRO2=RPOW(I+1,J,K)
  181.       GO TO 25
  182.    20 KRO2=RPO3(I,J,K)
  183.       IF(HO2.GE.0.) KRO2=RPO3(I+1,J,K)
  184.    25 CONTINUE
  185.       KRW2=RPW(I,J,K)
  186.       IF(HW2.GE.0.) KRW2=RPW(I+1,J,K)
  187.       KRG2=RPG(I,J,K)
  188.       IF(HG2.GE.0.) KRG2=RPG(I+1,J,K)
  189.       MO2=4.0*KRO2/((BO(I+1,J,K)+BO(I,J,K)) * (MUO2+MUO))
  190.       MW2=4.0*KRW2/((BW(I+1,J,K)+BW(I,J,K)) * (MUW2+MUW))
  191.       MG2=4.0*KRG2/((BG(I+1,J,K)+BG(I,J,K)) * (MUG2+MUG))
  192. 125   AOE=TX(I+1,J,K)*MO2
  193.       AWE=TX(I+1,J,K)*MW2
  194.       AGE=TX(I+1,J,K)*MG2
  195.       IF(J.EQ.1)GO TO 135
  196.       IF(VP(I,J-1,K).LE.0.0) GO TO 135
  197. C      I,J-1,K BLOCK
  198.       P3=P(I,J-1,K)
  199.       BPT=PBOT(I,J-1,K)
  200.       IPVTR=IPVT(I,J-1,K)
  201.       CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P3,RSO3)
  202.       CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P3,MUO3)
  203.       CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P3,RSW3)
  204.       CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P3,MUW3)
  205.       CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P3,MUG3)
  206.       SO3=SO(I,J-1,K)
  207.       SW3=SW(I,J-1,K)
  208.       SG3=SG(I,J-1,K)
  209.       IROCKR=IROCK(I,J-1,K)
  210.       PCOW3=CAPOW(I,J-1,K)
  211.       PCGO3=CAPGO(I,J-1,K)
  212.       RO3=(RHOSCO(IPVTR) + RSO3*RHOSCG(IPVTR))/BO(I,J-1,K)
  213.       RW3=(RHOSCW(IPVTR) + RSW3*RHOSCG(IPVTR))/BW(I,J-1,K)
  214.       RG3=RHOSCG(IPVTR)/BG(I,J-1,K)
  215.       FACT=-D288*(EL(I,J-1,K)-EL(I,J,K))
  216.       GOW3=(RO3+RO)*FACT
  217.       GWW3=(RW3+RW)*FACT + PCOW-PCOW3
  218.       GGW3=(RG3+RG)*FACT + PCGO3-PCGO
  219.       P33=P3-PP
  220.       HO3=P33+GOW3
  221.       HW3=P33+GWW3
  222.       HG3=P33+GGW3
  223.       IF(ITHREE(IROCKB).EQ.1) GO TO 30
  224.       KRO3=RPOW(I,J,K)
  225.       IF(HO3.GE.0.) KRO3=RPOW(I,J-1,K)
  226.       GO TO 35
  227.    30 KRO3=RPO3(I,J,K)
  228.       IF(HO3.GE.0.) KRO3=RPO3(I,J-1,K)
  229.    35 CONTINUE
  230.       KRW3=RPW(I,J,K)
  231.       IF(HW3.GE.0.) KRW3=RPW(I,J-1,K)
  232.       KRG3=RPG(I,J,K)
  233.       IF(HG3.GE.0.) KRG3=RPG(I,J-1,K)
  234.       MO3=4.0*KRO3/((BO(I,J-1,K)+BO(I,J,K)) * (MUO3+MUO))
  235.       MW3=4.0*KRW3/((BW(I,J-1,K)+BW(I,J,K)) * (MUW3+MUW))
  236.       MG3=4.0*KRG3/((BG(I,J-1,K)+BG(I,J,K)) * (MUG3+MUG))
  237. 135   AOS=TY(I,J,K)*MO3
  238.       AWS=TY(I,J,K)*MW3
  239.       AGS=TY(I,J,K)*MG3
  240.       IF(J.EQ.JJ)GO TO 140
  241.       IF(VP(I,J+1,K).LE.0.0) GO TO 140
  242. C      I,J+1,K BLOCK
  243.       P4=P(I,J+1,K)
  244.       BPT=PBOT(I,J+1,K)
  245.       IPVTR=IPVT(I,J+1,K)
  246.       CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P4,RSO4)
  247.       CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P4,MUO4)
  248.       CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P4,RSW4)
  249.       CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P4,MUW4)
  250.       CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P4,MUG4)
  251.       SO4=SO(I,J+1,K)
  252.       SW4=SW(I,J+1,K)
  253.       SG4=SG(I,J+1,K)
  254.       IROCKR=IROCK(I,J+1,K)
  255.       PCOW4=CAPOW(I,J+1,K)
  256.       PCGO4=CAPGO(I,J+1,K)
  257.       RO4=(RHOSCO(IPVTR) + RSO4*RHOSCG(IPVTR))/BO(I,J+1,K)
  258.       RW4=(RHOSCW(IPVTR) + RSW4*RHOSCG(IPVTR))/BW(I,J+1,K)
  259.       RG4=RHOSCG(IPVTR)/BG(I,J+1,K)
  260.       FACT=-D288*(EL(I,J+1,K)-EL(I,J,K))
  261.       GOW4=(RO4+RO)*FACT
  262.       GWW4=(RW4+RW)*FACT + PCOW-PCOW4
  263.       GGW4=(RG4+RG)*FACT + PCGO4-PCGO
  264.       P44=P4-PP
  265.       HO4=P44+GOW4
  266.       HW4=P44+GWW4
  267.       HG4=P44+GGW4
  268.       IF(ITHREE(IROCKB).EQ.1) GO TO 40
  269.       KRO4=RPOW(I,J,K)
  270.       IF(HO4.GE.0.) KRO4=RPOW(I,J+1,K)
  271.       GO TO 45
  272.    40 KRO4=RPO3(I,J,K)
  273.       IF(HO4.GE.0.) KRO4=RPOW(I,J+1,K)
  274.    45 CONTINUE
  275.       KRW4=RPW(I,J,K)
  276.       IF(HW4.GE.0.) KRW4=RPW(I,J+1,K)
  277.       KRG4=RPG(I,J,K)
  278.       IF(HG4.GE.0.) KRG4=RPG(I,J+1,K)
  279.       MO4=4.0*KRO4/((BO(I,J+1,K)+BO(I,J,K)) * (MUO4+MUO))
  280.       MW4=4.0*KRW4/((BW(I,J+1,K)+BW(I,J,K)) * (MUW4+MUW))
  281.       MG4=4.0*KRG4/((BG(I,J+1,K)+BG(I,J,K)) * (MUG4+MUG))
  282. 140   AON=TY(I,J+1,K)*MO4
  283.       AWN=TY(I,J+1,K)*MW4
  284.       AGN=TY(I,J+1,K)*MG4
  285.       IF(K.EQ.1)GO TO 145
  286.       IF(VP(I,J,K-1).LE.0.0) GO TO 145
  287. C      I,J,K-1 BLOCK
  288.       P5=P(I,J,K-1)
  289.       BPT=PBOT(I,J,K-1)
  290.       IPVTR=IPVT(I,J,K-1)
  291.       CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P5,RSO5)
  292.       CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P5,MUO5)
  293.       CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P5,RSW5)
  294.       CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P5,MUW5)
  295.       CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P5,MUG5)
  296.       SO5=SO(I,J,K-1)
  297.       SW5=SW(I,J,K-1)
  298.       SG5=SG(I,J,K-1)
  299.       IROCKR=IROCK(I,J,K-1)
  300.       PCOW5=CAPOW(I,J,K-1)
  301.       PCGO5=CAPGO(I,J,K-1)
  302.       RO5=(RHOSCO(IPVTR) + RSO5*RHOSCG(IPVTR))/BO(I,J,K-1)
  303.       RW5=(RHOSCW(IPVTR) + RSW5*RHOSCG(IPVTR))/BW(I,J,K-1)
  304.       RG5=RHOSCG(IPVTR)/BG(I,J,K-1)
  305.       FACT=-D288*(EL(I,J,K-1)-EL(I,J,K))
  306.       GOW5=(RO5+RO)*FACT
  307.       GWW5=(RW5+RW)*FACT + PCOW-PCOW5
  308.       GGW5=(RG5+RG)*FACT + PCGO5-PCGO
  309.       P55=P5-PP
  310.       HO5=P55+GOW5
  311.       HW5=P55+GWW5
  312.       HG5=P55+GGW5
  313.       IF(ITHREE(IROCKB).EQ.1) GO TO 50
  314.       KRO5=RPOW(I,J,K)
  315.       IF(HO5.GE.0.) KRO5=RPOW(I,J,K-1)
  316.       GO TO 55
  317.    50 KRO5=RPO3(I,J,K)
  318.       IF(HO5.GE.0.) KRO5=RPO3(I,J,K-1)
  319.    55 CONTINUE
  320.       KRW5=RPW(I,J,K)
  321.       IF(HW5.GE.0.) KRW5=RPW(I,J,K-1)
  322.       KRG5=RPG(I,J,K)
  323.       IF(HG5.GE.0.) KRG5=RPG(I,J,K-1)
  324.       MO5=4.0*KRO5/((BO(I,J,K-1)+BO(I,J,K)) * (MUO5+MUO))
  325.       MW5=4.0*KRW5/((BW(I,J,K-1)+BW(I,J,K)) * (MUW5+MUW))
  326.       MG5=4.0*KRG5/((BG(I,J,K-1)+BG(I,J,K)) * (MUG5+MUG))
  327. 145   AOT=TZ(I,J,K)*MO5
  328.       AWT=TZ(I,J,K)*MW5
  329.       AGT=TZ(I,J,K)*MG5
  330.       IF(K.EQ.KK)GO TO 150
  331.       IF(VP(I,J,K+1).LE.0.0) GO TO 150
  332. C      I,J,K+1 BLOCK
  333.       P6=P(I,J,K+1)
  334.       BPT=PBOT(I,J,K+1)
  335.       IPVTR=IPVT(I,J,K+1)
  336.       CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P6,RSO6)
  337.       CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P6,MUO6)
  338.       CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P6,RSW6)
  339.       CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P6,MUW6)
  340.       CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P6,MUG6)
  341.       SO6=SO(I,J,K+1)
  342.       SW6=SW(I,J,K+1)
  343.       SG6=SG(I,J,K+1)
  344.       IROCKR=IROCK(I,J,K+1)
  345.       PCOW6=CAPOW(I,J,K+1)
  346.       PCGO6=CAPGO(I,J,K+1)
  347.       RO6=(RHOSCO(IPVTR) + RSO6*RHOSCG(IPVTR))/BO(I,J,K+1)
  348.       RW6=(RHOSCW(IPVTR) + RSW6*RHOSCG(IPVTR))/BW(I,J,K+1)
  349.       RG6=RHOSCG(IPVTR)/BG(I,J,K+1)
  350.       FACT=-D288*(EL(I,J,K+1)-EL(I,J,K))
  351.       GOW6=(RO6+RO)*FACT
  352.       GWW6=(RW6+RW)*FACT + PCOW-PCOW6
  353.       GGW6=(RG6+RG)*FACT + PCGO6-PCGO
  354.       P66=P6-PP
  355.       HO6=P66+GOW6
  356.       HW6=P66+GWW6
  357.       HG6=P66+GGW6
  358.       IF(ITHREE(IROCKB).EQ.1) GO TO 60
  359.       KRO6=RPOW(I,J,K)
  360.       IF(HO6.GE.0.) KRO6=RPOW(I,J,K+1)
  361.       GO TO 65
  362.    60 KRO6=RPO3(I,J,K)
  363.       IF(HO6.GE.0.) KRO6=RPO3(I,J,K+1)
  364.    65 CONTINUE
  365.       KRW6=RPW(I,J,K)
  366.       IF(HW6.GE.0.) KRW6=RPW(I,J,K+1)
  367.       KRG6=RPG(I,J,K)
  368.       IF(HG6.GE.0.) KRG6=RPG(I,J,K+1)
  369.       MO6=4.0*KRO6/((BO(I,J,K+1)+BO(I,J,K)) * (MUO6+MUO))
  370.       MW6=4.0*KRW6/((BW(I,J,K+1)+BW(I,J,K)) * (MUW6+MUW))
  371.       MG6=4.0*KRG6/((BG(I,J,K+1)+BG(I,J,K)) * (MUG6+MUG))
  372. 150   AOB=TZ(I,J,K+1)*MO6
  373.       AWB=TZ(I,J,K+1)*MW6
  374.       AGB=TZ(I,J,K+1)*MG6
  375.       RSO1A=0.5*(RSO1+RSO)
  376.       RSO2A=0.5*(RSO2+RSO)
  377.       RSO3A=0.5*(RSO3+RSO)
  378.       RSO4A=0.5*(RSO4+RSO)
  379.       RSO5A=0.5*(RSO5+RSO)
  380.       RSO6A=0.5*(RSO6+RSO)
  381.       RSW1A=0.5*(RSW1+RSW)
  382.       RSW2A=0.5*(RSW2+RSW)
  383.       RSW3A=0.5*(RSW3+RSW)
  384.       RSW4A=0.5*(RSW4+RSW)
  385.       RSW5A=0.5*(RSW5+RSW)
  386.       RSW6A=0.5*(RSW6+RSW)
  387.       AO1=AOW*GOW1
  388.       AO2=AOE*GOW2
  389.       AO3=AOS*GOW3
  390.       AO4=AON*GOW4
  391.       AO5=AOT*GOW5
  392.       AO6=AOB*GOW6
  393.       AW1=AWW*GWW1
  394.       AW2=AWE*GWW2
  395.       AW3=AWS*GWW3
  396.       AW4=AWN*GWW4
  397.       AW5=AWT*GWW5
  398.       AW6=AWB*GWW6
  399.       GOWT(I,J,K)= AO1 + AO2 + AO3 + AO4 + AO5 + AO6
  400.       GWWT(I,J,K)= AW1 + AW2 + AW3 + AW4 + AW5 + AW6
  401.       GGWT(I,J,K)=AGW*GGW1+AGE*GGW2+AGS*GGW3+AGN*GGW4+AGT*GGW5+AGB*GGW6
  402.      &+RSO1A*AO1+RSO2A*AO2+RSO3A*AO3+RSO4A*AO4+RSO5A*AO5+RSO6A*AO6
  403.      &+RSW1A*AW1+RSW2A*AW2+RSW3A*AW3+RSW4A*AW4+RSW5A*AW5+RSW6A*AW6
  404.       QOWG(I,J,K)=(BO(I,J,K)-BG(I,J,K)*RSO)*(-GOWT(I,J,K)+QO(I,J,K)) +
  405.      &   (BW(I,J,K)-BG(I,J,K)*RSW)*(-GWWT(I,J,K)+QW(I,J,K)) +
  406.      &    BG(I,J,K)*(-GGWT(I,J,K)+QG(I,J,K))
  407.       AW(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO1-RSO)) * AOW +
  408.      & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW1-RSW)) * AWW +
  409.      &  BG(I,J,K)*AGW
  410.       AE(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO2-RSO)) * AOE +
  411.      & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW2-RSW)) * AWE +
  412.      &  BG(I,J,K)*AGE
  413.       AS(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO3-RSO)) * AOS +
  414.      & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW3-RSW)) * AWS +
  415.      &  BG(I,J,K)*AGS
  416.       AN(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO4-RSO)) * AON +
  417.      & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW4-RSW)) * AWN +
  418.      &  BG(I,J,K)*AGN
  419.       AT(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO5-RSO)) * AOT +
  420.      & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW5-RSW)) * AWT +
  421.      &  BG(I,J,K)*AGT
  422.       AB(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO6-RSO)) * AOB +
  423.      & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW6-RSW)) * AWB +
  424.      &  BG(I,J,K)*AGB
  425.       OW(I,J,K)=AOW
  426.       OE(I,J,K)=AOE
  427.       OS(I,J,K)=AOS
  428.       ON(I,J,K)=AON
  429.       OT(I,J,K)=AOT
  430.       OB(I,J,K)=AOB
  431.       WW(I,J,K)=AWW
  432.       WE(I,J,K)=AWE
  433.       WS(I,J,K)=AWS
  434.       WN(I,J,K)=AWN
  435.       WT(I,J,K)=AWT
  436.       WB(I,J,K)=AWB
  437.       IF(KCOFF.EQ.0)GO TO 200
  438.       WRITE(IOCODE,33)
  439.       WRITE(IOCODE,2) I,J,K,MO1,MO2,MO3,MO4,MO5,MO6
  440.       WRITE(IOCODE,2) I,J,K,MW1,MW2,MW3,MW4,MW5,MW6
  441.       WRITE(IOCODE,2) I,J,K,MG1,MG2,MG3,MG4,MG5,MG6
  442.       WRITE(IOCODE,2) I,J,K,AOW,AOE,AOS,AON,AOT,AOB,BO(I,J,K),RSO
  443.       WRITE(IOCODE,2) I,J,K,AWW,AWE,AWS,AWN,AWT,AWB,BW(I,J,K),RSW
  444.       WRITE(IOCODE,2) I,J,K,AGW,AGE,AGS,AGN,AGT,AGB,BG(I,J,K)
  445.       WRITE(IOCODE,2) I,J,K,GOWT(I,J,K),QO(I,J,K),GWWT(I,J,K),QW(I,J,K),
  446.      &     GGWT(I,J,K),QG(I,J,K),QOWG(I,J,K)
  447. 200   CONTINUE
  448. C**** CALCULATE MAIN DIAGONAL AND RHS VECTOR
  449.       DO 300 K=1,KK
  450.       DO 300 J=1,JJ
  451.       DO 300 I=1,II
  452.       SUM(I,J,K)=AW(I,J,K)+AE(I,J,K)+AS(I,J,K)+AN(I,J,K)+
  453.      &    AT(I,J,K)+AB(I,J,K)
  454.       GAM(I,J,K)=VP(I,J,K)*CT(I,J,K)*DIV1
  455.       E(I,J,K)=-SUM(I,J,K) - GAM(I,J,K)
  456.       B(I,J,K)= QOWG(I,J,K) - GAM(I,J,K)*PN(I,J,K)
  457. 300   CONTINUE
  458. C*** CALC. COEF. FOR 0 PV BLOCKS.
  459. C*** ASSUMES NO PRESSURE CHANGE WITH TIME.
  460.       DO 350 K=1,KK
  461.       DO 350 J=1,JJ
  462.       DO 350 I=1,II
  463.       IF(VP(I,J,K).GT.0.0) GO TO 350
  464.       E(I,J,K)=-1.0
  465.       B(I,J,K)=-PN(I,J,K)
  466.       AW(I,J,K)=0.0
  467.       AE(I,J,K)=0.0
  468.       AS(I,J,K)=0.0
  469.       AN(I,J,K)=0.0
  470.       AT(I,J,K)=0.0
  471.       AB(I,J,K)=0.0
  472.   350 CONTINUE
  473.       IF(KSM1.EQ.0)RETURN
  474.       IF(N.NE.1.AND.N.NE.NN.AND.N.NE.KSM)RETURN
  475.       WRITE(IOCODE,4)
  476.       DO 404 K=1,KK
  477.       DO 404 J=1,JJ
  478.       DO 404 I=1,II
  479.       WRITE(IOCODE,2) I,J,K,AT(I,J,K),AS(I,J,K),AW(I,J,K),E(I,J,K),
  480.      & AE(I,J,K)   ,AN(I,J,K),AB(I,J,K),B(I,J,K)
  481. 404   CONTINUE
  482. 4     FORMAT(//T3,'NODE      AT(I,J,K)     AS(I,J,K)       AW(I,J,K)',
  483.      &'      E(I,J,K)     AE(I,J,K)       AN(I,J,K)      AB(I,J,K)',
  484.      &'      B(I,J,K)'/)
  485. 2     FORMAT(1X,3I3,8E15.6)
  486. 33    FORMAT(//)
  487.       RETURN
  488.       END
  489.