      MODULE GWFSPGMODULE
        INTEGER,SAVE,POINTER  ::NSEEP,MXSPG,NSPGVL,ISPGCB,IPRSPG
        INTEGER,SAVE,POINTER  ::NPSPG,ISPGPB,NNPSPG
        CHARACTER(LEN=16),SAVE, DIMENSION(:),     POINTER   ::SPGAUX
        REAL,             SAVE, DIMENSION(:,:),   POINTER   ::SEEP
        REAL,             SAVE, DIMENSION(:,:,:), POINTER   ::QSPG
      TYPE GWFSPGTYPE
        INTEGER,POINTER  ::NSEEP,MXSPG,NSPGVL,ISPGCB,IPRSPG
        INTEGER,POINTER  ::NPSPG,ISPGPB,NNPSPG
        CHARACTER(LEN=16), DIMENSION(:),     POINTER     ::SPGAUX
        REAL,              DIMENSION(:,:),   POINTER     ::SEEP
        REAL,              DIMENSION(:,:,:), POINTER     ::QSPG
      END TYPE
      TYPE(GWFSPGTYPE), SAVE:: GWFSPGDAT(10)
      END MODULE GWFSPGMODULE



      SUBROUTINE GWF2SPG7AR(IN,IGRID)
C     ******************************************************************
C     ALLOCATE ARRAY STORAGE FOR SEEPS AND READ PARAMETER DEFINITIONS
C     ******************************************************************
C
C     SPECIFICATIONS:
C     ------------------------------------------------------------------
      USE GLOBAL,       ONLY:IOUT,NCOL,NROW,NLAY,IFREFM
      USE GWFSPGMODULE, ONLY:NSEEP,MXSPG,NSPGVL,ISPGCB,IPRSPG,NPSPG,
     1                       ISPGPB,NNPSPG,SPGAUX,SEEP,QSPG
      CHARACTER*200 LINE
C     ------------------------------------------------------------------
      ALLOCATE(NSEEP,MXSPG,NSPGVL,ISPGCB,IPRSPG)
      ALLOCATE(NPSPG,ISPGPB,NNPSPG)
C
C1------IDENTIFY PACKAGE AND INITIALIZE NSEEP.
      WRITE(IOUT,1)IN
    1 FORMAT(1X,/1X,'SPG -- SEEPAGE PACKAGE, VERSION 7, 7/4/2021',
     1' INPUT READ FROM UNIT ',I4)
      NSEEP=0
      NNPSPG=0
C
C2------READ MAXIMUM NUMBER OF SEEPS AND UNIT OR FLAG FOR
C2------CELL-BY-CELL FLOW TERMS.
      CALL URDCOM(IN,IOUT,LINE)
      CALL UPARLSTAL(IN,IOUT,LINE,NPSPG,MXPD)
      IF(IFREFM.EQ.0) THEN
         READ(LINE,'(2I10)') MXACTS,ISPGCB
         LLOC=21
      ELSE
         LLOC=1
         CALL URWORD(LINE,LLOC,ISTART,ISTOP,2,MXACTS,R,IOUT,IN)
         CALL URWORD(LINE,LLOC,ISTART,ISTOP,2,ISPGCB,R,IOUT,IN)
      END IF
      WRITE(IOUT,3) MXACTS
    3 FORMAT(1X,'MAXIMUM OF ',I6,' ACTIVE SEEPAGES AT ONE TIME')
      IF(ISPGCB.LT.0) WRITE(IOUT,7)
    7 FORMAT(1X,'CELL-BY-CELL FLOWS WILL BE PRINTED WHEN ICBCFL NOT 0')
         IF(ISPGCB.GT.0) WRITE(IOUT,8) ISPGCB
    8 FORMAT(1X,'CELL-BY-CELL FLOWS WILL BE SAVED ON UNIT ',I4)
C
C3------READ AUXILIARY VARIABLES AND CBC ALLOCATION OPTION.
      ALLOCATE (SPGAUX(20))
      NAUX=0
      IPRSPG=1
   10 CALL URWORD(LINE,LLOC,ISTART,ISTOP,1,N,R,IOUT,IN)
      IF(LINE(ISTART:ISTOP).EQ.'AUXILIARY' .OR.
     1        LINE(ISTART:ISTOP).EQ.'AUX') THEN
         CALL URWORD(LINE,LLOC,ISTART,ISTOP,1,N,R,IOUT,IN)
         IF(NAUX.LT.20) THEN
            NAUX=NAUX+1
            SPGAUX(NAUX)=LINE(ISTART:ISTOP)
            WRITE(IOUT,12) SPGAUX(NAUX)
   12       FORMAT(1X,'AUXILIARY SEEPAGE VARIABLE: ',A)
         END IF
         GO TO 10
      ELSE IF(LINE(ISTART:ISTOP).EQ.'NOPRINT') THEN
         WRITE(IOUT,13)
   13    FORMAT(1X,'LISTS OF SEEPAGE CELLS WILL NOT BE PRINTED')
         IPRSPG = 0
         GO TO 10
      END IF
C3A-----THERE ARE FIVE INPUT DATA VALUES PLUS ONE LOCATION FOR
C3A-----CELL-BY-CELL FLOW.
      NSPGVL=6+NAUX
C
C4------ALLOCATE SPACE FOR SEEPAGE ARRAYs.
      ISPGPB=MXACTS+1
      MXSPG=MXACTS+MXPD
      ALLOCATE (SEEP(NSPGVL,MXSPG))
      ALLOCATE (QSPG(NCOL,NROW,NLAY))
C
C5------READ NAMED PARAMETERS.
      WRITE(IOUT,1000) NPSPG
 1000 FORMAT(1X,//1X,I5,' Seepage parameters')
      IF(NPSPG.GT.0) THEN
        LSTSUM=ISPGPB
        DO 120 K=1,NPSPG
          LSTBEG=LSTSUM
          CALL UPARLSTRP(LSTSUM,MXSPG,IN,IOUT,IP,'SPG','SPG',1,
     &                   NUMINST)
          NLST=LSTSUM-LSTBEG
          IF(NUMINST.EQ.0) THEN
C5A-----READ PARAMETER WITHOUT INSTANCES
            CALL ULSTRD(NLST,SEEP,LSTBEG,NSPGVL,MXSPG,1,IN,IOUT,
     &      'SEEP. NO.  LAYER   ROW   COL     SEEP. EL.  STRESS FACTOR',
     &        SPGAUX,20,NAUX,IFREFM,NCOL,NROW,NLAY,5,5,IPRSPG)
          ELSE
C5B-----READ INSTANCES
            NINLST=NLST/NUMINST
            DO 110 I=1,NUMINST
            CALL UINSRP(I,IN,IOUT,IP,IPRSPG)
            CALL ULSTRD(NINLST,SEEP,LSTBEG,NSPGVL,MXSPG,1,IN,IOUT,
     &      'SEEP. NO.  LAYER   ROW   COL     SEEP. EL.  STRESS FACTOR',
     &        SPGAUX,20,NAUX,IFREFM,NCOL,NROW,NLAY,5,5,IPRSPG)
            LSTBEG=LSTBEG+NINLST
  110       CONTINUE
          END IF
  120   CONTINUE
      END IF
C
C6------RETURN
      CALL SGWF2SPG7PSV(IGRID)
      RETURN
      END
      SUBROUTINE GWF2SPG7RP(IN,IGRID)
C     ******************************************************************
C     READ SEEPAGE HEAD, CONDUCTANCE AND BOTTOM ELEVATION
C     ******************************************************************
C
C     SPECIFICATIONS:
C     ------------------------------------------------------------------
      USE GLOBAL,       ONLY:IOUT,NCOL,NROW,NLAY,IFREFM
      USE GWFSPGMODULE, ONLY:NSEEP,MXSPG,NSPGVL,IPRSPG,NPSPG,
     1                       ISPGPB,NNPSPG,SPGAUX,SEEP
C     ------------------------------------------------------------------
      CALL SGWF2SPG7PNT(IGRID)
C
C1------READ ITMP (NUMBER OF SEEPAGES OR FLAG TO REUSE DATA) AND
C1------NUMBER OF PARAMETERS.
      IF(NPSPG.GT.0) THEN
         IF(IFREFM.EQ.0) THEN
            READ(IN,'(2I10)') ITMP,NP
         ELSE
            READ(IN,*) ITMP,NP
         END IF
      ELSE
         NP=0
         IF(IFREFM.EQ.0) THEN
            READ(IN,'(I10)') ITMP
         ELSE
            READ(IN,*) ITMP
         END IF
      END IF
C
C------CALCULATE SOME CONSTANTS
      NAUX=NSPGVL-6
      IOUTU = IOUT
      IF(IPRSPG.EQ.0) IOUTU=-IOUT
C
C2------DETERMINE THE NUMBER OF NON-PARAMETER SEEPAGES.
      IF(ITMP.LT.0) THEN
         WRITE(IOUT,7)
    7    FORMAT(1X,/1X,
     1        'REUSING NON-PARAMETER SEEPAGES FROM LAST STRESS PERIOD')
      ELSE
         NNPSPG=ITMP
      END IF
C
C3------IF THERE ARE NEW NON-PARAMETER SEEPAGES, READ THEM.
      MXACTS=ISPGPB-1
      IF(ITMP.GT.0) THEN
         IF(NNPSPG.GT.MXACTS) THEN
            WRITE(IOUT,99) NNPSPG,MXACTS
   99       FORMAT(1X,/1X,'THE NUMBER OF ACTIVE SEEPAGES (',I6,
     1                     ') IS GREATER THAN MXACTS(',I6,')')
            CALL USTOP(' ')
         END IF
         CALL ULSTRD(NNPSPG,SEEP,1,NSPGVL,MXSPG,1,IN,IOUT,
     1     'SEEP. NO.  LAYER   ROW   COL     SEEP. EL.  CONDUCTANCE',
     2     SPGAUX,20,NAUX,IFREFM,NCOL,NROW,NLAY,5,5,IPRSPG)
      END IF
      NSEEP=NNPSPG
C
C1C-----IF THERE ARE ACTIVE SPG PARAMETERS, READ THEM AND SUBSTITUTE
      CALL PRESET('SPG')
      IF(NP.GT.0) THEN
         NREAD=NSPGVL-1
         DO 30 N=1,NP
         CALL UPARLSTSUB(IN,'SPG',IOUTU,'SPG',SEEP,NSPGVL,MXSPG,NREAD,
     1                MXACTS,NSEEP,5,5,
     2     'SEEP. NO.  LAYER   ROW   COL     SEEP. EL.  CONDUCTANCE',
     3            SPGAUX,20,NAUX)
   30    CONTINUE
      END IF
C
C3------PRINT NUMBER OF SEEPAGES IN CURRENT STRESS PERIOD.
      WRITE (IOUT,101) NSEEP
  101 FORMAT(1X,/1X,I6,' SEEPAGES')
C
C8------RETURN.
  260 RETURN
      END
      SUBROUTINE GWF2SPG7FM(KPER,IGRID)
C     ******************************************************************
C     SETS GROUNDWATER LEVEL EQUAL TO SEEPAGE LEVEL IF CELL IS SEEPING
C     ******************************************************************
C
C        SPECIFICATIONS:
C     ------------------------------------------------------------------
      USE GLOBAL,       ONLY:NCOL,NROW,NLAY,IBOUND,HNEW,CV,CC,CR,IUNIT
      USE GWFSPGMODULE, ONLY:NSEEP,SEEP,QSPG
C
      DOUBLE PRECISION HHNEW,EEL
C     ------------------------------------------------------------------
      CALL SGWF2SPG7PNT(IGRID)
C
      ZERO=0.0
C1------INITIALIZE QSPG.
      DO 50 IL=1,NLAY
      DO 50 IR=1,NROW
      DO 50 IC=1,NCOL
      QSPG(IC,IR,IL)=ZERO
 50   CONTINUE
C
C2------IF NSEEP<=0 THERE ARE NO SEEPAGES. RETURN.
      IF(NSEEP.LE.0) RETURN
C
C3------COLLECT FLUXES IN EACH CELL, FROM ALL ACTIVE STRESSES.
C
C3B-----ADD RECHARGE FLUX TO QSPG.
      IF(IUNIT(8).GT.0) CALL SGWF2SPG7H()
C
C3L-----ADD BCF FLUX TO QSPG.
      IF(IUNIT(1).GT.0) CALL SGWF2SPG7R(KPER)
C
C4------PROCESS EACH CELL IN THE SEEPAGE LIST TO CHECK FLOW.
C4------CELL DISCHARGES GROUNDWATER IF Q<0,
C4------IS INTERMEDIATE IF STRESS FLUXES>Q>0 AND
C4------IS INFILTRATING IF Q>STRESS FLUXES.
C4------THEREFOR SUBSTRACT EXTERNAL (NON SEEPAGE)
C4------FLUXES, QSPG, FROM Q.
C4------NOTE WE CHECK FLOW FOR ALL SEEPAGE CELLS BEFORE CHECKING
C4------HEAD (C5) TO AVOID USING NEW HEAD WHEN CHECKING FLOW.
      DO 101 L=1,NSEEP
C
C4A-----GET COLUMN, ROW AND LAYER OF CELL CONTAINING SEEPAGE.
      IL=SEEP(1,L)
      IR=SEEP(2,L)
      IC=SEEP(3,L)
C
C4B-----IF SEEPAGE IS NOT ACTIVE SKIP CELL.
      IF(IBOUND(IC,IR,IL).NE.-777) GOTO 101
C
C4C-----GET HEAD AT THE CELL.
      HHNEW=HNEW(IC,IR,IL)
C
C4D-----CALCULATE FIRST VERTICAL FLUX TO SEEPAGE CELL
C4D-----SEEPAGE CELL IS IN A ONE-LAYER SYSTEM.
      FLUXV=ZERO
      IF(NLAY.GT.1) THEN
C
C4D1----SEEPAGE CELL IS IN FIRST LAYER OF A MULTI-LAYER SYSTEM.
         IF(IL.EQ.1) THEN
            FLUXV=CV(IC,IR,IL)*HHNEW-CV(IC,IR,IL)*HNEW(IC,IR,IL+1)
C
C4D2----SEEPAGE CELL IS BETWEEN FIRST LAYER AND
C4D3----BOTTOM LAYER OF A MULTI-LAYER SYSTEM.
         ELSE IF(IL.LT.NLAY.AND.IL.GT.1) THEN
            FLUXV=CV(IC,IR,IL-1)*HHNEW-
     1           CV(IC,IR,IL-1)*HNEW(IC,IR,IL-1) +
     1           CV(IC,IR,IL)*HHNEW-CV(IC,IR,IL)*HNEW(IC,IR,IL+1)
C
C4D4----SEEPAGE CELL IS IN LOWEST LAYER OF A MULTI-LAYER SYSTEM.
         ELSE
            FLUXV=CV(IC,IR,IL-1)*HHNEW-
     1           CV(IC,IR,IL-1)*HNEW(IC,IR,IL-1)
         END IF
      ENDIF
C
C4E-----FORMULATE FLUX TO SEEPAGE CELL.
      Q=FLUXV
      IF(IR.LT.NROW) Q=Q+CC(IC,IR,IL)*HHNEW-
     1     CC(IC,IR,IL)*HNEW(IC,IR+1,IL)
      IF(IR.GT.1) Q=Q+CC(IC,IR-1,IL)*HHNEW-
     1     CC(IC,IR-1,IL)*HNEW(IC,IR-1,IL)
      IF(IC.LT.NCOL) Q=Q+CR(IC,IR,IL)*HHNEW-
     1     CR(IC,IR,IL)*HNEW(IC+1,IR,IL)
      IF(IC.GT.1) Q=Q+CR(IC-1,IR,IL)*HHNEW-
     1     CR(IC-1,IR,IL)*HNEW(IC-1,IR,IL)
      Q=Q+QSPG(IC,IR,IL)
C
C4F-----IF FLUX IS LARGER THAN ZERO CELL MUST BE TURNED BACK TO INTERNAL.
C4F-----HERE WE ONLY TAG IT TO AVOID DEALING WITH IT IN C5.
      IF(Q.GT.ZERO) IBOUND(IC,IR,IL)=-778
 101  CONTINUE

C5------PROCESS EACH CELL IN THE SEEPAGE LIST TO CHECK HEAD.
      DO 102 L=1,NSEEP
C
C5A-----GET COLUMN, ROW AND LAYER OF CELL CONTAINING SEEPAGE.
      IL=SEEP(1,L)
      IR=SEEP(2,L)
      IC=SEEP(3,L)
C
C5B-----GET SEEPAGE PARAMETERS FROM SEEPAGE LIST.
      EL=SEEP(4,L)
      EEL=EL
C
C5C-----GET HEAD AT THE CELL.
      HHNEW=HNEW(IC,IR,IL)
C
C5D-----IF CELL IS INTERNAL AND HEAD HIGHER THAN SEEPAGE LEVEL
C5D-----TURN IT TO CONSTANT HEAD.
      IF(IBOUND(IC,IR,IL).GT.0 .AND. HHNEW.GT.EEL) THEN
         IBOUND(IC,IR,IL)=-777
         HNEW(IC,IR,IL)=EEL
C5E-----TURN PREVIOUSLY TAGGED CELL (C4F) INTO INTERNAL CELL.
      ELSEIF(IBOUND(IC,IR,IL).EQ.-778) THEN
         IBOUND(IC,IR,IL)=1
      ENDIF
 102  CONTINUE
C
C16-----RETURN.
      RETURN
      END
      SUBROUTINE GWF2SPG7BD(KSTP,KPER,IGRID)
C     ******************************************************************
C     CALCULATE VOLUMETRIC BUDGET FOR SEEPAGES
C     ******************************************************************
C
C        SPECIFICATIONS:
C     ------------------------------------------------------------------
      USE GLOBAL,      ONLY:IOUT,NCOL,NROW,NLAY,IBOUND,HNEW,BUFF,CV,CC,
     1                      CR,IUNIT
      USE GWFBASMODULE,ONLY:MSUM,ICBCFL,IAUXSV,DELT,PERTIM,TOTIM,
     1                      VBVL,VBNM
      USE GWFSPGMODULE,ONLY:NSEEP,ISPGCB,SEEP,NSPGVL,SPGAUX,QSPG
C
      CHARACTER*16 TEXT
      DOUBLE PRECISION HHNEW,EEL,CCSPG,CEL,RATOUT,QQ
C
      DATA TEXT /'        SEEPAGES'/
C     ------------------------------------------------------------------
      CALL SGWF2SPG7PNT(IGRID)
C
C1------INITIALIZE CELL-BY-CELL FLOW TERM FLAG (IBD) AND
C1------ACCUMULATOR (RATOUT).
      ZERO=0.
      RATOUT=ZERO
      IBD=0
      IF(ISPGCB.LT.0 .AND. ICBCFL.NE.0) IBD=-1
      IF(ISPGCB.GT.0) IBD=ICBCFL
      IBDLBL=0
C
C2------IF CELL-BY-CELL FLOWS WILL BE SAVED AS A LIST, WRITE HEADER.
      IF(IBD.EQ.2) THEN
         NAUX=NSPGVL-6
         IF(IAUXSV.EQ.0) NAUX=0
         CALL UBDSV4(KSTP,KPER,TEXT,NAUX,SPGAUX,ISPGCB,NCOL,NROW,NLAY,
     1          NSEEP,IOUT,DELT,PERTIM,TOTIM,IBOUND)
      END IF
C
C3------CLEAR THE BUFFER AND QSPG.
      DO 50 IL=1,NLAY
      DO 50 IR=1,NROW
      DO 50 IC=1,NCOL
      BUFF(IC,IR,IL)=ZERO
      QSPG(IC,IR,IL)=ZERO
50    CONTINUE
C
C4------IF THERE ARE NO SEEPAGES THEN DO NOT ACCUMULATE SEEPAGE FLOW.
      IF(NSEEP.LE.0) GO TO 200
C
C5------COLLECT FLUXES IN EACH CELL, FROM ALL ACTIVE STRESSES.
C
C5B-----ADD RECHARGE FLUX TO QSPG.
      IF(IUNIT(8).GT.0) CALL SGWF2SPG7H()
C
C5L-----ADD BCF FLUX TO QSPG.
      IF(IUNIT(1).GT.0) CALL SGWF2SPG7R(KPER)
C
C6------LOOP THROUGH EACH SEEPAGE CALCULATING FLOW.
      DO 100 L=1,NSEEP
      Q=ZERO
C
C6A-----GET LAYER, ROW & COLUMN OF CELL CONTAINING SEEPAGE.
      IL=SEEP(1,L)
      IR=SEEP(2,L)
      IC=SEEP(3,L)
C
C6B-----IF SEEPAGE IS NOT ACTIVE SKIP CALCULATION.
      IF(IBOUND(IC,IR,IL).NE.-777) GO TO 60
C
C6C-----SET TEMPORARY CONSTANT HEAD CELL BACK TO INTERNAL CELL.
      IBOUND(IC,IR,IL)=1
C
C6C-----GET HEAD AT THE CELL.
      HHNEW=HNEW(IC,IR,IL)
C
C7A----SEEPAGE CELL IS IN A ONE-LAYER SYSTEM.
      FLUXV=ZERO
      IF(NLAY.GT.1) THEN
C
C7B----SEEPAGE CELL IS IN FIRST LAYER OF A MULTI-LAYER SYSTEM.
         IF(IL.EQ.1) THEN
            FLUXV=CV(IC,IR,IL)*HHNEW-CV(IC,IR,IL)*HNEW(IC,IR,IL+1)
C
C7C----SEEPAGE CELL IS BETWEEN FIRST LAYER AND
C7C----BOTTOM LAYER OF A MULTI-LAYER SYSTEM.
         ELSE IF(IL.LT.NLAY.AND.IL.GT.1) THEN
            FLUXV=CV(IC,IR,IL-1)*HHNEW-
     1           CV(IC,IR,IL-1)*HNEW(IC,IR,IL-1) +
     1           CV(IC,IR,IL)*HHNEW-CV(IC,IR,IL)*HNEW(IC,IR,IL+1)
C
C7D----SEEPAGE CELL IS IN LOWEST LAYER OF A MULTI-LAYER SYSTEM.
         ELSE
            FLUXV=CV(IC,IR,IL-1)*HHNEW-
     1           CV(IC,IR,IL-1)*HNEW(IC,IR,IL-1)
         END IF
      ENDIF
C
C7E----FORMULATE FLUX TO SEEPAGE CELL.
      Q=FLUXV
      IF(IR.LT.NROW) Q=Q+CC(IC,IR,IL)*HHNEW-
     1     CC(IC,IR,IL)*HNEW(IC,IR+1,IL)
      IF(IR.GT.1) Q=Q+CC(IC,IR-1,IL)*HHNEW-
     1     CC(IC,IR-1,IL)*HNEW(IC,IR-1,IL)
      IF(IC.LT.NCOL) Q=Q+CR(IC,IR,IL)*HHNEW-
     1     CR(IC,IR,IL)*HNEW(IC+1,IR,IL)
      IF(IC.GT.1) Q=Q+CR(IC-1,IR,IL)*HHNEW-
     1     CR(IC-1,IR,IL)*HNEW(IC-1,IR,IL)
      Q=Q+QSPG(IC,IR,IL)
C
C8----ACCUMULATE FLOW OUT OF CELL IF THERE IS GROUNDWATER DISCHARGE.
      IF(Q.LT.ZERO) THEN
         RATOUT=RATOUT-Q
      ELSE
         WRITE(IOUT,59) Q,IL,IR,IC
 59      FORMAT(1X,'WARNING THERE IS AN UNACCOUNTED FLUX OF (',
     1        1PG15.6,') IN SEEPAGE CELL (',I3,',',I3,',',I3,').',
     2        'MAYBE YOUR CONVERGENCE CRITERIA IS TOO BIG (CHECK '
     3        'TOTAL WATER BALANCE) OR YOU USE AN UNSUPPORTED ',
     4        'PACKAGE.')
      ENDIF
C
C9-----PRINT THE INDIVIDUAL RATES IF REQUESTED(ISPGCB<0).
 60   IF(IBD.LT.0) THEN
         IF(IBDLBL.EQ.0) WRITE(IOUT,61) TEXT,KPER,KSTP
   61    FORMAT(1X,/1X,A,'   PERIOD ',I4,'   STEP ',I3)
         WRITE(IOUT,62) L,IL,IR,IC,Q
   62    FORMAT(1X,'SEEPAGE ',I6,'   LAYER ',I3,'   ROW ',I5,
     1       '   COL ',I5,'   RATE ',1PG15.6)
         IBDLBL=1
      END IF
C
C10-----ADD Q TO BUFFER.
      BUFF(IC,IR,IL)=BUFF(IC,IR,IL)+Q
C
C11-----IF SAVING CELL-BY-CELL FLOWS IN A LIST, WRITE FLOW.  ALSO
C11-----COPY FLOW TO SEEP.
   99 IF(IBD.EQ.2) CALL UBDSVB(ISPGCB,NCOL,NROW,IC,IR,IL,Q,
     1                  SEEP(:,L),NSPGVL,NAUX,6,IBOUND,NLAY)
      SEEP(NSPGVL,L)=Q
  100 CONTINUE
C
C12------IF CELL-BY-CELL FLOW WILL BE SAVED AS A 3-D ARRAY,
C12------CALL UBUDSV TO SAVE THEM.
      IF(IBD.EQ.1) CALL UBUDSV(KSTP,KPER,TEXT,ISPGCB,BUFF,NCOL,NROW,
     1                          NLAY,IOUT)
C
C13------MOVE RATES,VOLUMES & LABELS INTO ARRAYS FOR PRINTING.
  200 ROUT=RATOUT
      VBVL(3,MSUM)=ZERO
      VBVL(4,MSUM)=ROUT
      VBVL(2,MSUM)=VBVL(2,MSUM)+ROUT*DELT
      VBNM(MSUM)=TEXT
C
C14------INCREMENT BUDGET TERM COUNTER.
      MSUM=MSUM+1
C
C15------RETURN.
      RETURN
      END
      SUBROUTINE SGWF2SPG7H()
C     ******************************************************************
C     ADD RECHARGE TO QSPG
C     ******************************************************************
C
C        SPECIFICATIONS:
C     ------------------------------------------------------------------
      USE GLOBAL,      ONLY:NCOL,NROW,NLAY,IBOUND
      USE GWFSPGMODULE,ONLY:QSPG
      USE GWFRCHMODULE,ONLY:NRCHOP,IRCH,RECH
C     ------------------------------------------------------------------
C
C1------IF NRCHOP IS 1 RECHARGE IS IN TOP LAYER.
      IF(NRCHOP.NE.1) GO TO 15
      DO 10 IR=1,NROW
      DO 10 IC=1,NCOL
C
C2A-----IF CELL IS EXTERNAL THERE IS NO RECHARGE INTO IT.
         IF(IBOUND(IC,IR,1).LE.0 .AND. IBOUND(IC,IR,1).NE.-777) GO TO 10
C
C2B-----ADD RECHARGE RATE TO QSPG.
         QSPG(IC,IR,1)=QSPG(IC,IR,1)-RECH(IC,IR)
 10   CONTINUE
      GO TO 100
C
C3------IF OPTION IS 2 THEN RECHARGE IS IN INDICATOR ARRAY.
   15 IF(NRCHOP.NE.2) GO TO 25
      DO 20 IR=1,NROW
      DO 20 IC=1,NCOL
C
C3A-----LAYER INDEX IS IN INDICATOR ARRAY.
         IL=IRCH(IC,IR)
C
C3B-----IF THE CELL IS EXTERNAL THERE IS NO RECHARGE INTO IT.
         IF(IBOUND(IC,IR,IL).LE.0 .AND. IBOUND(IC,IR,IL).NE.-777)
     1      GO TO 20
C
C3C-----ADD RECHARGE RATE TO QSPG.
         QSPG(IC,IR,IL)=QSPG(IC,IR,IL)-RECH(IC,IR)
 20   CONTINUE
      GO TO 100
C
C4------IF OPTION IS 3 RECHARGE IS INTO HIGHEST INTERNAL CELL.
   25 IF(NRCHOP.NE.3) GO TO 100
      DO 30 IR=1,NROW
      DO 30 IC=1,NCOL
      DO 28 IL=1,NLAY
C
C4A-----IF CELL IS CONSTANT HEAD MOVE ON TO NEXT HORIZONTAL LOCATION.
         IF(IBOUND(IC,IR,IL).LT.0 .AND. IBOUND(IC,IR,IL).NE.-777)
     1      GO TO 30
C
C4B-----IF CELL IS INACTIVE MOVE DOWN A LAYER.
         IF (IBOUND(IC,IR,IL).EQ.0) GO TO 28
C
C4C-----SET LAYER INDEX ARRAY EQUAL TO ACTIVE LAYER NUMBER.
         QSPG(IC,IR,IL)=QSPG(IC,IR,IL)-RECH(IC,IR)
         GO TO 30
 28   CONTINUE
 30   CONTINUE
 100  CONTINUE
C
C5------RETURN
      RETURN
      END
      SUBROUTINE SGWF2SPG7R(KPER)
C     ******************************************************************
C     ADD STORAGE AND LEAKAGE CORRECTION TO QSPG
C     ******************************************************************
C
C        SPECIFICATIONS:
C     ------------------------------------------------------------------
      USE GLOBAL,      ONLY:NCOL,NROW,NLAY,IBOUND,HNEW,HOLD,CV,BOTM,
     1     LBOTM,ISSFLG
      USE GWFSPGMODULE,ONLY:QSPG
      USE GWFBASMODULE,ONLY:DELT
      USE GWFBCFMODULE,ONLY:LAYCON,SC1,SC2
C     ------------------------------------------------------------------
C
      ISS=ISSFLG(KPER)
C
      KT=0
      ONE=1.
      TLED=ONE/DELT
C
C1------IF THE SIMULATION IS TRANSIENT DETERMINE FLUX FROM STORAGE
      IF(ISS.NE.0) GO TO 201
      KT=0
      DO 200 K=1,NLAY
C
C2------SEE IF THIS LAYER IS CONVERTIBLE OR NON-CONVERTIBLE.
      IF(LAYCON(K).EQ.3 .OR. LAYCON(K).EQ.2) GO TO 150
C3------NON-CONVERTIBLE LAYER, SO USE PRIMARY STORAGE.
      DO 140 I=1,NROW
      DO 140 J=1,NCOL
      IF(IBOUND(J,I,K).LE.0 .AND. IBOUND(J,I,K).NE.-777) GO TO 140
      RHO=SC1(J,I,K)*TLED
      QSPG(J,I,K)=QSPG(J,I,K)-(RHO*HOLD(J,I,K)-RHO*HNEW(J,I,K))
  140 CONTINUE
      GO TO 200
C
C4------A CONVERTIBLE LAYER, SO CHECK OLD AND NEW HEADS TO DETERMINE
C4------WHEN TO USE PRIMARY AND SECONDARY STORAGE
  150 KT=KT+1
      DO 180 I=1,NROW
      DO 180 J=1,NCOL
C
C5A-----IF THE CELL IS EXTERNAL THEN SKIP IT.
      IF(IBOUND(J,I,K).LE.0 .AND. IBOUND(J,I,K).NE.-777) GO TO 180
      TP=BOTM(J,I,LBOTM(K)-1)
      RHO2=SC2(J,I,KT)*TLED
      RHO1=SC1(J,I,K)*TLED
C
C5B-----FIND STORAGE FACTOR AT START OF TIME STEP.
      SOLD=RHO2
      IF(HOLD(J,I,K).GT.TP) SOLD=RHO1
C
C5C-----FIND STORAGE FACTOR AT END OF TIME STEP.
      HTMP=HNEW(J,I,K)
      SNEW=RHO2
      IF(HTMP.GT.TP) SNEW=RHO1
C
C5D-----ADD FLUX TO QSPG.
      QSPG(J,I,K)=QSPG(J,I,K)-(SOLD*(HOLD(J,I,K)-TP)-SNEW*(HTMP-TP))
C
  180 CONTINUE
C
  200 CONTINUE
C
C6------FOR EACH LAYER DETERMINE IF CORRECTION TERMS ARE NEEDED FOR
C6------FLOW DOWN INTO PARTIALLY SATURATED LAYERS.
  201 DO 300 K=1,NLAY
C
C7------SEE IF CORRECTION IS NEEDED FOR LEAKAGE FROM ABOVE.
      IF(LAYCON(K).NE.3 .AND. LAYCON(K).NE.2) GO TO 250
      IF(K.EQ.1) GO TO 250
C
C7A-----FOR EACH CELL MAKE THE CORRECTION IF NEEDED.
      DO 220 I=1,NROW
      DO 220 J=1,NCOL
C
C7B-----IF THE CELL IS EXTERNAL(IBOUND<=0) THEN SKIP IT.
      IF(IBOUND(J,I,K).LE.0 .AND. IBOUND(J,I,K).NE.-777) GO TO 220
      HTMP=HNEW(J,I,K)
C
C7C-----IF HEAD IS ABOVE TOP THEN CORRECTION NOT NEEDED
      IF(HTMP.GE.BOTM(J,I,LBOTM(K)-1)) GO TO 220
C
C7D-----WITH HEAD BELOW TOP ADD CORRECTION TERMS TO QSPG.
      QSPG(J,I,K)=QSPG(J,I,K)+CV(J,I,K-1)*(BOTM(J,I,LBOTM(K)-1)-HTMP)
  220 CONTINUE
C
C8------SEE IF THIS LAYER MAY NEED CORRECTION FOR LEAKAGE TO BELOW.
  250 IF(K.EQ.NLAY) GO TO 300
      IF(LAYCON(K+1).NE.3 .AND. LAYCON(K+1).NE.2) GO TO 300
C
C8A-----FOR EACH CELL MAKE THE CORRECTION IF NEEDED.
      DO 280 I=1,NROW
      DO 280 J=1,NCOL
C
C8B-----IF CELL IS EXTERNAL (IBOUND<=0) THEN SKIP IT.
      IF(IBOUND(J,I,K).LE.0 .AND. IBOUND(J,I,K).NE.-777) GO TO 280
C
C8C-----IF HEAD IN THE LOWER CELL IS LESS THAN TOP ADD CORRECTION
C8C-----TERM TO QSPG.
      HTMP=HNEW(J,I,K+1)
      IF(HTMP.LT.BOTM(J,I,LBOTM(K+1)-1)) QSPG(J,I,K)=QSPG(J,I,K)
     1     -CV(J,I,K)*(BOTM(J,I,LBOTM(K+1)-1)-HTMP)
  280 CONTINUE
  300 CONTINUE
C
C9------RETURN
      RETURN
      END
      SUBROUTINE GWF2SPG7DA(IGRID)
C  Deallocate SPG MEMORY
      USE GWFSPGMODULE
C
        CALL SGWF2SPG7PNT(IGRID)
        DEALLOCATE(NSEEP)
        DEALLOCATE(MXSPG)
        DEALLOCATE(NSPGVL)
        DEALLOCATE(ISPGCB)
        DEALLOCATE(IPRSPG)
        DEALLOCATE(NPSPG)
        DEALLOCATE(ISPGPB)
        DEALLOCATE(NNPSPG)
        DEALLOCATE(SPGAUX)
        DEALLOCATE(SEEP)
        DEALLOCATE(QSPG)
C
      RETURN
      END
      SUBROUTINE SGWF2SPG7PNT(IGRID)
C  Change SPG data to a different grid.
      USE GWFSPGMODULE
C
        NSEEP=>GWFSPGDAT(IGRID)%NSEEP
        MXSPG=>GWFSPGDAT(IGRID)%MXSPG
        NSPGVL=>GWFSPGDAT(IGRID)%NSPGVL
        ISPGCB=>GWFSPGDAT(IGRID)%ISPGCB
        IPRSPG=>GWFSPGDAT(IGRID)%IPRSPG
        NPSPG=>GWFSPGDAT(IGRID)%NPSPG
        ISPGPB=>GWFSPGDAT(IGRID)%ISPGPB
        NNPSPG=>GWFSPGDAT(IGRID)%NNPSPG
        SPGAUX=>GWFSPGDAT(IGRID)%SPGAUX
        SEEP=>GWFSPGDAT(IGRID)%SEEP
        QSPG=>GWFSPGDAT(IGRID)%QSPG
C
      RETURN
      END
      SUBROUTINE SGWF2SPG7PSV(IGRID)
C  Save SPG data for a grid.
      USE GWFSPGMODULE
C
        GWFSPGDAT(IGRID)%NSEEP=>NSEEP
        GWFSPGDAT(IGRID)%MXSPG=>MXSPG
        GWFSPGDAT(IGRID)%NSPGVL=>NSPGVL
        GWFSPGDAT(IGRID)%ISPGCB=>ISPGCB
        GWFSPGDAT(IGRID)%IPRSPG=>IPRSPG
        GWFSPGDAT(IGRID)%NPSPG=>NPSPG
        GWFSPGDAT(IGRID)%ISPGPB=>ISPGPB
        GWFSPGDAT(IGRID)%NNPSPG=>NNPSPG
        GWFSPGDAT(IGRID)%SPGAUX=>SPGAUX
        GWFSPGDAT(IGRID)%SEEP=>SEEP
        GWFSPGDAT(IGRID)%QSPG=>QSPG
C
      RETURN
      END
