1      SUBROUTINE FHB1AL(ISUM,LENX,LCFLLC,LCBDTM,LCFLRT,LCBDFV,LCBDHV,
2     &                  LCHDLC,LCSBHD,NBDTIM,NFLW,NHED,IN,IOUT,IFHBCB,
3     &                  NFHBX1,NFHBX2,IFHBD3,IFHBD4,IFHBD5,IFHBSS,ISS)
4C
5C-----VERSION 0000 10JAN1997 FHB1AL
6C     ******************************************************************
7C     ALLOCATE ARRAY STORAGE FOR FLOW AND HEAD BOUNDARY PACKAGE
8C     ******************************************************************
9C
10C        SPECIFICATIONS:
11C     ------------------------------------------------------------------
12      COMMON /FHBCOM/ FHBXNM(10),FHBXWT(10)
13      CHARACTER*16 FHBXNM
14      CHARACTER*80 LINE
15C     ------------------------------------------------------------------
16C
17C1------IDENTIFY PACKAGE
18      WRITE(IOUT,1)IN
19    1 FORMAT(1H0,'FHB1 -- SPECIFIED FLOW PACKAGE, VERSION 1,12/3/96',
20     &' INPUT READ FROM',I3)
21C
22C2------READ NUMBER OF TIMES, NUMBER OF SPECIFIED-FLOW CELLS AND
23C2------UNIT OR FLAG FOR CELL-BY-CELL FLOW TERMS, NUMBER OF
24C2------AUXILIARY VARIABLES.
25      READ(IN,*) NBDTIM,NFLW,NHED,IFHBSS,IFHBCB,NFHBX1,NFHBX2
26C
27C3------PRINT NBDTIM, STOP IF NO TIMES ARE TO BE SPECIFIED
28      IF(NFLW.LT.1.AND.NHED.LT.1) THEN
29         WRITE(IOUT,4)
30 4       FORMAT(1X,'SPECIFIED FLOW AND HEAD BOUNDARY OPTION ',
31     &   'CANCELLED.',/,1X,'NO BOUNDARY CELLS WERE SPECIFIED.')
32         IN=0
33         RETURN
34      ENDIF
35      IF(NBDTIM.LT.1) THEN
36         WRITE(IOUT,6)
37 6       FORMAT(1X, 'SIMULATION ABORTING.  NOT ENOUGH TIMES ',
38     &   'SPECIFIED FOR FHB1 PACKAGE.')
39         STOP
40      ELSE IF(NBDTIM.EQ.1) THEN
41         WRITE(IOUT,8)
42 8       FORMAT(1X,' SPECIFIED FLOW AND HEAD VALUES WILL REMAIN ',
43     &   'CONSTANT FOR ENTIRE SIMULATION.')
44      ELSE
45         WRITE(IOUT,10) NBDTIM
46 10      FORMAT(1H ,'TOTAL OF',I5,' TIMES WILL BE USED TO DEFINE ',
47     &   'VARIATIONS IN FLOW AND HEAD.')
48      ENDIF
49C
50C4------PRINT NFLW AND NHED AND STEADY-STATE OPTION
51      WRITE(IOUT,12) NFLW
52 12   FORMAT(1H ,'FLOW WILL BE SPECIFIED AT A TOTAL OF',I5,' CELLS.')
53      WRITE(IOUT,14) NHED
54 14   FORMAT(1H ,'HEAD WILL BE SPECIFIED AT A TOTAL OF',I5,' CELLS.')
55      IF(ISS.EQ.0) THEN
56         WRITE(IOUT,15)
57 15      FORMAT(1H ,'FHB STEADY-STATE OPTION FLAG WILL BE IGNORED,'/,
58     &          1H ,'SIMULATION IS TRANSIENT.')
59      ELSE
60      IF(IFHBSS.EQ.0) THEN
61         WRITE(IOUT,16)
62 16      FORMAT(1H ,'FLOW, HEAD, AND AUX VARIABLES AT TIME=0 WILL BE ',
63     &        /,1H ,'USED IN STEADY-STATE SIMULATIONS.')
64      ELSE
65         WRITE(IOUT,18)
66 18      FORMAT(1H ,'FLOW, HEAD, AND AUX VARIABLES WILL BE ',
67     &      'INTERPOLATED',/,1H ,'IN STEADY-STATE SIMULATIONS.')
68      ENDIF
69      ENDIF
70C
71C5------IF CELL-BY-CELL FLOW TERMS ARE TO BE SAVED THEN PRINT UNIT #
72      IF(IFHBCB.GT.0) WRITE(IOUT,20) IFHBCB
73 20   FORMAT(1X,'CELL-BY-CELL FLOWS WILL BE RECORDED ON UNIT',I3)
74      IF(IFHBCB.LT.0) WRITE(IOUT,24)
75 24   FORMAT(1X,'CELL-BY-CELL FLOWS WILL BE PRINTED WHEN ICBCFL NOT 0')
76C
77C6------READ AUXILIARY VARIABLES
78      IF(NFHBX1.GT.5.OR.NFHBX2.GT.5) THEN
79         WRITE(IOUT,*) ' ABORTING. A MAXIMUM OF 5 AUXILIARY VARIABLES',
80     &   ' CAN BE DEFINED BY FHB.'
81         STOP
82      ENDIF
83      WRITE(IOUT,26) NFHBX1
84 26   FORMAT(1X,I2,' AUXILIARY VARIABLES FOR SPECIFIED-FLOW CELLS WILL',
85     & /,'  BE DEFINED BY FHB FOR USE BY OTHER PACKAGES.')
86      IF(NFHBX1.LT.1) GO TO 38
87      WRITE(IOUT,28)
88 28   FORMAT('       NAME      WEIGHTING FACTOR',/,1X,32('-'))
89      DO 30 NX=1,NFHBX1
90      READ(IN,'(A)') LINE
91      LLOC=1
92      CALL URWORD(LINE,LLOC,ISTART,ISTOP,1,N,R,IOUT,IN)
93      FHBXNM(NX)=LINE(ISTART:ISTOP)
94      CALL URWORD(LINE,LLOC,ISTART,ISTOP,3,N,FHBXWT(NX),IOUT,IN)
95      WRITE(IOUT,29) FHBXNM(NX),FHBXWT(NX)
96 29   FORMAT(1X,A16,F11.2)
97      IF(FHBXWT(NX).LT.0.0.OR.FHBXWT(NX).GT.1.0) THEN
98      WRITE(IOUT,*) ' Aborting. Weights for Auxiliary variables cannot'
99      WRITE(IOUT,*) ' be less than 0.0 or greater than 1.0.'
100      STOP
101      ENDIF
102 30   CONTINUE
103 38   WRITE(IOUT,126) NFHBX2
104 126  FORMAT(1X,I2,' AUXILIARY VARIABLES FOR SPECIFIED-HEAD CELLS WILL',
105     & /,'  BE DEFINED BY FHB FOR USE BY OTHER PACKAGES.')
106      IF(NFHBX2.LT.1) GO TO 200
107      WRITE(IOUT,28)
108      DO 130 NX=1,NFHBX2
109      READ(IN,'(A)') LINE
110      LLOC=1
111      CALL URWORD(LINE,LLOC,ISTART,ISTOP,1,N,R,IOUT,IN)
112      FHBXNM(5+NX)=LINE(ISTART:ISTOP)
113      CALL URWORD(LINE,LLOC,ISTART,ISTOP,3,N,FHBXWT(5+NX),IOUT,IN)
114      WRITE(IOUT,129) FHBXNM(5+NX),FHBXWT(5+NX)
115 129  FORMAT(1X,A16,F11.2)
116      IF(FHBXWT(5+NX).LT.0.0.OR.FHBXWT(5+NX).GT.1.0) THEN
117      WRITE(IOUT,*) ' Aborting. Weights for Auxiliary variables cannot'
118      WRITE(IOUT,*) ' be less than 0.0 or greater than 1.0.'
119      STOP
120      ENDIF
121 130  CONTINUE
122C7------ALLOCATE SPACE FOR ARRAYS BDTIM, IFLLOC, FLWRAT, BDFV,
123C7------IHDLOC, SBHED, AND BDHV
124 200  IFHBD3=NBDTIM*(1+NFHBX1)
125      IFHBD4=2+NFHBX1
126      IFHBD5=NBDTIM*(1+NFHBX2)
127      LCBDTM=ISUM
128      ISUM=ISUM+NBDTIM
129      LCFLLC=ISUM
130      ISUM=ISUM+NFLW*4
131      LCFLRT=ISUM
132      ISUM=ISUM+NFLW*IFHBD3
133      LCBDFV=ISUM
134      ISUM=ISUM+NFLW*IFHBD4
135      LCHDLC=ISUM
136      ISUM=ISUM+NHED*4
137      LCSBHD=ISUM
138      ISUM=ISUM+NHED*IFHBD5
139      LCBDHV=ISUM
140      ISUM=ISUM+NHED*NFHBX2
141      ISP=ISUM-LCBDTM
142C
143C8------PRINT NUMBER OF SPACES IN X ARRAY USED BY FLOW PACKAGE.
144      WRITE(IOUT,210) ISP
145 210  FORMAT(1X,I8,' ELEMENTS IN X ARRAY ARE USED BY FHB1')
146      ISUM1=ISUM-1
147      WRITE(IOUT,220) ISUM1,LENX
148 220  FORMAT(1X,I8,' ELEMENTS OF X ARRAY USED OUT OF',I8)
149C
150C9------IF THERE ISN'T ENOUGH SPACE IN THE X ARRAY THEN PRINT
151C9------A WARNING MESSAGE.
152      IF(ISUM1.GT.LENX) WRITE(IOUT,230)
153 230  FORMAT(1X,'   ***X ARRAY MUST BE DIMENSIONED LARGER***')
154C10-----RETURN
155      RETURN
156      END
157      SUBROUTINE FHB1RP(IBOUND,NROW,NCOL,NLAY,IFLLOC,BDTIM,NBDTIM,
158     &             FLWRAT,NFLW,NHED,IHDLOC,SBHED,IN,IOUT,
159     &             NFHBX1,NFHBX2,IFHBD3,IFHBD5)
160C
161C
162C-----VERSION 0000 10JAN1997 FHB1RP
163C     ******************************************************************
164C     READ TIMES, CELL LOCATIONS, RATES, AND HEADS FOR FLOW AND HEAD
165C     BOUNDARY PACKAGE
166C     ******************************************************************
167C
168C        SPECIFICATIONS:
169C     ------------------------------------------------------------------
170      COMMON /FHBCOM/ FHBXNM(10),FHBXWT(10)
171      CHARACTER*16 FHBXNM
172      CHARACTER*1 DSH1
173      DIMENSION IBOUND(NCOL,NROW,NLAY),BDTIM(NBDTIM),IFLLOC(4,NFLW),
174     &      FLWRAT(IFHBD3,NFLW),IHDLOC(4,NHED), SBHED(IFHBD5,NHED)
175      DATA DSH1/'-'/
176C     ------------------------------------------------------------------
177C
178C1------READ TIMES AT WHICH SPECIFIED FLOW AND HEAD VALUES WILL BE READ
179      READ(IN,*) IFHBUN,CNSTM,IFHBPT
180      WRITE(IOUT,10) IFHBUN,CNSTM
181 10   FORMAT(1X,'TIMES FOR SPECIFIED-FLOW AND HEAD VALUES WILL BE READ',
182     & ' ON UNIT',I4,' AND',/,
183     &' MULTIPLIED BY',G12.4,'.')
184      READ(IFHBUN,*) (BDTIM(L),L=1,NBDTIM)
185      DO 12 L=1,NBDTIM
186      BDTIM(L)=BDTIM(L)*CNSTM
187 12   CONTINUE
188C
189C2------IF DESIRED, PRINT TABLE OF TIMES
190      IF(IFHBPT.GT.0) THEN
191         WRITE(IOUT,20) NBDTIM
192 20      FORMAT(1X,I5,' TIMES FOR SPECIFYING FLOWS AND HEADS:')
193         WRITE(IOUT,22) (L,L=1,NBDTIM)
194 22      FORMAT(16X,I8,4I12)
195         ND=MIN0(60,NBDTIM*12)
196         WRITE(IOUT,24) (DSH1,M=1,ND)
197 24      FORMAT(17X,60A1)
198         WRITE(IOUT,26) (BDTIM(L),L=1,NBDTIM)
199 26      FORMAT(17X,5G12.4)
200      ENDIF
201C
202C3------MAKE SURE THAT FIRST TIME IS ZERO AND THAT TIMES INCREASE
203      ICHK1=0
204      ICHK2=0
205      IF(BDTIM(1).NE.0.0) THEN
206         WRITE(IOUT,30)
207 30      FORMAT(1X,'STARTING TIME FOR SPECIFIED FLOWS AND HEADS MUST',
208     &   ' BE ZERO. ABORTING.')
209         ICHK1=1
210      ENDIF
211      DO 40 L=2,NBDTIM
212      IF(BDTIM(L).LT.BDTIM(L-1)) THEN
213         WRITE(IOUT,32)
214 32      FORMAT(1X,'TIMES FOR SPECIFIED FLOWS MUST INCREASE.',
215     &   '  ABORTING.')
216         ICHK2=1
217         GO TO 42
218      ENDIF
219 40   CONTINUE
220 42   IF(ICHK1.EQ.1.OR.ICHK2.EQ.1) STOP
221C
222C4A-----READ CELL INDICIES AND SPECIFIED-FLOW RATES
223      IF(NFLW.LT.1) GO TO 70
224      READ(IN,*) IFHBUN,CNSTM,IFHBPT
225      WRITE(IOUT,50) IFHBUN,CNSTM
226 50   FORMAT(/,1X,'CELL INDICIES AND SPECIFIED-FLOW RATES ',
227     & 'WILL BE READ ON UNIT',I4,'. RATES WILL',/,
228     & 1X,'BE MULTIPLIED BY',G12.4,'.')
229      IF(IFHBPT.GT.0) THEN
230         WRITE(IOUT,52)
231 52      FORMAT(1X,'LAYER  ROW  COL IAUX              FLOW RATES')
232         ND=MIN0(79,19+NBDTIM*12)
233         WRITE(IOUT,54) (DSH1,M=1,ND)
234 54      FORMAT(1X,78A1)
235      ENDIF
236      DO 59 N=1,NFLW
237      READ(IFHBUN,*) (IFLLOC(I,N),I=1,4),(FLWRAT(L,N),L=1,NBDTIM)
238      DO 56 L=1,NBDTIM
239      FLWRAT(L,N)=FLWRAT(L,N)*CNSTM
240 56   CONTINUE
241C
242C4B-----IF DESIRED, PRINT TABLE OF SPECIFIED-FLOW CELL LOCATIONS
243C4B-----AND RATES
244      IF(IFHBPT.GT.0) THEN
245         WRITE(IOUT,58) (IFLLOC(I,N),I=1,4),(FLWRAT(L,N),L=1,NBDTIM)
246 58      FORMAT(1X,I4,3I5,5G12.4,/,(20X,5G12.4))
247      ENDIF
248 59   CONTINUE
249C
250C5A------READ VALUES OF AUXILIARY VARIABLES FOR SPECIFIED-FLOW CELLS
251      IF(NFHBX1.LT.1) GO TO 70
252      DO 69 NX=1,NFHBX1
253      NS=NBDTIM*NX
254      READ(IN,*) IFHBUN,CNSTM,IFHBPT
255      WRITE(IOUT,61) FHBXNM(NX),IFHBUN,CNSTM
256 61   FORMAT(/,1X,A16,
257     & 'FOR SPECIFIED-FLOW CELLS WILL BE READ ON UNIT',I4,'.',/,
258     &  ' VALUES WILL BE MULTIPLIED BY',G12.4,'.')
259      IF(IFHBPT.GT.0) THEN
260         WRITE(IOUT,62) FHBXNM(NX)
261 62      FORMAT(1X,'LAYER  ROW  COL IAUX  ',A16)
262         WRITE(IOUT,54) (DSH1,M=1,ND)
263      ENDIF
264      DO 68 N=1,NFLW
265      READ(IFHBUN,*) (FLWRAT(NS+L,N),L=1,NBDTIM)
266      DO 66 L=1,NBDTIM
267      FLWRAT(NS+L,N)=FLWRAT(NS+L,N)*CNSTM
268 66   CONTINUE
269C
270C5B------IF DESIRED, PRINT TABLE OF AUXILIARY VARIABLE VALUES AT
271C5B------SPECIFIED-FLOW CELL LOCATIONS
272      IF(IFHBPT.GT.0) THEN
273         WRITE(IOUT,58) (IFLLOC(I,N),I=1,4),
274     &                        (FLWRAT(NS+L,N),L=1,NBDTIM)
275 67      FORMAT(1X,I4,2I6,5G12.4,/,(17X,5G12.4))
276      ENDIF
277 68   CONTINUE
278 69   CONTINUE
279C
280C6------READ CELL INDICIES AND SPECIFIED-HEAD VALUES
281 70   IF(NHED.LT.1) GO TO 300
282      READ(IN,*) IFHBUN,CNSTM,IFHBPT
283      WRITE(IOUT,71) IFHBUN,CNSTM
284 71   FORMAT(/,1X,'CELL INDICIES AND SPECIFIED-HEAD VALUES ',
285     & 'WILL BE READ ON UNIT',I4,'. HEAD VALUES',/,
286     & 1X,'WILL BE MULTIPLIED BY',G12.4,'.')
287      IF(IFHBPT.GT.0) THEN
288         WRITE(IOUT,72)
289 72      FORMAT(1X,'LAYER  ROW  COL IAUX             HEAD VALUES')
290         ND=MIN0(79,19+NBDTIM*12)
291         WRITE(IOUT,74) (DSH1,M=1,ND)
292 74      FORMAT(1X,79A1)
293      ENDIF
294      DO 80 N=1,NHED
295      READ(IFHBUN,*) (IHDLOC(I,N),I=1,4),(SBHED(L,N),L=1,NBDTIM)
296      DO 75 L=1,NBDTIM
297      SBHED(L,N)=SBHED(L,N)*CNSTM
298 75   CONTINUE
299C
300C7------AT SPECIFIED-HEAD LOCATIONS, SET IBOUND TO NEGATIVE NUMBER.
301C7------IGNORE SPECIFIED-HEAD CONDITIONS AT CELLS WHERE IBOUND IS ZERO
302      K=IHDLOC(1,N)
303      I=IHDLOC(2,N)
304      J=IHDLOC(3,N)
305      IF(IBOUND(J,I,K).NE.0) THEN
306         IBOUND(J,I,K)=-IABS(IBOUND(J,I,K))
307      ELSE
308         WRITE(IOUT,76) (IHDLOC(I,N),I=1,3)
309 76      FORMAT(1X,'SPECIFIED-HEAD VALUE IGNORED AT ROW',I5,', COLUMN',
310     &   I5,', AND LAYER',I5,'.')
311      ENDIF
312C
313C8------IF DESIRED, PRINT TABLE OF SPECIFIED-FLOW CELL LOCATIONS
314C8------AND RATES
315      IF(IFHBPT.GT.0) THEN
316         IF(IBOUND(J,I,K).NE.0)
317     &   WRITE(IOUT,58) (IHDLOC(I,N),I=1,4),(SBHED(L,N),L=1,NBDTIM)
318      ENDIF
319 80   CONTINUE
320C
321C9A------READ VALUES OF AUXILIARY VARIABLES FOR SPECIFIED-HEAD CELLS
322      IF(NFHBX2.LT.1) GO TO 300
323      DO 169 NX=1,NFHBX2
324      NS=NBDTIM*NX
325      READ(IN,*) IFHBUN,CNSTM,IFHBPT
326      WRITE(IOUT,161) FHBXNM(5+NX),IFHBUN,CNSTM
327 161  FORMAT(/,1X,A16,
328     & 'FOR SPECIFIED-HEAD CELLS WILL BE READ ON UNIT',I4,'.',/,
329     &  ' VALUES WILL BE MULTIPLIED BY',G12.4,'.')
330      IF(IFHBPT.GT.0) THEN
331         WRITE(IOUT,62) FHBXNM(5+NX)
332         WRITE(IOUT,54) (DSH1,M=1,ND)
333      ENDIF
334      DO 168 N=1,NHED
335      READ(IFHBUN,*) (SBHED(NS+L,N),L=1,NBDTIM)
336      DO 166 L=1,NBDTIM
337      SBHED(NS+L,N)=SBHED(NS+L,N)*CNSTM
338 166  CONTINUE
339C
340C9B------IF DESIRED, PRINT TABLE OF AUXILIARY VARIABLE VALUES AT
341C9B------SPECIFIED-HEAD CELL LOCATIONS
342      IF(IFHBPT.GT.0) THEN
343         WRITE(IOUT,58) (IHDLOC(I,N),I=1,4),
344     &                        (SBHED(NS+L,N),L=1,NBDTIM)
345      ENDIF
346 168  CONTINUE
347 169  CONTINUE
348C
349C10-----RETURN
350 300  RETURN
351      END
352      SUBROUTINE FHB1AD(HNEW,HOLD,NCOL,NROW,NLAY,ISS,TOTIM,DELT,BDTIM,
353     & NBDTIM,FLWRAT,BDFV,BDHV,NFLW,SBHED,IHDLOC,NHED,
354     & NFHBX1,NFHBX2,IFHBD3,IFHBD4,IFHBD5,IFHBSS)
355C
356C------VERSION 0000 10JAN1997 FHB1AD
357C     ******************************************************************
358C     COMPUTE SPECIFIED FLOWS AND HEADS AT CURRENT TIME STEP
359C     ******************************************************************
360C
361C        SPECIFICATIONS:
362C     ------------------------------------------------------------------
363      COMMON /FHBCOM/ FHBXNM(10),FHBXWT(10)
364      CHARACTER*16 FHBXNM
365      DOUBLE PRECISION HNEW
366C
367      DIMENSION BDTIM(NBDTIM),FLWRAT(IFHBD3,NFLW),BDFV(IFHBD4,NFLW),
368     &          BDHV(NFHBX2,NHED),SBHED(IFHBD5,NHED),IHDLOC(4,NHED),
369     &          HNEW(NCOL,NROW,NLAY),HOLD(NCOL,NROW,NLAY)
370C     ------------------------------------------------------------------
371C
372C1------IF THIS IS A STEADY-STATE SIMULATION OR A TRANSIENT SIMULATION
373C1------WITH CONSTANT SPECIFIED FLOWS AND HEADS, SET VALUES AND RETURN
374      IF((ISS.NE.0.AND.IFHBSS.EQ.0).OR.NBDTIM.EQ.1) THEN
375         IF(NFLW.LT.1) GO TO 6
376         DO 5 NF=1,NFLW
377         BDFV(1,NF)=FLWRAT(1,NF)
378         IF(NFHBX1.LT.1) GO TO 5
379         DO 4 NX=1,NFHBX1
380         N1=2+NX
381         N2=1+NX*NBDTIM
382         BDFV(N1,NF)=FLWRAT(N2,NF)
383 4       CONTINUE
384 5       CONTINUE
385 6       IF(NHED.LT.1) RETURN
386         DO 10 NH=1,NHED
387         K=IHDLOC(1,NH)
388         I=IHDLOC(2,NH)
389         J=IHDLOC(3,NH)
390         HNEW(J,I,K)=SBHED(1,NH)
391         HOLD(J,I,K)=SBHED(1,NH)
392         IF(NFHBX2.LT.1) GO TO 10
393         DO 8 NX=1,NFHBX2
394         N2=1+NX*NBDTIM
395         BDHV(NX,NF)=SBHED(N2,NF)
396 8       CONTINUE
397 10      CONTINUE
398         RETURN
399      ENDIF
400C
401C2------FIND ARRAY INDICES OF TIMES AROUND TIME AT START OF CURRENT
402C2------TIME STEP
403      IF(NFLW.LT.1) GO TO 200
404      T2=TOTIM
405      T1=TOTIM-DELT
406      DO 20 L=2,NBDTIM
407      IF(T1.LE.BDTIM(L)) THEN
408         IB1=L-1
409         IB2=L
410         GO TO 40
411      ENDIF
412 20   CONTINUE
413      IB1=NBDTIM-1
414      IB2=NBDTIM
415C
416C3------COMPUTE FACTOR FOR INTERPOLATION OR EXTRAPOLATION OF FLOW AT
417C3------START OF CURRENT TIME STEP
418 40   QFACT1=(T1-BDTIM(IB1))/(BDTIM(IB2)-BDTIM(IB1))
419C
420C4------FIND ARRAY INDICES OF TIMES AROUND TIME AT END OF CURRENT
421C4------TIME STEP
422      DO 60 L=IB2,NBDTIM
423      IF(T2.LE.BDTIM(L)) THEN
424         IB3=L-1
425         IB4=L
426         GO TO 70
427      ENDIF
428 60   CONTINUE
429      IB4=NBDTIM
430      IB3=NBDTIM-1
431C
432C5------COMPUTE FACTOR FOR INTERPOLATION OR EXTRAPOLATION OF FLOW AT
433C5------END OF CURRENT TIME STEP
434 70   QFACT2=(T2-BDTIM(IB3))/(BDTIM(IB4)-BDTIM(IB3))
435C
436C6------COMPUTE SPECIFIED FLOW RATES FOR THIS TIME STEP
437      NPI=IB4-IB2
438      DO 90 NF=1,NFLW
439      QA=FLWRAT(IB1,NF)
440      QB=FLWRAT(IB2,NF)
441      QC=FLWRAT(IB3,NF)
442      QD=FLWRAT(IB4,NF)
443      Q1=(QA+QFACT1*(QB-QA))
444      Q2=(QC+QFACT2*(QD-QC))
445      IF(NPI.EQ.0) THEN
446         BDFV(1,NF)=0.5*(Q1+Q2)
447      ELSE
448         TP=T1
449         QP=Q1
450         SUM1=0.0
451         DO 80 NI=IB2,IB3
452         QN=FLWRAT(NI,NF)
453         DDT=BDTIM(NI)-TP
454         SUM1=SUM1+DDT*0.5*(QN+QP)
455         TP=BDTIM(NI)
456         QP=QN
457 80      CONTINUE
458         DDT=T2-TP
459         SUM1=SUM1+DDT*0.5*(Q2+QP)
460         BDFV(1,NF)=SUM1/DELT
461      ENDIF
462 90   CONTINUE
463C
464C7-----COMPUTE VALUES OF AUXILIARY VARIABLES FOR SPECIFIED-FLOW
465C7-----CELLS FOR CURRENT TIME STEP
466      IF(NFHBX1.LT.1) GO TO 200
467      DO 190 NX=1,NFHBX1
468      N1=2+NX
469      N2=NX*NBDTIM
470      TT=TOTIM-(1.-FHBXWT(NX))*DELT
471      DO 120 L=2,NBDTIM
472      IF(TT.LE.BDTIM(L)) THEN
473         IB1=L-1
474         IB2=L
475         GO TO 140
476      ENDIF
477 120  CONTINUE
478      IB1=NBDTIM-1
479      IB2=NBDTIM
480 140  XFACT=(TT-BDTIM(IB1))/(BDTIM(IB2)-BDTIM(IB1))
481      DO 150 NF=1,NFLW
482      XX=FLWRAT(N2+IB1,NF)+XFACT*(FLWRAT(N2+IB2,NF)-FLWRAT(N2+IB1,NF))
483      BDFV(N1,NF)=XX
484 150  CONTINUE
485 190  CONTINUE
486C8------FIND ARRAY INDICES OF TIMES AROUND TIME AT END OF CURRENT
487C8------TIME STEP, COMPUTE FACTOR OF INTERPOLATION OR EXTRAPOLATION
488C8------OF HEAD
489 200  IF(NHED.LT.1) RETURN
490      TT=TOTIM
491      DO 220 L=2,NBDTIM
492      IF(TT.LE.BDTIM(L)) THEN
493         IB1=L-1
494         IB2=L
495         GO TO 240
496      ENDIF
497 220  CONTINUE
498      IB1=NBDTIM-1
499      IB2=NBDTIM
500 240  HFACT=(TT-BDTIM(IB1))/(BDTIM(IB2)-BDTIM(IB1))
501C
502C9------AT EACH SPECIFIED-HEAD LOCATION, INTERPOLATE OR EXTRAPOLATE
503C9------HEAD. SET HNEW AND HOLD EQUAL TO COMPUTED HEAD
504      DO 250 NH=1,NHED
505      K=IHDLOC(1,NH)
506      I=IHDLOC(2,NH)
507      J=IHDLOC(3,NH)
508      HH=SBHED(IB1,NH)+HFACT*(SBHED(IB2,NH)-SBHED(IB1,NH))
509      HNEW(J,I,K)=HH
510      HOLD(J,I,K)=HH
511 250  CONTINUE
512C
513C10----COMPUTE VALUES OF AUXILIARY VARIABLES FOR FOR SPECIFIED-HEAD
514C10----CELLS FOR CURRENT TIME STEP
515      IF(NFHBX2.LT.1) RETURN
516      DO 390 NX=1,NFHBX2
517      N2=NX*NBDTIM
518      TT=TOTIM-(1.-FHBXWT(5+NX))*DELT
519      DO 320 L=2,NBDTIM
520      IF(TT.LE.BDTIM(L)) THEN
521         IB1=L-1
522         IB2=L
523         GO TO 340
524      ENDIF
525 320  CONTINUE
526      IB1=NBDTIM-1
527      IB2=NBDTIM
528 340  XFACT=(TT-BDTIM(IB1))/(BDTIM(IB2)-BDTIM(IB1))
529      DO 350 NF=1,NHED
530      XX=SBHED(N2+IB1,NF)+XFACT*(SBHED(N2+IB2,NF)-SBHED(N2+IB1,NF))
531      BDHV(NX,NF)=XX
532 350  CONTINUE
533 390  CONTINUE
534C
535C11------RETURN
536      RETURN
537      END
538      SUBROUTINE FHB1FM(RHS,IBOUND,IFLLOC,BDFV,NFLW,NCOL,NROW,NLAY,
539     &                   IFHBD4)
540C
541C-----VERSION 0000 10JAN1997 FHB1FM
542C
543C     ******************************************************************
544C     SUBTRACT SPECIFIED Q FROM RHS
545C     ******************************************************************
546C
547C        SPECIFICATIONS:
548C     ------------------------------------------------------------------
549      DIMENSION RHS(NCOL,NROW,NLAY),IBOUND(NCOL,NROW,NLAY),
550     1            IFLLOC(4,NFLW),BDFV(IFHBD4,NFLW)
551C     ------------------------------------------------------------------
552C
553C1------PROCESS EACH SPECIFIED-FLOW LOCATION IN THE LIST.
554      IF(NFLW.LE.0) RETURN
555      DO 100 L=1,NFLW
556      IR=IFLLOC(2,L)
557      IC=IFLLOC(3,L)
558      IL=IFLLOC(1,L)
559      Q=BDFV(1,L)
560C
561C1A-----IF THE CELL IS INACTIVE THEN BYPASS PROCESSING.
562      IF(IBOUND(IC,IR,IL).LE.0) GO TO 100
563C
564C1B-----IF THE CELL IS VARIABLE HEAD THEN SUBTRACT Q FROM
565C       THE RHS ACCUMULATOR.
566      RHS(IC,IR,IL)=RHS(IC,IR,IL)-Q
567  100 CONTINUE
568C
569C2------RETURN
570      RETURN
571      END
572      SUBROUTINE FHB1BD(IFLLOC,BDFV,NFLW,VBNM,VBVL,MSUM,IBOUND,DELT,
573     &      NCOL,NROW,NLAY,KSTP,KPER,IFHBCB,ICBCFL,BUFF,IOUT,IFHBD4)
574C
575C-----VERSION 0000 10JAN1997 FHB1BD
576C     ******************************************************************
577C     CALCULATE VOLUMETRIC BUDGET FOR SPECIFIED FLOWS
578C     ******************************************************************
579C
580C        SPECIFICATIONS:
581C     ------------------------------------------------------------------
582      CHARACTER*16 VBNM(MSUM),TEXT
583      DOUBLE PRECISION RATIN,RATOUT,QQ
584      DIMENSION VBVL(4,MSUM),IBOUND(NCOL,NROW,NLAY),
585     1          BUFF(NCOL,NROW,NLAY),IFLLOC(4,NFLW),BDFV(IFHBD4,NFLW)
586C
587      DATA TEXT/' SPECIFIED FLOWS'/
588C     ------------------------------------------------------------------
589C
590C1------CLEAR RATIN AND RATOUT ACCUMULATORS.
591      ZERO=0.
592      RATIN=ZERO
593      RATOUT=ZERO
594      IBD=0
595      IF(IFHBCB.LT.0 .AND. ICBCFL.NE.0) IBD=-1
596      IF(IFHBCB.GT.0) IBD=ICBCFL
597C
598C2A----IF CELL-BY-CELL FLOWS WILL BE SAVED AS A LIST, WRITE HEADER.
599      IF(IBD.EQ.2) CALL UBDSV2(KSTP,KPER,TEXT,IFHBCB,NCOL,NROW,NLAY,
600     1          NFLW,IOUT,DELT,PERTIM,TOTIM,IBOUND)
601C
602C2B----CLEAR THE BUFFER.
603      DO 50 IL=1,NLAY
604      DO 50 IR=1,NROW
605      DO 50 IC=1,NCOL
606      BUFF(IC,IR,IL)=ZERO
607   50 CONTINUE
608C
609C3A----IF THERE ARE NO SPECIFIED-FLOW CELLS, DO NOT ACCUMULATE FLOW
610      IF(NFLW.EQ.0) GO TO 200
611C
612C3B-----PROCESS SPECIFIED-FLOW CELLS ONE AT A TIME.
613   60 DO 100 L=1,NFLW
614C
615C3C-----GET LAYER, ROW, AND COLUMN NUMBERS
616      IR=IFLLOC(2,L)
617      IC=IFLLOC(3,L)
618      IL=IFLLOC(1,L)
619      Q=ZERO
620C
621C3D-----IF THE CELL IS NO-FLOW OR CONSTANT-HEAD, IGNORE IT.
622      IF(IBOUND(IC,IR,IL).LE.0)GO TO 97
623C
624C4A-----GET FLOW RATE FROM SPECIFIED-FLOW LIST
625      Q=BDFV(1,L)
626      QQ=Q
627C
628C4B-----PRINT THE INDIVIDUAL RATES IF REQUESTED(IFHBCB<0).
629      IF(IBD.LT.0) THEN
630        WRITE(IOUT,900) TEXT,KPER,KSTP,L,IL,IR,IC,Q
631  900   FORMAT(1H0,4A4,'   PERIOD',I3,'   STEP',I3,' SEQ NO',I4,
632     1    '   LAYER',I3,'   ROW ',I4,'   COL',I4,'   RATE',G15.7)
633      ENDIF
634C
635C4C-----ADD FLOW RATE TO BUFFER.
636      BUFF(IC,IR,IL)=BUFF(IC,IR,IL)+Q
637C
638C5A-----SEE IF FLOW RATE IS NEGATIVE, ZERO, OR POSITIVE.
639      IF(Q) 90,97,80
640C
641C5B-----FLOW RATE IS POSITIVE (RECHARGE). ADD IT TO RATIN.
642   80 RATIN=RATIN+QQ
643      GO TO 97
644C
645C5C-----FLOW RATE IS NEGATIVE(DISCHARGE). ADD IT TO RATOUT.
646   90 RATOUT=RATOUT-QQ
647C
648C6------IF CELL-BY-CELL FLOWS ARE BEING SAVED AS A LIST, WRITE FLOW.
649C6------RETURN THE ACTUAL FLOW IN THE BDFV ARRAY.
650   97 IF(IBD.EQ.2) CALL UBDSVA(IFHBCB,NCOL,NROW,IC,IR,IL,Q,IBOUND,NLAY)
651      BDFV(2,L)=Q
652  100 CONTINUE
653C
654C7------IF CELL-BY-CELL FLOWS WILL BE SAVED AS A 3-D ARRAY,
655C7------CALL UBUDSV TO SAVE THEM
656      IF(IBD.EQ.1) CALL UBUDSV(KSTP,KPER,TEXT,IFHBCB,BUFF,NCOL,NROW,
657     1                          NLAY,IOUT)
658C
659C8------MOVE RATES, VOLUMES & LABELS INTO ARRAYS FOR PRINTING.
660 200  RIN=RATIN
661      ROUT=RATOUT
662      VBVL(3,MSUM)=RIN
663      VBVL(4,MSUM)=ROUT
664      VBVL(1,MSUM)=VBVL(1,MSUM)+RIN*DELT
665      VBVL(2,MSUM)=VBVL(2,MSUM)+ROUT*DELT
666      VBNM(MSUM)=TEXT
667C
668C9------INCREMENT BUDGET TERM COUNTER (MSUM).
669      MSUM=MSUM+1
670C
671C10-----RETURN
672      RETURN
673      END
674
675
676
677
678