34 character(len=LENFTYPE) ::
ftype =
'ZDG'
35 character(len=16) ::
text =
' ZDG'
39 integer(I4B),
dimension(:),
pointer,
contiguous :: idcxs => null()
40 real(dp),
dimension(:),
pointer,
contiguous :: width => null()
41 real(dp),
dimension(:),
pointer,
contiguous :: slope => null()
42 real(dp),
dimension(:),
pointer,
contiguous :: rough => null()
43 real(dp),
pointer :: unitconv => null()
77 subroutine zdg_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, &
78 mempath, dis, cxs, unitconv)
80 class(
bndtype),
pointer :: packobj
81 integer(I4B),
intent(in) :: id
82 integer(I4B),
intent(in) :: ibcnum
83 integer(I4B),
intent(in) :: inunit
84 integer(I4B),
intent(in) :: iout
85 character(len=*),
intent(in) :: namemodel
86 character(len=*),
intent(in) :: pakname
87 character(len=*),
intent(in) :: mempath
90 real(dp),
intent(in) :: unitconv
99 call packobj%set_names(ibcnum, namemodel, pakname,
ftype, mempath)
103 call zdgobj%allocate_scalars()
106 call packobj%pack_initialize()
108 packobj%inunit = inunit
111 packobj%ibcnum = ibcnum
126 zdgobj%unitconv = unitconv
145 call this%BndExtType%allocate_scalars()
148 call mem_allocate(this%unitconv,
'UNITCONV', this%memoryPath)
151 this%unitconv =
dzero
167 integer(I4B),
dimension(:),
pointer,
contiguous,
optional :: nodelist
168 real(DP),
dimension(:, :),
pointer,
contiguous,
optional :: auxvar
172 call this%BndExtType%allocate_arrays(nodelist, auxvar)
175 call mem_setptr(this%idcxs,
'IDCXS', this%input_mempath)
176 call mem_setptr(this%width,
'WIDTH', this%input_mempath)
177 call mem_setptr(this%slope,
'SLOPE', this%input_mempath)
178 call mem_setptr(this%rough,
'ROUGH', this%input_mempath)
181 call mem_checkin(this%idcxs,
'IDCXS', this%memoryPath, &
182 'IDCXS', this%input_mempath)
183 call mem_checkin(this%width,
'WIDTH', this%memoryPath, &
184 'WIDTH', this%input_mempath)
185 call mem_checkin(this%slope,
'SLOPE', this%memoryPath, &
186 'SLOPE', this%input_mempath)
187 call mem_checkin(this%rough,
'ROUGH', this%memoryPath, &
188 'ROUGH', this%input_mempath)
206 call this%BndExtType%bnd_da()
238 call this%BndExtType%source_options()
244 call this%log_zdg_options(found)
262 write (this%iout,
'(/1x,a)')
'PROCESSING '//trim(adjustl(this%text)) &
270 write (this%iout,
'(1x,a)') &
271 'END OF '//trim(adjustl(this%text))//
' OPTIONS'
286 if (this%iper /=
kper)
return
289 call this%BndExtType%bnd_rp()
292 if (this%iprpak /= 0)
then
293 call this%write_list()
312 integer(I4B) :: i, node
315 real(DP) :: absdhdxsq
322 if (this%nbound == 0)
return
325 do i = 1, this%nbound
327 node = this%nodelist(i)
328 if (this%ibound(node) <= 0)
then
335 absdhdxsq = this%slope(i)**
dhalf
336 depth = this%xnew(node) - this%dis%bot(node)
340 cond = this%get_cond(i, depth, absdhdxsq, this%unitconv)
341 q = -cond * this%slope(i)
345 cond = this%get_cond(i, depth + eps, absdhdxsq, this%unitconv)
346 qeps = -cond * this%slope(i)
349 derv = (qeps - q) / eps
353 this%rhs(i) = -q + derv * this%xnew(node)
368 function get_cond(this, i, depth, absdhdxsq, unitconv)
result(c)
372 integer(I4B),
intent(in) :: i
373 real(dp),
intent(in) :: depth
374 real(dp),
intent(in) :: absdhdxsq
375 real(dp),
intent(in) :: unitconv
377 integer(I4B) :: idcxs
385 real(dp) :: conveyance
387 idcxs = this%idcxs(i)
388 width = this%width(i)
389 rough = this%rough(i)
390 slope = this%slope(i)
391 roughc = this%cxs%get_roughness(idcxs, width, depth, rough, &
393 a = this%cxs%get_area(idcxs, width, depth)
394 r = this%cxs%get_hydraulic_radius(idcxs, width, depth, area=a)
397 conveyance = this%cxs%get_conveyance(idcxs, width, depth, rough)
398 c = conveyance / absdhdxsq
408 subroutine zdg_fc(this, rhs, ia, idxglo, matrix_sln)
411 real(DP),
dimension(:),
intent(inout) :: rhs
412 integer(I4B),
dimension(:),
intent(in) :: ia
413 integer(I4B),
dimension(:),
intent(in) :: idxglo
421 if (this%imover == 1)
then
422 call this%pakmvrobj%fc()
426 do i = 1, this%nbound
428 rhs(n) = rhs(n) + this%rhs(i)
430 call matrix_sln%add_value_pos(idxglo(ipos), this%hcof(i))
434 if (this%imover == 1 .and. this%rhs(i) >
dzero)
then
435 call this%pakmvrobj%accumulate_qformvr(i, this%rhs(i))
454 this%listlabel = trim(this%filtyp)//
' NO.'
455 if (this%dis%ndim == 3)
then
456 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
457 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'ROW'
458 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'COL'
459 elseif (this%dis%ndim == 2)
then
460 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
461 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'CELL2D'
463 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'NODE'
465 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'FLOW RATE'
466 if (this%inamedbound == 1)
then
467 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'BOUNDARY NAME'
507 call this%obs%StoreObsType(
'zdg', .true., indx)
512 call this%obs%StoreObsType(
'to-mvr', .true., indx)
535 call this%obs%obs_bd_clear()
538 do i = 1, this%obs%npakobs
539 obsrv => this%obs%pakobs(i)%obsrv
540 if (obsrv%BndFound)
then
541 do n = 1, obsrv%indxbnds_count
543 jj = obsrv%indxbnds(n)
544 select case (obsrv%ObsTypeId)
546 if (this%imover == 1)
then
547 v = this%pakmvrobj%get_qtomvr(jj)
555 errmsg =
'Unrecognized observation type: '//trim(obsrv%ObsTypeId)
558 call this%obs%SaveOneSimval(obsrv, v)
561 call this%obs%SaveOneSimval(obsrv,
dnodata)
581 integer(I4B) :: i, nlinks
585 nlinks = this%TsManager%boundtslinks%Count()
588 if (
associated(tslink))
then
589 if (tslink%JCol == 1)
then
610 integer(I4B),
intent(in) :: col
611 integer(I4B),
intent(in) :: row
617 bndval = this%idcxs(row)
619 bndval = this%width(row)
621 bndval = this%slope(row)
623 bndval = this%rough(row)
625 errmsg =
'Programming error. ZDG bound value requested column '&
626 &
'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 dtwothirds
real constant 2/3
real(dp), parameter dnodata
real no data constant
real(dp), parameter dem1
real constant 1e-1
real(dp), parameter dhalf
real constant 1/2
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.
real(dp) function, public get_perturbation(x)
Calculate a numerical perturbation given the value of x.
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 ZDG package methods.
subroutine zdg_allocate_scalars(this)
@ brief Allocate scalars
subroutine define_listlabel(this)
@ brief Define the list label for the package
logical function zdg_obs_supported(this)
Determine if observations are supported.
subroutine zdg_cf(this)
@ brief Formulate the package hcof and rhs terms.
subroutine zdg_fc(this, rhs, ia, idxglo, matrix_sln)
@ brief Copy hcof and rhs terms into solution.
real(dp) function zdg_bound_value(this, col, row)
@ brief Return a bound value
subroutine log_zdg_options(this, found)
@ brief Log SWF specific package options
real(dp) function get_cond(this, i, depth, absdhdxsq, unitconv)
Calculate conductance-like term.
subroutine, public zdg_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, mempath, dis, cxs, unitconv)
@ brief Create a new package object
subroutine zdg_rp_ts(this)
Assign time series links for the package.
character(len=16) text
package flow text string
subroutine zdg_da(this)
@ brief Deallocate package memory
subroutine zdg_options(this)
@ brief Source additional options for package
subroutine zdg_rp(this)
@ brief SWF read and prepare
character(len=lenftype) ftype
package ftype
subroutine zdg_allocate_arrays(this, nodelist, auxvar)
@ brief Allocate arrays
subroutine zdg_bd_obs(this)
Save observations for the package.
subroutine zdg_df_obs(this)
Define the observation types available in 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.