MODFLOW 6  version 6.7.0.dev2
USGS Modular Hydrologic Model
gwtsrcmodule Module Reference

Data Types

type  gwtsrctype
 

Functions/Subroutines

subroutine, public src_create (packobj, id, ibcnum, inunit, iout, namemodel, depvartype, pakname, fmi)
 Create a source loading package. More...
 
subroutine src_options (this, option, found)
 Set additional options specific to the GwtSrcType. More...
 
subroutine src_da (this)
 Deallocate memory. More...
 
subroutine src_allocate_scalars (this)
 Allocate scalars. More...
 
subroutine src_allocate_arrays (this, nodelist, auxvar)
 Allocate arrays. More...
 
subroutine src_rp (this)
 
subroutine set_nodesontop (this)
 Store nodelist in nodesontop. More...
 
subroutine src_cf (this)
 Formulate the HCOF and RHS terms. More...
 
subroutine src_fc (this, rhs, ia, idxglo, matrix_sln)
 Add matrix terms related to specified mass source loading. More...
 
subroutine define_listlabel (this)
 Define list labels. More...
 
logical function src_obs_supported (this)
 Support function for specified mass source loading observations. More...
 
subroutine src_df_obs (this)
 Define observations. More...
 
subroutine src_rp_ts (this)
 Procedure related to time series. More...
 

Variables

character(len=lenftype) ftype = 'SRC'
 
character(len=16) text = ' SRC'
 

Function/Subroutine Documentation

◆ define_listlabel()

subroutine gwtsrcmodule::define_listlabel ( class(gwtsrctype), intent(inout)  this)
private

Define the list heading that is written to iout when PRINT_INPUT option is used.

Definition at line 309 of file gwt-src.f90.

310  ! -- dummy
311  class(GwtSrcType), intent(inout) :: this
312  ! -- local
313  !
314  ! -- create the header list label
315  this%listlabel = trim(this%filtyp)//' NO.'
316  if (this%dis%ndim == 3) then
317  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'LAYER'
318  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'ROW'
319  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'COL'
320  elseif (this%dis%ndim == 2) then
321  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'LAYER'
322  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'CELL2D'
323  else
324  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'NODE'
325  end if
326  write (this%listlabel, '(a, a16)') trim(this%listlabel), 'STRESS RATE'
327  if (this%inamedbound == 1) then
328  write (this%listlabel, '(a, a16)') trim(this%listlabel), 'BOUNDARY NAME'
329  end if

◆ set_nodesontop()

subroutine gwtsrcmodule::set_nodesontop ( class(gwtsrctype), intent(inout)  this)
private

Definition at line 209 of file gwt-src.f90.

210  implicit none
211  ! -- dummy
212  class(GwtSrcType), intent(inout) :: this
213  ! -- local
214  integer(I4B) :: n
215  ! !
216  ! ! -- allocate if necessary
217  ! if (.not. associated(this%nodesontop)) then
218  ! allocate (this%nodesontop(this%maxbound))
219  ! end if
220  !
221  ! -- copy nodelist into nodesontop
222  do n = 1, this%nbound
223  this%nodesontop(n) = this%nodelist(n)
224  end do

◆ src_allocate_arrays()

subroutine gwtsrcmodule::src_allocate_arrays ( class(gwtsrctype this,
integer(i4b), dimension(:), optional, pointer, contiguous  nodelist,
real(dp), dimension(:, :), optional, pointer, contiguous  auxvar 
)

Allocate scalars specific to this source loading package

Parameters
nodelistpackage nodelist
auxvarpackage aux variable array

Definition at line 166 of file gwt-src.f90.

168  ! -- dummy
169  class(GwtSrcType) :: this
170  integer(I4B), dimension(:), pointer, contiguous, optional :: nodelist !< package nodelist
171  real(DP), dimension(:, :), pointer, contiguous, optional :: auxvar !< package aux variable array
172  ! local
173  integer(I4B) :: n
174  !
175  ! -- call standard BndType allocate scalars
176  call this%BndType%allocate_arrays(nodelist, auxvar)
177  !
178  ! -- allocate the object and assign values to object variables
179  if (this%highest_sat) then
180  call mem_allocate(this%nodesontop, this%maxbound, 'NODESONTOP', &
181  this%memoryPath)
182  end if
183  !
184  ! -- Set values
185  if (this%highest_sat) then
186  do n = 1, this%maxbound
187  this%nodesontop(n) = 0
188  end do
189  end if
190 

◆ src_allocate_scalars()

subroutine gwtsrcmodule::src_allocate_scalars ( class(gwtsrctype this)

Allocate scalars specific to this source loading package

Definition at line 146 of file gwt-src.f90.

148  ! -- dummy
149  class(GwtSrcType) :: this
150  !
151  ! -- call standard BndType allocate scalars
152  call this%BndType%allocate_scalars()
153  !
154  ! -- allocate the object and assign values to object variables
155  call mem_allocate(this%highest_sat, 'HIGHEST_SAT', this%memoryPath)
156  !
157  ! -- Set values
158  this%highest_sat = .false.
159 

◆ src_cf()

subroutine gwtsrcmodule::src_cf ( class(gwtsrctype this)
private

This subroutine:

  • calculates hcof and rhs terms
  • skip if no sources

Definition at line 233 of file gwt-src.f90.

234  ! -- dummy
235  class(GwtSrcType) :: this
236  ! -- local
237  integer(I4B) :: i, node
238  real(DP) :: q
239  !
240  ! -- Return if no sources
241  if (this%nbound == 0) return
242  !
243  ! -- Calculate hcof and rhs for each source entry
244  do i = 1, this%nbound
245  !
246  ! -- Find the node number
247  if (this%highest_sat) then
248  node = this%nodesontop(i)
249  else
250  node = this%nodelist(i)
251  end if
252  !
253  ! -- reset nodelist to highest active
254  if (this%highest_sat) then
255  if (this%fmi%gwfsat(node) == 0) &
256  call this%dis%highest_saturated(node, this%fmi%gwfsat)
257  this%nodelist(i) = node
258  end if
259 
260  this%hcof(i) = dzero
261  if (this%ibound(node) <= 0) then
262  this%rhs(i) = dzero
263  cycle
264  end if
265  q = this%bound(1, i)
266  this%rhs(i) = -q
267  end do

◆ src_create()

subroutine, public gwtsrcmodule::src_create ( class(bndtype), pointer  packobj,
integer(i4b), intent(in)  id,
integer(i4b), intent(in)  ibcnum,
integer(i4b), intent(in)  inunit,
integer(i4b), intent(in)  iout,
character(len=*), intent(in)  namemodel,
character(len=lenvarname), intent(in)  depvartype,
character(len=*), intent(in)  pakname,
type(tspfmitype), intent(in), target  fmi 
)

This subroutine points bndobj to the newly created package

Definition at line 53 of file gwt-src.f90.

55  ! -- dummy
56  class(BndType), pointer :: packobj
57  integer(I4B), intent(in) :: id
58  integer(I4B), intent(in) :: ibcnum
59  integer(I4B), intent(in) :: inunit
60  integer(I4B), intent(in) :: iout
61  character(len=*), intent(in) :: namemodel
62  character(len=*), intent(in) :: pakname
63  character(len=LENVARNAME), intent(in) :: depvartype
64  type(TspFmiType), intent(in), target :: fmi
65  ! -- local
66  type(GwtSrcType), pointer :: srcobj
67  !
68  ! -- allocate the object and assign values to object variables
69  allocate (srcobj)
70  packobj => srcobj
71  !
72  ! -- create name and memory path
73  call packobj%set_names(ibcnum, namemodel, pakname, ftype)
74  packobj%text = text
75  !
76  ! -- allocate scalars
77  call srcobj%allocate_scalars()
78  !
79  ! -- initialize package
80  call packobj%pack_initialize()
81 
82  packobj%inunit = inunit
83  packobj%iout = iout
84  packobj%id = id
85  packobj%ibcnum = ibcnum
86  packobj%ncolbnd = 1
87  packobj%iscloc = 1
88  !
89  ! -- Store the appropriate label based on the dependent variable
90  srcobj%depvartype = depvartype
91 
92  srcobj%fmi => fmi
93 
Here is the caller graph for this function:

◆ src_da()

subroutine gwtsrcmodule::src_da ( class(gwtsrctype this)
private

Definition at line 124 of file gwt-src.f90.

125  ! -- modules
127  ! -- dummy
128  class(GwtSrcType) :: this
129  !
130  ! -- Deallocate parent package
131  call this%BndType%bnd_da()
132  !
133  ! -- arrays
134  if (this%highest_sat) then
135  call mem_deallocate(this%nodesontop, "NODESONTOP", this%memoryPath)
136  end if
137  !
138  ! -- scalars
139  call mem_deallocate(this%highest_sat)

◆ src_df_obs()

subroutine gwtsrcmodule::src_df_obs ( class(gwtsrctype this)
private

This subroutine:

  • stores observation types supported by SRC package.
  • overrides BndTypebnd_df_obs

Definition at line 353 of file gwt-src.f90.

354  implicit none
355  ! -- dummy
356  class(GwtSrcType) :: this
357  ! -- local
358  integer(I4B) :: indx
359  !
360  call this%obs%StoreObsType('src', .true., indx)
361  this%obs%obsData(indx)%ProcessIdPtr => defaultobsidprocessor
362  !
363  ! -- Store obs type and assign procedure pointer
364  ! for to-mvr observation type.
365  call this%obs%StoreObsType('to-mvr', .true., indx)
366  this%obs%obsData(indx)%ProcessIdPtr => defaultobsidprocessor
Here is the call graph for this function:

◆ src_fc()

subroutine gwtsrcmodule::src_fc ( class(gwtsrctype this,
real(dp), dimension(:), intent(inout)  rhs,
integer(i4b), dimension(:), intent(in)  ia,
integer(i4b), dimension(:), intent(in)  idxglo,
class(matrixbasetype), pointer  matrix_sln 
)
private

Copy rhs and hcof into solution rhs and amat

Definition at line 274 of file gwt-src.f90.

275  ! -- dummy
276  class(GwtSrcType) :: this
277  real(DP), dimension(:), intent(inout) :: rhs
278  integer(I4B), dimension(:), intent(in) :: ia
279  integer(I4B), dimension(:), intent(in) :: idxglo
280  class(MatrixBaseType), pointer :: matrix_sln
281  ! -- local
282  integer(I4B) :: i, n, ipos
283  !
284  ! -- pakmvrobj fc
285  if (this%imover == 1) then
286  call this%pakmvrobj%fc()
287  end if
288  !
289  ! -- Copy package rhs and hcof into solution rhs and amat
290  do i = 1, this%nbound
291  n = this%nodelist(i)
292  rhs(n) = rhs(n) + this%rhs(i)
293  ipos = ia(n)
294  call matrix_sln%add_value_pos(idxglo(ipos), this%hcof(i))
295  !
296  ! -- If mover is active and mass is being withdrawn,
297  ! store available mass (as positive value).
298  if (this%imover == 1 .and. this%rhs(i) > dzero) then
299  call this%pakmvrobj%accumulate_qformvr(i, this%rhs(i))
300  end if
301  end do

◆ src_obs_supported()

logical function gwtsrcmodule::src_obs_supported ( class(gwtsrctype this)
private

This function:

  • returns true because SRC package supports observations.
  • overrides BndTypebnd_obs_supported()

Definition at line 339 of file gwt-src.f90.

340  implicit none
341  ! -- dummy
342  class(GwtSrcType) :: this
343  !
344  src_obs_supported = .true.

◆ src_options()

subroutine gwtsrcmodule::src_options ( class(gwtsrctype), intent(inout)  this,
character(len=*), intent(inout)  option,
logical, intent(inout)  found 
)
private

This routine overrides BndTypebnd_options

Definition at line 100 of file gwt-src.f90.

101  ! -- dummy
102  class(GwtSrcType), intent(inout) :: this
103  character(len=*), intent(inout) :: option
104  logical, intent(inout) :: found
105  ! -- local
106  ! -- formats
107 
108  found = .true.
109  select case (option)
110  case ('HIGHEST_SATURATED')
111  this%highest_sat = .true.
112  write (this%iout, '(4x,a)') &
113  'Mass source loading rate will be applied to the highest cell at or below &
114  &the specified cellid with a non-zero saturation.'
115  case default
116  !
117  ! -- No options found
118  found = .false.
119  end select

◆ src_rp()

subroutine gwtsrcmodule::src_rp ( class(gwtsrctype), intent(inout)  this)
Parameters
[in,out]thisBndType object

Definition at line 193 of file gwt-src.f90.

194  ! -- modules
195  ! -- dummy
196  class(GwtSrcType), intent(inout) :: this !< BndType object
197 
198  ! call standard BndType rp
199  call this%BndType%bnd_rp()
200 
201  if (this%highest_sat) call this%set_nodesontop()
202 
203  return
204 

◆ src_rp_ts()

subroutine gwtsrcmodule::src_rp_ts ( class(gwtsrctype), intent(inout)  this)
private

Assign tsLinkText appropriately for all time series in use by package. In the SRC package only the SMASSRATE variable can be controlled by time series.

Definition at line 375 of file gwt-src.f90.

376  ! -- dummy
377  class(GwtSrcType), intent(inout) :: this
378  ! -- local
379  integer(I4B) :: i, nlinks
380  type(TimeSeriesLinkType), pointer :: tslink => null()
381  !
382  nlinks = this%TsManager%boundtslinks%Count()
383  do i = 1, nlinks
384  tslink => gettimeserieslinkfromlist(this%TsManager%boundtslinks, i)
385  if (associated(tslink)) then
386  if (tslink%JCol == 1) then
387  tslink%Text = 'SMASSRATE'
388  end if
389  end if
390  end do
Here is the call graph for this function:

Variable Documentation

◆ ftype

character(len=lenftype) gwtsrcmodule::ftype = 'SRC'
private

Definition at line 18 of file gwt-src.f90.

18  character(len=LENFTYPE) :: ftype = 'SRC'

◆ text

character(len=16) gwtsrcmodule::text = ' SRC'
private

Definition at line 19 of file gwt-src.f90.

19  character(len=16) :: text = ' SRC'