18 character(len=LENFTYPE) ::
ftype =
'GHB'
19 character(len=LENPACKAGENAME) ::
text =
' GHB'
22 real(dp),
dimension(:),
pointer,
contiguous :: bhead => null()
23 real(dp),
dimension(:),
pointer,
contiguous :: cond => null()
46 subroutine ghb_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, &
49 class(
bndtype),
pointer :: packobj
50 integer(I4B),
intent(in) :: id
51 integer(I4B),
intent(in) :: ibcnum
52 integer(I4B),
intent(in) :: inunit
53 integer(I4B),
intent(in) :: iout
54 character(len=*),
intent(in) :: namemodel
55 character(len=*),
intent(in) :: pakname
56 character(len=*),
intent(in) :: mempath
58 type(
ghbtype),
pointer :: ghbobj
65 call packobj%set_names(ibcnum, namemodel, pakname,
ftype, mempath)
69 call packobj%allocate_scalars()
72 call packobj%pack_initialize()
74 packobj%inunit = inunit
77 packobj%ibcnum = ibcnum
93 call this%BndExtType%bnd_da()
111 class(
ghbtype),
intent(inout) :: this
116 call this%BndExtType%source_options()
119 call mem_set_value(this%imover,
'MOVER', this%input_mempath, found%mover)
122 call this%log_ghb_options(found)
134 class(
ghbtype),
intent(inout) :: this
138 write (this%iout,
'(/1x,a)')
'PROCESSING '//trim(adjustl(this%text)) &
141 if (found%mover)
then
142 write (this%iout,
'(4x,A)')
'MOVER OPTION ENABLED'
146 write (this%iout,
'(1x,a)') &
147 'END OF '//trim(adjustl(this%text))//
' OPTIONS'
160 integer(I4B),
dimension(:),
pointer,
contiguous,
optional :: nodelist
161 real(DP),
dimension(:, :),
pointer,
contiguous,
optional :: auxvar
164 call this%BndExtType%allocate_arrays(nodelist, auxvar)
167 call mem_setptr(this%bhead,
'BHEAD', this%input_mempath)
168 call mem_setptr(this%cond,
'COND', this%input_mempath)
171 call mem_checkin(this%bhead,
'BHEAD', this%memoryPath, &
172 'BHEAD', this%input_mempath)
173 call mem_checkin(this%cond,
'COND', this%memoryPath, &
174 'COND', this%input_mempath)
186 class(
ghbtype),
intent(inout) :: this
188 if (this%iper /=
kper)
return
191 call this%BndExtType%bnd_rp()
194 if (this%ivsc == 1)
then
195 call this%ghb_store_user_cond()
199 if (this%iprpak /= 0)
then
200 call this%write_list()
214 class(
ghbtype),
intent(inout) :: this
216 character(len=LINELENGTH) :: errmsg
221 character(len=*),
parameter :: fmtghberr = &
222 "('GHB BOUNDARY (',i0,') HEAD (',f10.3,') IS LESS THAN CELL &
223 &BOTTOM (',f10.3,')')"
224 character(len=*),
parameter :: fmtcondmulterr = &
225 "('GHB BOUNDARY (',i0,') CONDUCTANCE MULTIPLIER (',g10.3,') IS &
227 character(len=*),
parameter :: fmtconderr = &
228 "('GHB BOUNDARY (',i0,') CONDUCTANCE (',g10.3,') IS LESS THAN &
232 do i = 1, this%nbound
233 node = this%nodelist(i)
234 bt = this%dis%bot(node)
236 if (this%bhead(i) < bt .and. this%icelltype(node) /= 0)
then
237 write (errmsg, fmt=fmtghberr) i, this%bhead(i), bt
240 if (this%iauxmultcol > 0)
then
241 if (this%auxvar(this%iauxmultcol, i) < dzero)
then
242 write (errmsg, fmt=fmtcondmulterr) &
243 i, this%auxvar(this%iauxmultcol, i)
247 if (this%cond(i) < dzero)
then
248 write (errmsg, fmt=fmtconderr) i, this%cond(i)
270 integer(I4B) :: i, node
273 if (this%nbound .eq. 0)
return
276 do i = 1, this%nbound
277 node = this%nodelist(i)
278 if (this%ibound(node) .le. 0)
then
283 this%hcof(i) = -this%cond_mult(i)
284 this%rhs(i) = -this%cond_mult(i) * this%bhead(i)
293 subroutine ghb_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
301 integer(I4B) :: i, n, ipos
302 real(DP) :: cond, bhead, qghb
305 if (this%imover == 1)
then
306 call this%pakmvrobj%fc()
310 do i = 1, this%nbound
312 rhs(n) = rhs(n) + this%rhs(i)
314 call matrix_sln%add_value_pos(idxglo(ipos), this%hcof(i))
318 bhead = this%bhead(i)
319 if (this%imover == 1 .and. this%xnew(n) > bhead)
then
320 cond = this%cond_mult(i)
321 qghb = cond * (this%xnew(n) - bhead)
322 call this%pakmvrobj%accumulate_qformvr(i, qghb)
335 class(
ghbtype),
intent(inout) :: this
338 this%listlabel = trim(this%filtyp)//
' NO.'
339 if (this%dis%ndim == 3)
then
340 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
341 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'ROW'
342 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'COL'
343 elseif (this%dis%ndim == 2)
then
344 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'LAYER'
345 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'CELL2D'
347 write (this%listlabel,
'(a, a7)') trim(this%listlabel),
'NODE'
349 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'STAGE'
350 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'CONDUCTANCE'
351 if (this%inamedbound == 1)
then
352 write (this%listlabel,
'(a, a16)') trim(this%listlabel),
'BOUNDARY NAME'
387 call this%obs%StoreObsType(
'ghb', .true., indx)
392 call this%obs%StoreObsType(
'to-mvr', .true., indx)
403 class(
ghbtype),
intent(inout) :: this
408 do n = 1, this%nbound
409 this%condinput(n) = this%cond_mult(n)
422 class(
ghbtype),
intent(inout) :: this
423 integer(I4B),
intent(in) :: row
427 if (this%iauxmultcol > 0)
then
428 cond = this%cond(row) * this%auxvar(this%iauxmultcol, row)
430 cond = this%cond(row)
443 class(
ghbtype),
intent(inout) :: this
444 integer(I4B),
intent(in) :: col
445 integer(I4B),
intent(in) :: row
451 bndval = this%bhead(row)
453 bndval = this%cond_mult(row)
455 errmsg =
'Programming error. GHB bound value requested column '&
456 &
'outside range of ncolbnd (2).'
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
integer(i4b), parameter lenpackagename
maximum length of the package name
integer(i4b), parameter lenftype
maximum length of a package type (DIS, WEL, OC, etc.)
real(dp), parameter dzero
real constant zero
subroutine ghb_cf(this)
Formulate the HCOF and RHS terms.
subroutine ghb_allocate_arrays(this, nodelist, auxvar)
Allocate arrays.
subroutine log_ghb_options(this, found)
Log options specific to GhbType.
subroutine ghb_ck(this)
Check ghb boundary condition data.
subroutine ghb_da(this)
Deallocate memory.
subroutine ghb_fc(this, rhs, ia, idxglo, matrix_sln)
Copy rhs and hcof into solution rhs and amat.
subroutine ghb_df_obs(this)
Store observation type supported by GHB package.
character(len=lenpackagename) text
real(dp) function ghb_bound_value(this, col, row)
Return requested boundary value.
real(dp) function cond_mult(this, row)
Apply multiplier to GHB conductance if option AUXMULTCOL is used.
character(len=lenftype) ftype
subroutine ghb_store_user_cond(this)
Store user-specified conductance for GHB boundary type.
logical function ghb_obs_supported(this)
Return true because GHB package supports observations.
subroutine ghb_rp(this)
Read and prepare.
subroutine, public ghb_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, mempath)
Create a New Ghb Package and point bndobj to the new package.
subroutine ghb_options(this)
Set options specific to GhbType.
subroutine define_listlabel(this)
Define the list heading that is written to iout when PRINT_INPUT option is used.
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 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.
integer(i4b) function, public count_errors()
Return number of errors.
subroutine, public store_error_filename(filename, terminate)
Store the erroring file name.
subroutine, public store_error_unit(iunit, terminate)
Store the file unit number.
This module contains simulation variables.
character(len=maxcharlen) errmsg
error message string
integer(i4b), pointer, public kper
current stress period number
This class is used to store a single deferred-length character string. It was designed to work in an ...