31 character(len=LENFTYPE) ::
ftype =
'FLW'
32 character(len=16) ::
text =
' FLW'
35 real(dp),
dimension(:),
pointer,
contiguous :: q => null()
63 subroutine flw_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, &
66 class(
bndtype),
pointer :: packobj
67 integer(I4B),
intent(in) :: id
68 integer(I4B),
intent(in) :: ibcnum
69 integer(I4B),
intent(in) :: inunit
70 integer(I4B),
intent(in) :: iout
71 character(len=*),
intent(in) :: namemodel
72 character(len=*),
intent(in) :: pakname
73 character(len=*),
intent(in) :: mempath
82 call packobj%set_names(ibcnum, namemodel, pakname,
ftype, mempath)
86 call flwobj%allocate_scalars()
89 call packobj%pack_initialize()
91 packobj%inunit = inunit
94 packobj%ibcnum = ibcnum
97 packobj%ictMemPath =
''
116 call this%BndExtType%allocate_scalars()
138 integer(I4B),
dimension(:),
pointer,
contiguous,
optional :: nodelist
139 real(DP),
dimension(:, :),
pointer,
contiguous,
optional :: auxvar
143 call this%BndExtType%allocate_arrays(nodelist, auxvar)
146 call mem_setptr(this%q,
'Q', this%input_mempath)
150 'Q', this%input_mempath)
168 call this%BndExtType%bnd_da()
194 call this%BndExtType%source_options()
200 call this%log_flw_options(found)
218 write (this%iout,
'(/1x,a)')
'PROCESSING '//trim(adjustl(this%text)) &
226 write (this%iout,
'(1x,a)') &
227 'END OF '//trim(adjustl(this%text))//
' OPTIONS'
242 if (this%iper /=
kper)
return
245 call this%BndExtType%bnd_rp()
248 if (this%iprpak /= 0)
then
249 call this%write_list()
266 integer(I4B) :: i, node
270 if (this%nbound == 0)
return
273 do i = 1, this%nbound
274 node = this%nodelist(i)
276 if (this%ibound(node) <= 0)
then
293 subroutine flw_fc(this, rhs, ia, idxglo, matrix_sln)
296 real(DP),
dimension(:),
intent(inout) :: rhs
297 integer(I4B),
dimension(:),
intent(in) :: ia
298 integer(I4B),
dimension(:),
intent(in) :: idxglo
306 if (this%imover == 1)
then
307 call this%pakmvrobj%fc()
311 do i = 1, this%nbound
313 rhs(n) = rhs(n) + this%rhs(i)
315 call matrix_sln%add_value_pos(idxglo(ipos), this%hcof(i))
319 if (this%imover == 1 .and. this%rhs(i) >
dzero)
then
320 call this%pakmvrobj%accumulate_qformvr(i, this%rhs(i))
339 this%listlabel = trim(this%filtyp)//
' NO.'
340 if (this%dis%ndim == 3)
then
341 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
342 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'ROW'
343 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'COL'
344 elseif (this%dis%ndim == 2)
then
345 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
346 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'CELL2D'
348 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'NODE'
350 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'FLOW RATE'
351 if (this%inamedbound == 1)
then
352 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'BOUNDARY NAME'
392 call this%obs%StoreObsType(
'flw', .true., indx)
397 call this%obs%StoreObsType(
'to-mvr', .true., indx)
420 call this%obs%obs_bd_clear()
423 do i = 1, this%obs%npakobs
424 obsrv => this%obs%pakobs(i)%obsrv
425 if (obsrv%BndFound)
then
426 do n = 1, obsrv%indxbnds_count
428 jj = obsrv%indxbnds(n)
429 select case (obsrv%ObsTypeId)
431 if (this%imover == 1)
then
432 v = this%pakmvrobj%get_qtomvr(jj)
440 errmsg =
'Unrecognized observation type: '//trim(obsrv%ObsTypeId)
443 call this%obs%SaveOneSimval(obsrv, v)
446 call this%obs%SaveOneSimval(obsrv,
dnodata)
466 integer(I4B) :: i, nlinks
470 nlinks = this%TsManager%boundtslinks%Count()
473 if (
associated(tslink))
then
474 if (tslink%JCol == 1)
then
489 integer(I4B),
intent(in) :: row
493 if (this%iauxmultcol > 0)
then
494 q = this%q(row) * this%auxvar(this%iauxmultcol, row)
514 integer(I4B),
intent(in) :: col
515 integer(I4B),
intent(in) :: row
521 bndval = this%q_mult(row)
523 errmsg =
'Programming error. FLW bound value requested column '&
524 &
'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
This module contains the FLW package methods.
subroutine flw_df_obs(this)
Define the observation types available in the package.
character(len=lenftype) ftype
package ftype
subroutine flw_cf(this)
@ brief Formulate the package hcof and rhs terms.
subroutine, public flw_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, mempath)
@ brief Create a new package object
subroutine flw_da(this)
@ brief Deallocate package memory
real(dp) function flw_bound_value(this, col, row)
@ brief Return a bound value
logical function flw_obs_supported(this)
Determine if observations are supported.
subroutine flw_allocate_scalars(this)
@ brief Allocate scalars
character(len=16) text
package flow text string
subroutine flw_rp(this)
@ brief SWF read and prepare
subroutine flw_bd_obs(this)
Save observations for the package.
subroutine define_listlabel(this)
@ brief Define the list label for the package
subroutine log_flw_options(this, found)
@ brief Log SWF specific package options
subroutine flw_options(this)
@ brief Source additional options for package
real(dp) function q_mult(this, row)
subroutine flw_allocate_arrays(this, nodelist, auxvar)
@ brief Allocate arrays
subroutine flw_fc(this, rhs, ia, idxglo, matrix_sln)
@ brief Copy hcof and rhs terms into solution.
subroutine flw_rp_ts(this)
Assign time series links for the package.
integer(i4b), pointer, public kper
current stress period number
type(timeserieslinktype) function, pointer, public gettimeserieslinkfromlist(list, indx)
Get time series link from a list.