36 character(len=LENFTYPE) ::
ftype =
'WEL'
37 character(len=16) ::
text =
' WEL'
40 real(dp),
dimension(:),
pointer,
contiguous :: q => null()
41 integer(I4B),
pointer :: iflowred => null()
42 real(dp),
pointer :: flowred => null()
43 integer(I4B),
pointer :: ioutafrcsv => null()
73 subroutine wel_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, &
76 class(
bndtype),
pointer :: packobj
77 integer(I4B),
intent(in) :: id
78 integer(I4B),
intent(in) :: ibcnum
79 integer(I4B),
intent(in) :: inunit
80 integer(I4B),
intent(in) :: iout
81 character(len=*),
intent(in) :: namemodel
82 character(len=*),
intent(in) :: pakname
83 character(len=*),
intent(in) :: mempath
85 type(
weltype),
pointer :: welobj
92 call packobj%set_names(ibcnum, namemodel, pakname,
ftype, mempath)
96 call welobj%allocate_scalars()
99 call packobj%pack_initialize()
101 packobj%inunit = inunit
104 packobj%ibcnum = ibcnum
123 call this%BndExtType%bnd_da()
148 call this%BndExtType%allocate_scalars()
151 call mem_allocate(this%iflowred,
'IFLOWRED', this%memoryPath)
152 call mem_allocate(this%flowred,
'FLOWRED', this%memoryPath)
153 call mem_allocate(this%ioutafrcsv,
'IOUTAFRCSV', this%memoryPath)
174 integer(I4B),
dimension(:),
pointer,
contiguous,
optional :: nodelist
175 real(DP),
dimension(:, :),
pointer,
contiguous,
optional :: auxvar
179 call this%BndExtType%allocate_arrays(nodelist, auxvar)
182 call mem_setptr(this%q,
'Q', this%input_mempath)
186 'Q', this%input_mempath)
203 class(
weltype),
intent(inout) :: this
205 character(len=LINELENGTH) :: fname
208 character(len=*),
parameter :: fmtflowred = &
209 &
"(4x, 'AUTOMATIC FLOW REDUCTION OF WELLS IMPLEMENTED.')"
210 character(len=*),
parameter :: fmtflowredv = &
211 &
"(4x, 'AUTOMATIC FLOW REDUCTION FRACTION (',g15.7,').')"
214 call this%BndExtType%source_options()
217 call mem_set_value(this%flowred,
'FLOWRED', this%input_mempath, found%flowred)
218 call mem_set_value(fname,
'AFRCSVFILE', this%input_mempath, found%afrcsvfile)
219 call mem_set_value(this%imover,
'MOVER', this%input_mempath, found%mover)
221 if (found%flowred)
then
225 if (this%flowred <=
dzero)
then
227 else if (this%flowred >
done)
then
232 if (found%afrcsvfile)
then
233 call this%wel_afr_csv_init(fname)
236 if (found%mover)
then
241 call this%log_wel_options(found)
253 class(
weltype),
intent(inout) :: this
257 character(len=*),
parameter :: fmtflowred = &
258 &
"(4x, 'AUTOMATIC FLOW REDUCTION OF WELLS IMPLEMENTED.')"
259 character(len=*),
parameter :: fmtflowredv = &
260 &
"(4x, 'AUTOMATIC FLOW REDUCTION FRACTION (',g15.7,').')"
263 write (this%iout,
'(/1x,a)')
'PROCESSING '//trim(adjustl(this%text)) &
266 if (found%flowred)
then
267 if (this%iflowred > 0) &
268 write (this%iout, fmtflowred)
269 write (this%iout, fmtflowredv) this%flowred
272 if (found%afrcsvfile)
then
276 if (found%mover)
then
277 write (this%iout,
'(4x,A)')
'MOVER OPTION ENABLED'
281 write (this%iout,
'(1x,a)') &
282 'END OF '//trim(adjustl(this%text))//
' OPTIONS'
294 class(
weltype),
intent(inout) :: this
297 if (this%iper /=
kper)
return
300 call this%BndExtType%bnd_rp()
303 if (this%iprpak /= 0)
then
304 call this%write_list()
321 integer(I4B) :: i, node, ict
329 if (this%nbound == 0)
return
332 do i = 1, this%nbound
333 node = this%nodelist(i)
335 if (this%ibound(node) <= 0)
then
340 if (this%iflowred /= 0 .and. q <
dzero)
then
341 ict = this%icelltype(node)
343 tp = this%dis%top(node)
344 bt = this%dis%bot(node)
346 tp = bt + this%flowred * thick
363 subroutine wel_fc(this, rhs, ia, idxglo, matrix_sln)
366 real(DP),
dimension(:),
intent(inout) :: rhs
367 integer(I4B),
dimension(:),
intent(in) :: ia
368 integer(I4B),
dimension(:),
intent(in) :: idxglo
376 if (this%imover == 1)
then
377 call this%pakmvrobj%fc()
381 do i = 1, this%nbound
383 rhs(n) = rhs(n) + this%rhs(i)
385 call matrix_sln%add_value_pos(idxglo(ipos), this%hcof(i))
389 if (this%imover == 1 .and. this%rhs(i) >
dzero)
then
390 call this%pakmvrobj%accumulate_qformvr(i, this%rhs(i))
404 subroutine wel_fn(this, rhs, ia, idxglo, matrix_sln)
407 real(DP),
dimension(:),
intent(inout) :: rhs
408 integer(I4B),
dimension(:),
intent(in) :: ia
409 integer(I4B),
dimension(:),
intent(in) :: idxglo
423 do i = 1, this%nbound
424 node = this%nodelist(i)
427 if (this%ibound(node) <= 0)
then
432 ict = this%icelltype(node)
433 if (this%iflowred /= 0 .and. ict /= 0)
then
438 tp = this%dis%top(node)
439 bt = this%dis%bot(node)
441 tp = bt + this%flowred * thick
443 drterm = drterm * this%q_mult(i)
445 call matrix_sln%add_value_pos(idxglo(ipos), drterm)
446 rhs(node) = rhs(node) + drterm * this%xnew(node)
458 class(
weltype),
intent(inout) :: this
459 character(len=*),
intent(in) :: fname
461 character(len=*),
parameter :: fmtafrcsv = &
462 "(4x, 'AUTO FLOW REDUCE INFORMATION WILL BE SAVED TO FILE: ', a, /4x, &
463 &'OPENED ON UNIT: ', I0)"
466 call openfile(this%ioutafrcsv, this%iout, fname,
'CSV', &
467 filstat_opt=
'REPLACE')
468 write (this%iout, fmtafrcsv) trim(adjustl(fname)), &
470 write (this%ioutafrcsv,
'(a)') &
471 'time,period,step,boundnumber,cellnumber,rate-requested,&
472 &rate-actual,wel-reduction'
481 class(
weltype),
intent(inout) :: this
484 integer(I4B) :: nodereduced
485 integer(I4B) :: nodeuser
488 do i = 1, this%nbound
489 nodereduced = this%nodelist(i)
492 if (this%ibound(nodereduced) <= 0)
then
495 v = this%q_mult(i) + this%rhs(i)
497 nodeuser = this%dis%get_nodeuser(nodereduced)
498 write (this%ioutafrcsv,
'(*(G0,:,","))') &
499 totim,
kper,
kstp, i, nodeuser, this%q_mult(i), this%simvals(i), v
512 class(
weltype),
intent(inout) :: this
515 this%listlabel = trim(this%filtyp)//
' NO.'
516 if (this%dis%ndim == 3)
then
517 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
518 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'ROW'
519 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'COL'
520 elseif (this%dis%ndim == 2)
then
521 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
522 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'CELL2D'
524 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'NODE'
526 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'STRESS RATE'
527 if (this%inamedbound == 1)
then
528 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'BOUNDARY NAME'
568 call this%obs%StoreObsType(
'wel', .true., indx)
573 call this%obs%StoreObsType(
'to-mvr', .true., indx)
578 call this%obs%StoreObsType(
'wel-reduction', .true., indx)
601 call this%obs%obs_bd_clear()
604 do i = 1, this%obs%npakobs
605 obsrv => this%obs%pakobs(i)%obsrv
606 if (obsrv%BndFound)
then
607 do n = 1, obsrv%indxbnds_count
609 jj = obsrv%indxbnds(n)
610 select case (obsrv%ObsTypeId)
612 if (this%imover == 1)
then
613 v = this%pakmvrobj%get_qtomvr(jj)
620 case (
'WEL-REDUCTION')
621 if (this%iflowred > 0)
then
622 v = this%q_mult(jj) + this%rhs(jj)
625 errmsg =
'Unrecognized observation type: '//trim(obsrv%ObsTypeId)
628 call this%obs%SaveOneSimval(obsrv, v)
631 call this%obs%SaveOneSimval(obsrv,
dnodata)
636 if (this%ioutafrcsv > 0)
then
637 call this%wel_afr_csv_write()
648 class(
weltype),
intent(inout) :: this
649 integer(I4B),
intent(in) :: row
653 if (this%iauxmultcol > 0)
then
654 q = this%q(row) * this%auxvar(this%iauxmultcol, row)
673 class(
weltype),
intent(inout) :: this
674 integer(I4B),
intent(in) :: col
675 integer(I4B),
intent(in) :: row
681 bndval = this%q_mult(row)
683 errmsg =
'Programming error. WEL bound value requested column '&
684 &
'outside range of ncolbnd (1).'
This module contains block parser methods.
This module contains the extended boundary package.
This module contains the base boundary package.
This module contains simulation constants.
integer(i4b), parameter linelength
maximum length of a standard line
real(dp), parameter dnodata
real no data constant
real(dp), parameter dem1
real constant 1e-1
integer(i4b), parameter lenftype
maximum length of a package type (DIS, WEL, OC, etc.)
real(dp), parameter dzero
real constant zero
real(dp), parameter done
real constant 1
This module defines variable data types.
character(len=lenmempath) function create_mem_path(component, subcomponent, context)
returns the path to the memory object
This module contains the derived types ObserveType and ObsDataType.
This module contains the derived type ObsType.
subroutine, public defaultobsidprocessor(obsrv, dis, inunitobs, iout)
@ brief Process IDstring provided for each observation
This module contains simulation methods.
subroutine, public store_error(msg, terminate)
Store an error message.
subroutine, public store_error_filename(filename, terminate)
Store the erroring file name.
This module contains simulation variables.
character(len=maxcharlen) errmsg
error message string
real(dp) function sqsaturationderivative(top, bot, x, c1, c2)
@ brief sQSaturationDerivative
real(dp) function sqsaturation(top, bot, x, c1, c2)
@ brief sQSaturation
real(dp), pointer, public totim
time relative to start of simulation
integer(i4b), pointer, public kstp
current time step number
integer(i4b), pointer, public kper
current stress period number
This module contains the WEL package methods.
subroutine, public wel_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, mempath)
@ brief Create a new package object
subroutine wel_allocate_scalars(this)
@ brief Allocate scalars
subroutine wel_fc(this, rhs, ia, idxglo, matrix_sln)
@ brief Copy hcof and rhs terms into solution.
subroutine define_listlabel(this)
@ brief Define the list label for the package
subroutine wel_rp(this)
@ brief WEL read and prepare
subroutine wel_allocate_arrays(this, nodelist, auxvar)
@ brief Allocate arrays
subroutine wel_options(this)
@ brief Source additional options for package
subroutine wel_afr_csv_write(this)
Write out auto flow reductions only when & where they occur.
real(dp) function q_mult(this, row)
real(dp) function wel_bound_value(this, col, row)
@ brief Return a bound value
subroutine wel_da(this)
@ brief Deallocate package memory
subroutine wel_fn(this, rhs, ia, idxglo, matrix_sln)
@ brief Add Newton-Raphson terms for package into solution.
subroutine wel_afr_csv_init(this, fname)
Initialize the auto flow reduce csv output file.
character(len=lenftype) ftype
package ftype
subroutine wel_cf(this)
@ brief Formulate the package hcof and rhs terms.
subroutine wel_bd_obs(this)
Save observations for the package.
logical function wel_obs_supported(this)
Determine if observations are supported.
subroutine log_wel_options(this, found)
@ brief Log WEL specific package options
character(len=16) text
package flow text string
subroutine wel_df_obs(this)
Define the observation types available in the package.