$DO66
C................................................................SOLONE
      SUBROUTINE SOLONE(II,JJ,KK,DIV1,D288,
     &KSM,KSM1,N,NN,KCOFF)
C      MACHINE DEPENDENT INCLUDE STATEMENT
$INCLUDE:'PARAMS.FOR'
C      FLOW EQ COEF WITH SINGLE POINT UPSTREAM REL PERMS
      REAL KROT,KROGT,KRWT,KRGT,MUOT,MUWT,MUGT,KX,KY,KZ
      REAL MUO4,MUW4,MUG4,MUO5,MUW5,MUG5,MUO6,MUW6,MUG6
     &,KRO1,KRW1,KRG1,KRO2,KRW2,KRG2,KRO3,KRW3,KRG3
     &,MUO1,MUW1,MUG1,MUO2,MUW2,MUG2,MUO3,MUW3,MUG3
     &,KRO4,KRW4,KRG4,KRO5,KRW5,KRG5,KRO6,KRW6,KRG6
     &,MO1,MW1,MG1,MO2,MW2,MG2,MO3,MW3,MG3
     &,MO4,MW4,MG4,MO5,MW5,MG5,MO6,MW6,MG6
     &,MUO,MUW,MUG
      COMMON /BUBBLE/ PBO,VSLOPE(LP8),BSLOPE(LP8),RSLOPE(LP8),PMAXT,
     & IREPRS,MPGT(LP8),
     & RHOSCO(LP8),RHOSCW(LP8),RHOSCG(LP8),MSAT(LP7),MPOT(LP8),
     & MPWT(LP8),PBOT(LP1,LP2,LP3),PBOTN(LP1,LP2,LP3)
      COMMON /COEF/ AW(LP1,LP2,LP3),AE(LP1,LP2,LP3),AN(LP1,LP2,LP3),
     & AS(LP1,LP2,LP3),AB(LP1,LP2,LP3),AT(LP1,LP2,LP3),E(LP1,LP2,LP3),
     & B(LP1,LP2,LP3)
      COMMON /SARRAY/ PN(LP1,LP2,LP3),IOCODE,IDMAX,
     & SON(LP1,LP2,LP3),SWN(LP1,LP2,LP3),SGN(LP1,LP2,LP3),
     & A1(LP1,LP2,LP3),A2(LP1,LP2,LP3),A3(LP1,LP2,LP3),
     & SUM(LP1,LP2,LP3),GAM(LP1,LP2,LP3),QS(LP1,LP2,LP3)
      COMMON /SPARM/ KX(LP1,LP2,LP3),KY(LP1,LP2,LP3),KZ(LP1,LP2,LP3),
     & EL(LP1,LP2,LP3),TX(LP4,LP2,LP3),TY(LP1,LP5,LP3),TZ(LP1,LP2,LP6),
     & PDAT(LP1,LP2,LP3),PDATUM,GRAD
      COMMON /SPRTPS/ P(LP1,LP2,LP3),SO(LP1,LP2,LP3),SW(LP1,LP2,LP3),
     & SG(LP1,LP2,LP3)
      COMMON /SPVT/ SAT(LP7,LP9),KROT(LP7,LP9),KRWT(LP7,LP9),
     & BGT(LP7,LP9),
     & KRGT(LP7,LP9),ITHREE(LP7),RSOT(LP7,LP9),BWPT(LP7,LP9),
     & PCOWT(LP7,LP9),PCGOT(LP7,LP9),KROGT(LP7,LP9),SWR(LP7),
     & POT(LP7,LP9),MUOT(LP7,LP9),BOT(LP7,LP9),BOPT(LP7,LP9),
     & RSOPT(LP7,LP9),PWT(LP7,LP9),MUWT(LP7,LP9),BWT(LP7,LP9),
     & RSWT(LP7,LP9),RSWPT(LP7,LP9),PGT(LP7,LP9),MUGT(LP7,LP9),
     & BGPT(LP7,LP9),CRT(LP7,LP9),IPVT(LP1,LP2,LP3),IROCK(LP1,LP2,LP3),
     & NROCK,NPVT,PSIT(LP7,LP9),PRT(LP7,LP9),WOROCK(LP7),GOROCK(LP7)
      COMMON /SSOLN/ BO(LP1,LP2,LP3),BW(LP1,LP2,LP3),BG(LP1,LP2,LP3),
     & QO(LP1,LP2,LP3),QW(LP1,LP2,LP3),QG(LP1,LP2,LP3),
     & GOWT(LP1,LP2,LP3),GWWT(LP1,LP2,LP3),GGWT(LP1,LP2,LP3),
     & OW(LP4,LP2,LP3),OE(LP4,LP2,LP3),WW(LP4,LP2,LP3),WE(LP4,LP2,LP3),
     & OS(LP1,LP5,LP3),ON(LP1,LP5,LP3),WS(LP1,LP5,LP3),WN(LP1,LP5,LP3),
     & OT(LP1,LP2,LP6),OB(LP1,LP2,LP6),WT(LP1,LP2,LP6),WB(LP1,LP2,LP6),
     & QOWG(LP1,LP2,LP3),VP(LP1,LP2,LP3),CT(LP1,LP2,LP3)
      COMMON /VECTOR/ DX(LP1,LP2,LP3),DY(LP1,LP2,LP3),DZ(LP1,LP2,LP3),
     & DZNET(LP1,LP2,LP3),IQN1(LP11),IQN2(LP11),IQN3(LP11),IHEDIN(80)
      DIMENSION RPW(LP1,LP2,LP3),RPG(LP1,LP2,LP3),RPOW(LP1,LP2,LP3),
     & RPO3(LP1,LP2,LP3),CAPOW(LP1,LP2,LP3),CAPGO(LP1,LP2,LP3)
      EQUIVALENCE (RPW,DX), (RPG,DY), (RPOW,DZ), (RPO3,KX),
     & (CAPOW,KY), (CAPGO,KZ)
      DATA RSO1,RSO2,RSO3,RSO4,RSO5,RSO6/6*0.0/
      DATA RSW1,RSW2,RSW3,RSW4,RSW5,RSW6/6*0.0/
      DO 2000 K=1,KK
      DO 2000 J=1,JJ
      DO 2000 I=1,II
      RPW(I,J,K)=0.
      RPG(I,J,K)=0.
      RPOW(I,J,K)=0.
      RPO3(I,J,K)=0.
      CAPOW(I,J,K)=0.
      CAPGO(I,J,K)=0.
      IF(VP(I,J,K).LE.0.0) GO TO 2000
      SSO=SO(I,J,K)
      SSW=SW(I,J,K)
      SSG=SG(I,J,K)
      IROCKB=IROCK(I,J,K)
      CALL INTERP(IROCKB,SAT,KROT,MSAT(IROCKB),SSO,KRO1)
      RPOW(I,J,K)=KRO1
      CALL INTERP(IROCKB,SAT,KRWT,MSAT(IROCKB),SSW,KRW1)
      RPW(I,J,K)=KRW1
      CALL INTERP(IROCKB,SAT,KRGT,MSAT(IROCKB),SSG,KRG1)
      RPG(I,J,K)=KRG1
      CALL TRIKRO(IROCKB,SSO,SSW,KRO3)
      RPO3(I,J,K)=KRO3
      CALL INTERP(IROCKB,SAT,PCOWT,MSAT(IROCKB),SSW,PCOW)
      CAPOW(I,J,K)=PCOW
      CALL INTERP(IROCKB,SAT,PCGOT,MSAT(IROCKB),SSG,PCGO)
      CAPGO(I,J,K)=PCGO
 2000 CONTINUE
      DO 200 K=1,KK
      DO 200 J=1,JJ
      DO 200 I=1,II
      IF(VP(I,J,K).LE.0.0) GO TO 200
C      I,J,K BLOCK
      PP=P(I,J,K)
      BPT=PBOT(I,J,K)
      IPVTR=IPVT(I,J,K)
      CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),PP,RSO)
      CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),PP,MUO)
      CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),PP,RSW)
      CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),PP,MUW)
      CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),PP,MUG)
      SSO=SO(I,J,K)
      SSW=SW(I,J,K)
      SSG=SG(I,J,K)
      IROCKB=IROCK(I,J,K)
      PCOW=CAPOW(I,J,K)
      PCGO=CAPGO(I,J,K)
      RO=(RHOSCO(IPVTR) + RSO*RHOSCG(IPVTR))/BO(I,J,K)
      RW=(RHOSCW(IPVTR) + RSW*RHOSCG(IPVTR))/BW(I,J,K)
      RG=RHOSCG(IPVTR)/BG(I,J,K)
      IF(I.EQ.1)GO TO 115
      IF(VP(I-1,J,K).LE.0.0) GO TO 115
C      I-1,J,K BLOCK
      P1=P(I-1,J,K)
      BPT=PBOT(I-1,J,K)
      IPVTR=IPVT(I-1,J,K)
      CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P1,RSO1)
      CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P1,MUO1)
      CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P1,RSW1)
      CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P1,MUW1)
      CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P1,MUG1)
      SO1S=SO(I-1,J,K)
      SW1S=SW(I-1,J,K)
      SG1S=SG(I-1,J,K)
      IROCKR=IROCK(I-1,J,K)
      PCOW1=CAPOW(I-1,J,K)
      PCGO1=CAPGO(I-1,J,K)
      RO1=(RHOSCO(IPVTR) + RSO1*RHOSCG(IPVTR))/BO(I-1,J,K)
      RW1=(RHOSCW(IPVTR) + RSW1*RHOSCG(IPVTR))/BW(I-1,J,K)
      RG1=RHOSCG(IPVTR)/BG(I-1,J,K)
      FACT=-D288*(EL(I-1,J,K)-EL(I,J,K))
      GOW1=(RO1+RO)*FACT
      GWW1=(RW1+RW)*FACT + PCOW-PCOW1
      GGW1=(RG1+RG)*FACT + PCGO1-PCGO
      P11=P1-PP
      HO1=P11+GOW1
      HW1=P11+GWW1
      HG1=P11+GGW1
      IF(ITHREE(IROCKB).EQ.1) GO TO 10
      KRO1=RPOW(I,J,K)
      IF(HO1.GE.0.) KRO1=RPOW(I-1,J,K)
      GO TO 15
   10 KRO1=RPO3(I,J,K)
      IF(HO1.GE.0.) KRO1=RPO3(I-1,J,K)
   15 CONTINUE
      KRW1=RPW(I,J,K)
      IF(HW1.GE.0.) KRW1=RPW(I-1,J,K)
      KRG1=RPG(I,J,K)
      IF(HG1.GE.0.) KRG1=RPG(I-1,J,K)
      MO1=4.0*KRO1/((BO(I-1,J,K)+BO(I,J,K)) * (MUO1+MUO))
      MW1=4.0*KRW1/((BW(I-1,J,K)+BW(I,J,K)) * (MUW1+MUW))
      MG1=4.0*KRG1/((BG(I-1,J,K)+BG(I,J,K)) * (MUG1+MUG))
115   AOW=TX(I,J,K)*MO1
      AWW=TX(I,J,K)*MW1
      AGW=TX(I,J,K)*MG1
      IF(I.EQ.II)GO TO 125
      IF(VP(I+1,J,K).LE.0.0) GO TO 125
C      I+1,J,K BLOCK
      P2=P(I+1,J,K)
      BPT=PBOT(I+1,J,K)
      IPVTR=IPVT(I+1,J,K)
      CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P2,RSO2)
      CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P2,MUO2)
      CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P2,RSW2)
      CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P2,MUW2)
      CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P2,MUG2)
      SO2=SO(I+1,J,K)
      SW2=SW(I+1,J,K)
      SG2=SG(I+1,J,K)
      IROCKR=IROCK(I+1,J,K)
      PCOW2=CAPOW(I+1,J,K)
      PCGO2=CAPGO(I+1,J,K)
      RO2=(RHOSCO(IPVTR) + RSO2*RHOSCG(IPVTR))/BO(I+1,J,K)
      RW2=(RHOSCW(IPVTR) + RSW2*RHOSCG(IPVTR))/BW(I+1,J,K)
      RG2=RHOSCG(IPVTR)/BG(I+1,J,K)
      FACT=-D288*(EL(I+1,J,K)-EL(I,J,K))
      GOW2=(RO2+RO)*FACT
      GWW2=(RW2+RW)*FACT + PCOW-PCOW2
      GGW2=(RG2+RG)*FACT + PCGO2-PCGO
      P22=P2-PP
      HO2=P22+GOW2
      HW2=P22+GWW2
      HG2=P22+GGW2
      IF(ITHREE(IROCKB).EQ.1) GO TO 20
      KRO2=RPOW(I,J,K)
      IF(HO2.GE.0.) KRO2=RPOW(I+1,J,K)
      GO TO 25
   20 KRO2=RPO3(I,J,K)
      IF(HO2.GE.0.) KRO2=RPO3(I+1,J,K)
   25 CONTINUE
      KRW2=RPW(I,J,K)
      IF(HW2.GE.0.) KRW2=RPW(I+1,J,K)
      KRG2=RPG(I,J,K)
      IF(HG2.GE.0.) KRG2=RPG(I+1,J,K)
      MO2=4.0*KRO2/((BO(I+1,J,K)+BO(I,J,K)) * (MUO2+MUO))
      MW2=4.0*KRW2/((BW(I+1,J,K)+BW(I,J,K)) * (MUW2+MUW))
      MG2=4.0*KRG2/((BG(I+1,J,K)+BG(I,J,K)) * (MUG2+MUG))
125   AOE=TX(I+1,J,K)*MO2
      AWE=TX(I+1,J,K)*MW2
      AGE=TX(I+1,J,K)*MG2
      IF(J.EQ.1)GO TO 135
      IF(VP(I,J-1,K).LE.0.0) GO TO 135
C      I,J-1,K BLOCK
      P3=P(I,J-1,K)
      BPT=PBOT(I,J-1,K)
      IPVTR=IPVT(I,J-1,K)
      CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P3,RSO3)
      CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P3,MUO3)
      CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P3,RSW3)
      CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P3,MUW3)
      CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P3,MUG3)
      SO3=SO(I,J-1,K)
      SW3=SW(I,J-1,K)
      SG3=SG(I,J-1,K)
      IROCKR=IROCK(I,J-1,K)
      PCOW3=CAPOW(I,J-1,K)
      PCGO3=CAPGO(I,J-1,K)
      RO3=(RHOSCO(IPVTR) + RSO3*RHOSCG(IPVTR))/BO(I,J-1,K)
      RW3=(RHOSCW(IPVTR) + RSW3*RHOSCG(IPVTR))/BW(I,J-1,K)
      RG3=RHOSCG(IPVTR)/BG(I,J-1,K)
      FACT=-D288*(EL(I,J-1,K)-EL(I,J,K))
      GOW3=(RO3+RO)*FACT
      GWW3=(RW3+RW)*FACT + PCOW-PCOW3
      GGW3=(RG3+RG)*FACT + PCGO3-PCGO
      P33=P3-PP
      HO3=P33+GOW3
      HW3=P33+GWW3
      HG3=P33+GGW3
      IF(ITHREE(IROCKB).EQ.1) GO TO 30
      KRO3=RPOW(I,J,K)
      IF(HO3.GE.0.) KRO3=RPOW(I,J-1,K)
      GO TO 35
   30 KRO3=RPO3(I,J,K)
      IF(HO3.GE.0.) KRO3=RPO3(I,J-1,K)
   35 CONTINUE
      KRW3=RPW(I,J,K)
      IF(HW3.GE.0.) KRW3=RPW(I,J-1,K)
      KRG3=RPG(I,J,K)
      IF(HG3.GE.0.) KRG3=RPG(I,J-1,K)
      MO3=4.0*KRO3/((BO(I,J-1,K)+BO(I,J,K)) * (MUO3+MUO))
      MW3=4.0*KRW3/((BW(I,J-1,K)+BW(I,J,K)) * (MUW3+MUW))
      MG3=4.0*KRG3/((BG(I,J-1,K)+BG(I,J,K)) * (MUG3+MUG))
135   AOS=TY(I,J,K)*MO3
      AWS=TY(I,J,K)*MW3
      AGS=TY(I,J,K)*MG3
      IF(J.EQ.JJ)GO TO 140
      IF(VP(I,J+1,K).LE.0.0) GO TO 140
C      I,J+1,K BLOCK
      P4=P(I,J+1,K)
      BPT=PBOT(I,J+1,K)
      IPVTR=IPVT(I,J+1,K)
      CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P4,RSO4)
      CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P4,MUO4)
      CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P4,RSW4)
      CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P4,MUW4)
      CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P4,MUG4)
      SO4=SO(I,J+1,K)
      SW4=SW(I,J+1,K)
      SG4=SG(I,J+1,K)
      IROCKR=IROCK(I,J+1,K)
      PCOW4=CAPOW(I,J+1,K)
      PCGO4=CAPGO(I,J+1,K)
      RO4=(RHOSCO(IPVTR) + RSO4*RHOSCG(IPVTR))/BO(I,J+1,K)
      RW4=(RHOSCW(IPVTR) + RSW4*RHOSCG(IPVTR))/BW(I,J+1,K)
      RG4=RHOSCG(IPVTR)/BG(I,J+1,K)
      FACT=-D288*(EL(I,J+1,K)-EL(I,J,K))
      GOW4=(RO4+RO)*FACT
      GWW4=(RW4+RW)*FACT + PCOW-PCOW4
      GGW4=(RG4+RG)*FACT + PCGO4-PCGO
      P44=P4-PP
      HO4=P44+GOW4
      HW4=P44+GWW4
      HG4=P44+GGW4
      IF(ITHREE(IROCKB).EQ.1) GO TO 40
      KRO4=RPOW(I,J,K)
      IF(HO4.GE.0.) KRO4=RPOW(I,J+1,K)
      GO TO 45
   40 KRO4=RPO3(I,J,K)
      IF(HO4.GE.0.) KRO4=RPOW(I,J+1,K)
   45 CONTINUE
      KRW4=RPW(I,J,K)
      IF(HW4.GE.0.) KRW4=RPW(I,J+1,K)
      KRG4=RPG(I,J,K)
      IF(HG4.GE.0.) KRG4=RPG(I,J+1,K)
      MO4=4.0*KRO4/((BO(I,J+1,K)+BO(I,J,K)) * (MUO4+MUO))
      MW4=4.0*KRW4/((BW(I,J+1,K)+BW(I,J,K)) * (MUW4+MUW))
      MG4=4.0*KRG4/((BG(I,J+1,K)+BG(I,J,K)) * (MUG4+MUG))
140   AON=TY(I,J+1,K)*MO4
      AWN=TY(I,J+1,K)*MW4
      AGN=TY(I,J+1,K)*MG4
      IF(K.EQ.1)GO TO 145
      IF(VP(I,J,K-1).LE.0.0) GO TO 145
C      I,J,K-1 BLOCK
      P5=P(I,J,K-1)
      BPT=PBOT(I,J,K-1)
      IPVTR=IPVT(I,J,K-1)
      CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P5,RSO5)
      CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P5,MUO5)
      CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P5,RSW5)
      CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P5,MUW5)
      CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P5,MUG5)
      SO5=SO(I,J,K-1)
      SW5=SW(I,J,K-1)
      SG5=SG(I,J,K-1)
      IROCKR=IROCK(I,J,K-1)
      PCOW5=CAPOW(I,J,K-1)
      PCGO5=CAPGO(I,J,K-1)
      RO5=(RHOSCO(IPVTR) + RSO5*RHOSCG(IPVTR))/BO(I,J,K-1)
      RW5=(RHOSCW(IPVTR) + RSW5*RHOSCG(IPVTR))/BW(I,J,K-1)
      RG5=RHOSCG(IPVTR)/BG(I,J,K-1)
      FACT=-D288*(EL(I,J,K-1)-EL(I,J,K))
      GOW5=(RO5+RO)*FACT
      GWW5=(RW5+RW)*FACT + PCOW-PCOW5
      GGW5=(RG5+RG)*FACT + PCGO5-PCGO
      P55=P5-PP
      HO5=P55+GOW5
      HW5=P55+GWW5
      HG5=P55+GGW5
      IF(ITHREE(IROCKB).EQ.1) GO TO 50
      KRO5=RPOW(I,J,K)
      IF(HO5.GE.0.) KRO5=RPOW(I,J,K-1)
      GO TO 55
   50 KRO5=RPO3(I,J,K)
      IF(HO5.GE.0.) KRO5=RPO3(I,J,K-1)
   55 CONTINUE
      KRW5=RPW(I,J,K)
      IF(HW5.GE.0.) KRW5=RPW(I,J,K-1)
      KRG5=RPG(I,J,K)
      IF(HG5.GE.0.) KRG5=RPG(I,J,K-1)
      MO5=4.0*KRO5/((BO(I,J,K-1)+BO(I,J,K)) * (MUO5+MUO))
      MW5=4.0*KRW5/((BW(I,J,K-1)+BW(I,J,K)) * (MUW5+MUW))
      MG5=4.0*KRG5/((BG(I,J,K-1)+BG(I,J,K)) * (MUG5+MUG))
145   AOT=TZ(I,J,K)*MO5
      AWT=TZ(I,J,K)*MW5
      AGT=TZ(I,J,K)*MG5
      IF(K.EQ.KK)GO TO 150
      IF(VP(I,J,K+1).LE.0.0) GO TO 150
C      I,J,K+1 BLOCK
      P6=P(I,J,K+1)
      BPT=PBOT(I,J,K+1)
      IPVTR=IPVT(I,J,K+1)
      CALL INTPVT(IPVTR,BPT,RSLOPE(IPVTR),POT,RSOT,MPOT(IPVTR),P6,RSO6)
      CALL INTPVT(IPVTR,BPT,VSLOPE(IPVTR),POT,MUOT,MPOT(IPVTR),P6,MUO6)
      CALL INTERP(IPVTR,PWT,RSWT,MPWT(IPVTR),P6,RSW6)
      CALL INTERP(IPVTR,PWT,MUWT,MPWT(IPVTR),P6,MUW6)
      CALL INTERP(IPVTR,PGT,MUGT,MPGT(IPVTR),P6,MUG6)
      SO6=SO(I,J,K+1)
      SW6=SW(I,J,K+1)
      SG6=SG(I,J,K+1)
      IROCKR=IROCK(I,J,K+1)
      PCOW6=CAPOW(I,J,K+1)
      PCGO6=CAPGO(I,J,K+1)
      RO6=(RHOSCO(IPVTR) + RSO6*RHOSCG(IPVTR))/BO(I,J,K+1)
      RW6=(RHOSCW(IPVTR) + RSW6*RHOSCG(IPVTR))/BW(I,J,K+1)
      RG6=RHOSCG(IPVTR)/BG(I,J,K+1)
      FACT=-D288*(EL(I,J,K+1)-EL(I,J,K))
      GOW6=(RO6+RO)*FACT
      GWW6=(RW6+RW)*FACT + PCOW-PCOW6
      GGW6=(RG6+RG)*FACT + PCGO6-PCGO
      P66=P6-PP
      HO6=P66+GOW6
      HW6=P66+GWW6
      HG6=P66+GGW6
      IF(ITHREE(IROCKB).EQ.1) GO TO 60
      KRO6=RPOW(I,J,K)
      IF(HO6.GE.0.) KRO6=RPOW(I,J,K+1)
      GO TO 65
   60 KRO6=RPO3(I,J,K)
      IF(HO6.GE.0.) KRO6=RPO3(I,J,K+1)
   65 CONTINUE
      KRW6=RPW(I,J,K)
      IF(HW6.GE.0.) KRW6=RPW(I,J,K+1)
      KRG6=RPG(I,J,K)
      IF(HG6.GE.0.) KRG6=RPG(I,J,K+1)
      MO6=4.0*KRO6/((BO(I,J,K+1)+BO(I,J,K)) * (MUO6+MUO))
      MW6=4.0*KRW6/((BW(I,J,K+1)+BW(I,J,K)) * (MUW6+MUW))
      MG6=4.0*KRG6/((BG(I,J,K+1)+BG(I,J,K)) * (MUG6+MUG))
150   AOB=TZ(I,J,K+1)*MO6
      AWB=TZ(I,J,K+1)*MW6
      AGB=TZ(I,J,K+1)*MG6
      RSO1A=0.5*(RSO1+RSO)
      RSO2A=0.5*(RSO2+RSO)
      RSO3A=0.5*(RSO3+RSO)
      RSO4A=0.5*(RSO4+RSO)
      RSO5A=0.5*(RSO5+RSO)
      RSO6A=0.5*(RSO6+RSO)
      RSW1A=0.5*(RSW1+RSW)
      RSW2A=0.5*(RSW2+RSW)
      RSW3A=0.5*(RSW3+RSW)
      RSW4A=0.5*(RSW4+RSW)
      RSW5A=0.5*(RSW5+RSW)
      RSW6A=0.5*(RSW6+RSW)
      AO1=AOW*GOW1
      AO2=AOE*GOW2
      AO3=AOS*GOW3
      AO4=AON*GOW4
      AO5=AOT*GOW5
      AO6=AOB*GOW6
      AW1=AWW*GWW1
      AW2=AWE*GWW2
      AW3=AWS*GWW3
      AW4=AWN*GWW4
      AW5=AWT*GWW5
      AW6=AWB*GWW6
      GOWT(I,J,K)= AO1 + AO2 + AO3 + AO4 + AO5 + AO6
      GWWT(I,J,K)= AW1 + AW2 + AW3 + AW4 + AW5 + AW6
      GGWT(I,J,K)=AGW*GGW1+AGE*GGW2+AGS*GGW3+AGN*GGW4+AGT*GGW5+AGB*GGW6
     &+RSO1A*AO1+RSO2A*AO2+RSO3A*AO3+RSO4A*AO4+RSO5A*AO5+RSO6A*AO6
     &+RSW1A*AW1+RSW2A*AW2+RSW3A*AW3+RSW4A*AW4+RSW5A*AW5+RSW6A*AW6
      QOWG(I,J,K)=(BO(I,J,K)-BG(I,J,K)*RSO)*(-GOWT(I,J,K)+QO(I,J,K)) +
     &   (BW(I,J,K)-BG(I,J,K)*RSW)*(-GWWT(I,J,K)+QW(I,J,K)) +
     &    BG(I,J,K)*(-GGWT(I,J,K)+QG(I,J,K))
      AW(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO1-RSO)) * AOW +
     & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW1-RSW)) * AWW +
     &  BG(I,J,K)*AGW
      AE(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO2-RSO)) * AOE +
     & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW2-RSW)) * AWE +
     &  BG(I,J,K)*AGE
      AS(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO3-RSO)) * AOS +
     & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW3-RSW)) * AWS +
     &  BG(I,J,K)*AGS
      AN(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO4-RSO)) * AON +
     & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW4-RSW)) * AWN +
     &  BG(I,J,K)*AGN
      AT(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO5-RSO)) * AOT +
     & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW5-RSW)) * AWT +
     &  BG(I,J,K)*AGT
      AB(I,J,K)=(BO(I,J,K) + 0.5*BG(I,J,K)*(RSO6-RSO)) * AOB +
     & (BW(I,J,K) + 0.5*BG(I,J,K)*(RSW6-RSW)) * AWB +
     &  BG(I,J,K)*AGB
      OW(I,J,K)=AOW
      OE(I,J,K)=AOE
      OS(I,J,K)=AOS
      ON(I,J,K)=AON
      OT(I,J,K)=AOT
      OB(I,J,K)=AOB
      WW(I,J,K)=AWW
      WE(I,J,K)=AWE
      WS(I,J,K)=AWS
      WN(I,J,K)=AWN
      WT(I,J,K)=AWT
      WB(I,J,K)=AWB
      IF(KCOFF.EQ.0)GO TO 200
      WRITE(IOCODE,33)
      WRITE(IOCODE,2) I,J,K,MO1,MO2,MO3,MO4,MO5,MO6
      WRITE(IOCODE,2) I,J,K,MW1,MW2,MW3,MW4,MW5,MW6
      WRITE(IOCODE,2) I,J,K,MG1,MG2,MG3,MG4,MG5,MG6
      WRITE(IOCODE,2) I,J,K,AOW,AOE,AOS,AON,AOT,AOB,BO(I,J,K),RSO
      WRITE(IOCODE,2) I,J,K,AWW,AWE,AWS,AWN,AWT,AWB,BW(I,J,K),RSW
      WRITE(IOCODE,2) I,J,K,AGW,AGE,AGS,AGN,AGT,AGB,BG(I,J,K)
      WRITE(IOCODE,2) I,J,K,GOWT(I,J,K),QO(I,J,K),GWWT(I,J,K),QW(I,J,K),
     &     GGWT(I,J,K),QG(I,J,K),QOWG(I,J,K)
200   CONTINUE
C**** CALCULATE MAIN DIAGONAL AND RHS VECTOR
      DO 300 K=1,KK
      DO 300 J=1,JJ
      DO 300 I=1,II
      SUM(I,J,K)=AW(I,J,K)+AE(I,J,K)+AS(I,J,K)+AN(I,J,K)+
     &    AT(I,J,K)+AB(I,J,K)
      GAM(I,J,K)=VP(I,J,K)*CT(I,J,K)*DIV1
      E(I,J,K)=-SUM(I,J,K) - GAM(I,J,K)
      B(I,J,K)= QOWG(I,J,K) - GAM(I,J,K)*PN(I,J,K)
300   CONTINUE
C*** CALC. COEF. FOR 0 PV BLOCKS.
C*** ASSUMES NO PRESSURE CHANGE WITH TIME.
      DO 350 K=1,KK
      DO 350 J=1,JJ
      DO 350 I=1,II
      IF(VP(I,J,K).GT.0.0) GO TO 350
      E(I,J,K)=-1.0
      B(I,J,K)=-PN(I,J,K)
      AW(I,J,K)=0.0
      AE(I,J,K)=0.0
      AS(I,J,K)=0.0
      AN(I,J,K)=0.0
      AT(I,J,K)=0.0
      AB(I,J,K)=0.0
  350 CONTINUE
      IF(KSM1.EQ.0)RETURN
      IF(N.NE.1.AND.N.NE.NN.AND.N.NE.KSM)RETURN
      WRITE(IOCODE,4)
      DO 404 K=1,KK
      DO 404 J=1,JJ
      DO 404 I=1,II
      WRITE(IOCODE,2) I,J,K,AT(I,J,K),AS(I,J,K),AW(I,J,K),E(I,J,K),
     & AE(I,J,K)   ,AN(I,J,K),AB(I,J,K),B(I,J,K)
404   CONTINUE
4     FORMAT(//T3,'NODE      AT(I,J,K)     AS(I,J,K)       AW(I,J,K)',
     &'      E(I,J,K)     AE(I,J,K)       AN(I,J,K)      AB(I,J,K)',
     &'      B(I,J,K)'/)
2     FORMAT(1X,3I3,8E15.6)
33    FORMAT(//)
      RETURN
      END
