33 character(len=LENFTYPE) ::
ftype =
'CDB'
34 character(len=16) ::
text =
' CDB'
38 integer(I4B),
dimension(:),
pointer,
contiguous :: idcxs => null()
39 real(dp),
dimension(:),
pointer,
contiguous :: width => null()
40 real(dp),
pointer :: gravconv => null()
73 subroutine cdb_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, &
74 mempath, dis, cxs, lengthconv, timeconv)
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
86 real(dp),
intent(in) :: lengthconv
87 real(dp),
intent(in) :: timeconv
96 call packobj%set_names(ibcnum, namemodel, pakname,
ftype, mempath)
100 call cdbobj%allocate_scalars()
103 call packobj%pack_initialize()
105 packobj%inunit = inunit
108 packobj%ibcnum = ibcnum
120 cdbobj%gravconv =
dgravity * lengthconv * timeconv**2
139 call this%BndExtType%allocate_scalars()
142 call mem_allocate(this%gravconv,
'GRAVCONV', this%memoryPath)
145 this%gravconv =
dzero
161 integer(I4B),
dimension(:),
pointer,
contiguous,
optional :: nodelist
162 real(DP),
dimension(:, :),
pointer,
contiguous,
optional :: auxvar
166 call this%BndExtType%allocate_arrays(nodelist, auxvar)
169 call mem_setptr(this%idcxs,
'IDCXS', this%input_mempath)
170 call mem_setptr(this%width,
'WIDTH', this%input_mempath)
173 call mem_checkin(this%idcxs,
'IDCXS', this%memoryPath, &
174 'IDCXS', this%input_mempath)
175 call mem_checkin(this%width,
'WIDTH', this%memoryPath, &
176 'WIDTH', this%input_mempath)
194 call this%BndExtType%bnd_da()
224 call this%BndExtType%source_options()
230 call this%log_cdb_options(found)
248 write (this%iout,
'(/1x,a)')
'PROCESSING '//trim(adjustl(this%text)) &
256 write (this%iout,
'(1x,a)') &
257 'END OF '//trim(adjustl(this%text))//
' OPTIONS'
272 if (this%iper /=
kper)
return
275 call this%BndExtType%bnd_rp()
278 if (this%iprpak /= 0)
then
279 call this%write_list()
298 integer(I4B) :: i, node
306 if (this%nbound == 0)
return
309 do i = 1, this%nbound
311 node = this%nodelist(i)
312 if (this%ibound(node) <= 0)
then
319 depth = this%xnew(node) - this%dis%bot(node)
322 q = this%qcalc(i, depth)
326 qeps = this%qcalc(i, depth + eps)
329 derv = (qeps - q) / eps
333 this%rhs(i) = -q + derv * this%xnew(node)
342 function qcalc(this, i, depth)
result(q)
346 integer(I4B),
intent(in) :: i
347 real(dp),
intent(in) :: depth
351 integer(I4B) :: idcxs
356 idcxs = this%idcxs(i)
357 width = this%width(i)
358 a = this%cxs%get_area(idcxs, width, depth)
359 r = this%cxs%get_hydraulic_radius(idcxs, width, depth, area=a)
361 q = this%gravconv * a**
dtwo * r
377 subroutine cdb_fc(this, rhs, ia, idxglo, matrix_sln)
380 real(DP),
dimension(:),
intent(inout) :: rhs
381 integer(I4B),
dimension(:),
intent(in) :: ia
382 integer(I4B),
dimension(:),
intent(in) :: idxglo
390 if (this%imover == 1)
then
391 call this%pakmvrobj%fc()
395 do i = 1, this%nbound
397 rhs(n) = rhs(n) + this%rhs(i)
399 call matrix_sln%add_value_pos(idxglo(ipos), this%hcof(i))
403 if (this%imover == 1 .and. this%rhs(i) >
dzero)
then
404 call this%pakmvrobj%accumulate_qformvr(i, this%rhs(i))
423 this%listlabel = trim(this%filtyp)//
' NO.'
424 if (this%dis%ndim == 3)
then
425 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
426 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'ROW'
427 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'COL'
428 elseif (this%dis%ndim == 2)
then
429 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
430 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'CELL2D'
432 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'NODE'
434 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'FLOW RATE'
435 if (this%inamedbound == 1)
then
436 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'BOUNDARY NAME'
476 call this%obs%StoreObsType(
'cdb', .true., indx)
481 call this%obs%StoreObsType(
'to-mvr', .true., indx)
504 call this%obs%obs_bd_clear()
507 do i = 1, this%obs%npakobs
508 obsrv => this%obs%pakobs(i)%obsrv
509 if (obsrv%BndFound)
then
510 do n = 1, obsrv%indxbnds_count
512 jj = obsrv%indxbnds(n)
513 select case (obsrv%ObsTypeId)
515 if (this%imover == 1)
then
516 v = this%pakmvrobj%get_qtomvr(jj)
524 errmsg =
'Unrecognized observation type: '//trim(obsrv%ObsTypeId)
527 call this%obs%SaveOneSimval(obsrv, v)
530 call this%obs%SaveOneSimval(obsrv,
dnodata)
550 integer(I4B) :: i, nlinks
554 nlinks = this%TsManager%boundtslinks%Count()
557 if (
associated(tslink))
then
558 if (tslink%JCol == 1)
then
579 integer(I4B),
intent(in) :: col
580 integer(I4B),
intent(in) :: row
586 bndval = this%idcxs(row)
588 bndval = this%width(row)
590 errmsg =
'Programming error. CDB bound value requested column '&
591 &
'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.
real(dp), parameter dnodata
real no data constant
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 dgravity
real constant gravitational acceleration (m/(s s))
real(dp), parameter dzero
real constant zero
real(dp), parameter dprec
real constant machine precision
real(dp), parameter dtwo
real constant 2
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 CDB package methods.
subroutine cdb_allocate_arrays(this, nodelist, auxvar)
@ brief Allocate arrays
subroutine cdb_rp(this)
@ brief SWF read and prepare
subroutine, public cdb_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, mempath, dis, cxs, lengthconv, timeconv)
@ brief Create a new package object
real(dp) function qcalc(this, i, depth)
Calculate critical depth boundary flow.
logical function cdb_obs_supported(this)
Determine if observations are supported.
subroutine log_cdb_options(this, found)
@ brief Log SWF specific package options
subroutine cdb_bd_obs(this)
Save observations for the package.
character(len=16) text
package flow text string
character(len=lenftype) ftype
package ftype
subroutine cdb_rp_ts(this)
Assign time series links for the package.
subroutine cdb_da(this)
@ brief Deallocate package memory
subroutine cdb_allocate_scalars(this)
@ brief Allocate scalars
subroutine cdb_options(this)
@ brief Source additional options for package
subroutine cdb_cf(this)
@ brief Formulate the package hcof and rhs terms.
subroutine cdb_df_obs(this)
Define the observation types available in the package.
subroutine cdb_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
real(dp) function cdb_bound_value(this, col, row)
@ brief Return a bound value
integer(i4b), pointer, public kper
current stress period number
type(timeserieslinktype) function, pointer, public gettimeserieslinkfromlist(list, indx)
Get time series link from a list.