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