MODFLOW 6  version 6.7.0.dev3
USGS Modular Hydrologic Model
methoddismodule Module Reference

Data Types

type  methoddistype
 

Functions/Subroutines

subroutine, public create_method_dis (method)
 Create a new structured grid (DIS) tracking method. More...
 
subroutine deallocate (this)
 Destructor the tracking method. More...
 
subroutine load_cell (this, ic, cell)
 
subroutine load_dis (this, particle, next_level, submethod)
 Load the cell geometry and tracking method. More...
 
subroutine load_particle (this, cell, particle)
 Load cell properties into the particle, including. More...
 
subroutine update_flowja (this, cell, particle)
 Update cell-cell flows of particle mass. Every particle is currently assigned unit mass. More...
 
subroutine pass_dis (this, particle)
 Pass a particle to the next cell, if there is one. More...
 
subroutine apply_dis (this, particle, tmax)
 Apply the structured tracking method to a particle. More...
 
double precision function get_top (this, iatop)
 Returns a top elevation based on index iatop. More...
 
subroutine load_cell_defn (this, ic, defn)
 Loads cell definition from the grid. More...
 
subroutine load_cell_properties (this, ic, defn)
 
subroutine load_cell_neighbors (this, defn)
 Loads face neighbors to cell definition from the grid. Assumes cell index and number of vertices are already loaded. More...
 
subroutine load_cell_flows (this, defn)
 Load flows into the cell definition. These include face, boundary and net distributed flows. Assumes cell index and number of vertices are already loaded. More...
 
subroutine load_cell_face_flows (this, defn)
 
subroutine load_cell_boundary_flows (this, defn)
 Add boundary flows to the cell definition faceflow array. Assumes cell index and number of vertices are already loaded. More...
 

Function/Subroutine Documentation

◆ apply_dis()

subroutine methoddismodule::apply_dis ( class(methoddistype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle,
real(dp), intent(in)  tmax 
)

Definition at line 299 of file MethodDis.f90.

300  class(MethodDisType), intent(inout) :: this
301  type(ParticleType), pointer, intent(inout) :: particle
302  real(DP), intent(in) :: tmax
303 
304  call this%track(particle, 1, tmax)

◆ create_method_dis()

subroutine, public methoddismodule::create_method_dis ( type(methoddistype), pointer  method)

Definition at line 44 of file MethodDis.f90.

45  ! dummy
46  type(MethodDisType), pointer :: method
47  ! local
48  type(CellRectType), pointer :: cell
49 
50  allocate (method)
51  allocate (method%name)
52  call create_cell_rect(cell)
53  method%cell => cell
54  method%name = "dis"
55  method%delegates = .true.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ deallocate()

subroutine methoddismodule::deallocate ( class(methoddistype), intent(inout)  this)
private

Definition at line 59 of file MethodDis.f90.

60  class(MethodDisType), intent(inout) :: this
61  deallocate (this%name)

◆ get_top()

double precision function methoddismodule::get_top ( class(methoddistype), intent(inout)  this,
integer, intent(in)  iatop 
)
private

Definition at line 308 of file MethodDis.f90.

309  class(MethodDisType), intent(inout) :: this
310  integer, intent(in) :: iatop
311  double precision :: top
312 
313  if (iatop .lt. 0) then
314  top = this%fmi%dis%top(-iatop)
315  else
316  top = this%fmi%dis%bot(iatop)
317  end if

◆ load_cell()

subroutine methoddismodule::load_cell ( class(methoddistype), intent(inout)  this,
integer(i4b), intent(in)  ic,
type(cellrecttype), intent(inout)  cell 
)
private

Definition at line 64 of file MethodDis.f90.

65  ! dummy
66  class(MethodDisType), intent(inout) :: this
67  integer(I4B), intent(in) :: ic
68  type(CellRectType), intent(inout) :: cell
69  ! local
70  integer(I4B) :: icu
71  integer(I4B) :: irow
72  integer(I4B) :: jcol
73  integer(I4B) :: klay
74  real(DP) :: areax
75  real(DP) :: areay
76  real(DP) :: areaz
77  real(DP) :: dx
78  real(DP) :: dy
79  real(DP) :: dz
80  real(DP) :: factor
81  real(DP) :: term
82 
83  select type (dis => this%fmi%dis)
84  type is (distype)
85  icu = dis%get_nodeuser(ic)
86 
87  call get_ijk(icu, dis%nrow, dis%ncol, dis%nlay, &
88  irow, jcol, klay)
89 
90  dx = dis%delr(jcol)
91  dy = dis%delc(irow)
92  dz = cell%defn%top - cell%defn%bot
93 
94  cell%dx = dx
95  cell%dy = dy
96  cell%dz = dz
97  cell%sinrot = dzero
98  cell%cosrot = done
99  cell%xOrigin = cell%defn%polyvert(1, 1)
100  cell%yOrigin = cell%defn%polyvert(2, 1)
101  cell%zOrigin = cell%defn%bot
102  cell%ipvOrigin = 1
103 
104  factor = done / cell%defn%retfactor
105  factor = factor / cell%defn%porosity
106 
107  areaz = dx * dy
108  term = factor / areaz
109 
110  cell%vz1 = cell%defn%faceflow(6) * term
111  cell%vz2 = -cell%defn%faceflow(7) * term
112 
113  if (this%fmi%ibdgwfsat0(ic) == 0) then
114  cell%vx1 = dzero
115  cell%vx2 = dzero
116  cell%vy1 = dzero
117  cell%vy2 = dzero
118  cell%vz1 = dzero
119  cell%vz2 = dzero
120  return
121  end if
122 
123  areax = dy * dz
124  term = factor / areax
125  cell%vx1 = cell%defn%faceflow(1) * term
126  cell%vx2 = -cell%defn%faceflow(3) * term
127 
128  areay = dx * dz
129  term = factor / areay
130  cell%vy1 = cell%defn%faceflow(4) * term
131  cell%vy2 = -cell%defn%faceflow(2) * term
132 
133  end select
Here is the call graph for this function:

◆ load_cell_boundary_flows()

subroutine methoddismodule::load_cell_boundary_flows ( class(methoddistype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 488 of file MethodDis.f90.

489  ! dummy
490  class(MethodDisType), intent(inout) :: this
491  type(CellDefnType), pointer, intent(inout) :: defn
492  ! local
493  integer(I4B) :: ioffset
494 
495  ioffset = (defn%icell - 1) * this%fmi%max_faces
496  defn%faceflow(1) = defn%faceflow(1) + &
497  this%fmi%BoundaryFlows(ioffset + 1)
498  defn%faceflow(2) = defn%faceflow(2) + &
499  this%fmi%BoundaryFlows(ioffset + 2)
500  defn%faceflow(3) = defn%faceflow(3) + &
501  this%fmi%BoundaryFlows(ioffset + 3)
502  defn%faceflow(4) = defn%faceflow(4) + &
503  this%fmi%BoundaryFlows(ioffset + 4)
504  defn%faceflow(5) = defn%faceflow(1)
505  defn%faceflow(6) = defn%faceflow(6) + &
506  this%fmi%BoundaryFlows(ioffset + this%fmi%max_faces - 1)
507  defn%faceflow(7) = defn%faceflow(7) + &
508  this%fmi%BoundaryFlows(ioffset + this%fmi%max_faces)

◆ load_cell_defn()

subroutine methoddismodule::load_cell_defn ( class(methoddistype), intent(inout)  this,
integer(i4b), intent(in)  ic,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 321 of file MethodDis.f90.

322  ! dummy
323  class(MethodDisType), intent(inout) :: this
324  integer(I4B), intent(in) :: ic
325  type(CellDefnType), pointer, intent(inout) :: defn
326 
327  call this%load_cell_properties(ic, defn)
328  call this%fmi%dis%get_polyverts( &
329  defn%icell, &
330  defn%polyvert, &
331  closed=.true.)
332  call this%load_cell_neighbors(defn)
333  call this%load_cell_saturation_status(defn)
334  defn%ispv180(1:defn%npolyverts + 1) = .false.
335  call this%load_cell_flows(defn)
336 

◆ load_cell_face_flows()

subroutine methoddismodule::load_cell_face_flows ( class(methoddistype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 468 of file MethodDis.f90.

469  ! dummy
470  class(MethodDisType), intent(inout) :: this
471  type(CellDefnType), pointer, intent(inout) :: defn
472  ! local
473  integer(I4B) :: m, n, nfaces
474  real(DP) :: q
475 
476  nfaces = defn%npolyverts + 3
477  do m = 1, nfaces
478  n = defn%facenbr(m)
479  if (n > 0) then
480  q = this%fmi%gwfflowja(this%fmi%dis%con%ia(defn%icell) + n)
481  defn%faceflow(m) = defn%faceflow(m) + q
482  end if
483  end do

◆ load_cell_flows()

subroutine methoddismodule::load_cell_flows ( class(methoddistype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 441 of file MethodDis.f90.

442  class(MethodDisType), intent(inout) :: this
443  type(CellDefnType), pointer, intent(inout) :: defn
444 
445  ! Load face flows, including boundary flows. As with cell verts,
446  ! the face flow array wraps around. Top and bottom flows make up
447  ! the last two elements, respectively, for size npolyverts + 3.
448  ! If there is no flow through any face, set a no-exit-face flag.
449  defn%faceflow = dzero
450  call this%load_cell_boundary_flows(defn)
451  call this%load_cell_face_flows(defn)
452  call this%cap_cell_wt_flow(defn)
453  call this%load_cell_no_exit_face(defn)
454 
455  ! Add up net distributed flow
456  defn%distflow = this%fmi%SourceFlows(defn%icell) + &
457  this%fmi%SinkFlows(defn%icell) + &
458  this%fmi%StorageFlows(defn%icell)
459 
460  ! Set weak sink flag
461  if (this%fmi%SinkFlows(defn%icell) .ne. dzero) then
462  defn%iweaksink = 1
463  else
464  defn%iweaksink = 0
465  end if

◆ load_cell_neighbors()

subroutine methoddismodule::load_cell_neighbors ( class(methoddistype), intent(inout)  this,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 372 of file MethodDis.f90.

373  ! dummy
374  class(MethodDisType), intent(inout) :: this
375  type(CellDefnType), pointer, intent(inout) :: defn
376  ! local
377  integer(I4B) :: ic1
378  integer(I4B) :: ic2
379  integer(I4B) :: icu1
380  integer(I4B) :: icu2
381  integer(I4B) :: j1
382  integer(I4B) :: iloc
383  integer(I4B) :: ipos
384  integer(I4B) :: irow1
385  integer(I4B) :: irow2
386  integer(I4B) :: jcol1
387  integer(I4B) :: jcol2
388  integer(I4B) :: klay1
389  integer(I4B) :: klay2
390  integer(I4B) :: iedgeface
391 
392  select type (dis => this%fmi%dis)
393  type is (distype)
394  ! Load face neighbors
395  defn%facenbr = 0
396  ic1 = defn%icell
397  icu1 = dis%get_nodeuser(ic1)
398  call get_ijk(icu1, dis%nrow, dis%ncol, dis%nlay, &
399  irow1, jcol1, klay1)
400  call get_jk(icu1, dis%get_ncpl(), dis%nlay, j1, klay1)
401  do iloc = 1, dis%con%ia(ic1 + 1) - dis%con%ia(ic1) - 1
402  ipos = dis%con%ia(ic1) + iloc
403  ! mask could become relevant if PRT uses interface model
404  if (dis%con%mask(ipos) == 0) cycle
405  ic2 = dis%con%ja(ipos)
406  icu2 = dis%get_nodeuser(ic2)
407  call get_ijk(icu2, dis%nrow, dis%ncol, dis%nlay, &
408  irow2, jcol2, klay2)
409  if (klay2 == klay1) then
410  ! Edge (polygon) face neighbor
411  if (irow2 > irow1) then
412  ! Neighbor to the S
413  iedgeface = 4
414  else if (jcol2 > jcol1) then
415  ! Neighbor to the E
416  iedgeface = 3
417  else if (irow2 < irow1) then
418  ! Neighbor to the N
419  iedgeface = 2
420  else
421  ! Neighbor to the W
422  iedgeface = 1
423  end if
424  defn%facenbr(iedgeface) = int(iloc, 1)
425  else if (klay2 > klay1) then
426  ! Bottom face neighbor
427  defn%facenbr(defn%npolyverts + 2) = int(iloc, 1)
428  else
429  ! Top face neighbor
430  defn%facenbr(defn%npolyverts + 3) = int(iloc, 1)
431  end if
432  end do
433  end select
434  ! List of edge (polygon) faces wraps around
435  defn%facenbr(defn%npolyverts + 1) = defn%facenbr(1)
Here is the call graph for this function:

◆ load_cell_properties()

subroutine methoddismodule::load_cell_properties ( class(methoddistype), intent(inout)  this,
integer(i4b), intent(in)  ic,
type(celldefntype), intent(inout), pointer  defn 
)
private

Definition at line 339 of file MethodDis.f90.

340  ! dummy
341  class(MethodDisType), intent(inout) :: this
342  integer(I4B), intent(in) :: ic
343  type(CellDefnType), pointer, intent(inout) :: defn
344  ! local
345  integer(I4B) :: irow, icol, ilay, icu
346 
347  defn%icell = ic
348  defn%npolyverts = 4 ! rectangular cell always has 4 vertices
349  defn%iatop = get_iatop(this%fmi%dis%get_ncpl(), &
350  this%fmi%dis%get_nodeuser(ic))
351  defn%top = this%fmi%dis%bot(ic) + &
352  this%fmi%gwfsat(ic) * &
353  (this%fmi%dis%top(ic) - this%fmi%dis%bot(ic))
354  defn%bot = this%fmi%dis%bot(ic)
355  defn%sat = this%fmi%gwfsat(ic)
356  defn%porosity = this%porosity(ic)
357  defn%retfactor = this%retfactor(ic)
358  select type (dis => this%fmi%dis)
359  type is (distype)
360  icu = dis%get_nodeuser(ic)
361  call get_ijk(icu, dis%nrow, dis%ncol, dis%nlay, irow, icol, ilay)
362  defn%ilay = ilay
363  end select
364  defn%izone = this%izone(ic)
365  defn%can_be_rect = .true.
366  defn%can_be_quad = .false.
367 
Here is the call graph for this function:

◆ load_dis()

subroutine methoddismodule::load_dis ( class(methoddistype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle,
integer(i4b), intent(in)  next_level,
class(methodtype), intent(inout), pointer  submethod 
)
private

Definition at line 137 of file MethodDis.f90.

138  ! dummy
139  class(MethodDisType), intent(inout) :: this
140  type(ParticleType), pointer, intent(inout) :: particle
141  integer(I4B), intent(in) :: next_level
142  class(MethodType), pointer, intent(inout) :: submethod
143  ! local
144  integer(I4B) :: ic
145 
146  select type (cell => this%cell)
147  type is (cellrecttype)
148  ic = particle%itrdomain(next_level)
149  call this%load_cell_defn(ic, cell%defn)
150  call this%load_cell(ic, cell)
151  if (this%fmi%ibdgwfsat0(ic) == 0) then
152  call method_cell_ptb%init( &
153  fmi=this%fmi, &
154  cell=this%cell, &
155  events=this%events, &
156  tracktimes=this%tracktimes)
157  submethod => method_cell_ptb
158  else
159  call method_cell_plck%init( &
160  fmi=this%fmi, &
161  cell=this%cell, &
162  events=this%events, &
163  tracktimes=this%tracktimes)
164  submethod => method_cell_plck
165  end if
166  end select

◆ load_particle()

subroutine methoddismodule::load_particle ( class(methoddistype), intent(inout)  this,
type(cellrecttype), intent(inout), pointer  cell,
type(particletype), intent(inout), pointer  particle 
)
private

Definition at line 171 of file MethodDis.f90.

172  use particlemodule, only: term_boundary
173  use particleeventmodule, only: terminate
174  ! dummy
175  class(MethodDisType), intent(inout) :: this
176  type(CellRectType), pointer, intent(inout) :: cell
177  type(ParticleType), pointer, intent(inout) :: particle
178  ! local
179  integer(I4B) :: ic
180  integer(I4B) :: icu
181  integer(I4B) :: ilay
182  integer(I4B) :: irow
183  integer(I4B) :: icol
184  integer(I4B) :: inface
185  integer(I4B) :: idiag
186  integer(I4B) :: inbr
187  integer(I4B) :: ipos
188  real(DP) :: z
189  real(DP) :: zrel
190  real(DP) :: topfrom
191  real(DP) :: botfrom
192  real(DP) :: top
193  real(DP) :: bot
194  real(DP) :: sat
195 
196  select type (dis => this%fmi%dis)
197  type is (distype)
198  ! compute reduced/user node numbers and layer
199  inface = particle%iboundary(level_feature)
200  inbr = cell%defn%facenbr(inface)
201  idiag = this%fmi%dis%con%ia(cell%defn%icell)
202  ipos = idiag + inbr
203  ic = dis%con%ja(ipos)
204  icu = dis%get_nodeuser(ic)
205  call get_ijk(icu, dis%nrow, dis%ncol, dis%nlay, irow, icol, ilay)
206 
207  ! update node numbers and layer
208  particle%itrdomain(level_feature) = ic
209  particle%icu = icu
210  particle%ilay = ilay
211 
212  ! map/set particle entry face
213  if (inface .eq. 1) then
214  inface = 3
215  else if (inface .eq. 2) then
216  inface = 4
217  else if (inface .eq. 3) then
218  inface = 1
219  else if (inface .eq. 4) then
220  inface = 2
221  else if (inface .eq. 6) then
222  inface = 7
223  else if (inface .eq. 7) then
224  inface = 6
225  end if
226  particle%iboundary(level_feature) = inface
227 
228  ! map z between cells
229  z = particle%z
230  if (inface < 5) then
231  topfrom = cell%defn%top
232  botfrom = cell%defn%bot
233  zrel = (z - botfrom) / (topfrom - botfrom)
234  top = dis%top(ic)
235  bot = dis%bot(ic)
236  sat = this%fmi%gwfsat(ic)
237  z = bot + zrel * sat * (top - bot)
238  end if
239  particle%z = z
240  end select
@, public terminate
particle terminated
@ term_boundary
terminated at a boundary face
Definition: Particle.f90:30
Here is the call graph for this function:

◆ pass_dis()

subroutine methoddismodule::pass_dis ( class(methoddistype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle 
)
private

Definition at line 270 of file MethodDis.f90.

271  use particlemodule, only: term_boundary
272  use particleeventmodule, only: terminate
273  ! dummy
274  class(MethodDisType), intent(inout) :: this
275  type(ParticleType), pointer, intent(inout) :: particle
276  ! local
277  type(CellRectType), pointer :: cell
278  integer(I4B) :: iface
279  logical(LGP) :: no_neighbors
280 
281  select type (c => this%cell)
282  type is (cellrecttype)
283  cell => c
284  iface = particle%iboundary(level_feature)
285  no_neighbors = cell%defn%facenbr(iface) == 0
286 
287  ! todo AMP: reconsider when multiple models supported
288  if (no_neighbors) then
289  call this%terminate(particle, status=term_boundary)
290  return
291  end if
292 
293  call this%load_particle(cell, particle)
294  call this%update_flowja(cell, particle)
295  end select

◆ update_flowja()

subroutine methoddismodule::update_flowja ( class(methoddistype), intent(inout)  this,
type(cellrecttype), intent(inout), pointer  cell,
type(particletype), intent(inout), pointer  particle 
)

Definition at line 245 of file MethodDis.f90.

246  ! dummy
247  class(MethodDisType), intent(inout) :: this
248  type(CellRectType), pointer, intent(inout) :: cell
249  type(ParticleType), pointer, intent(inout) :: particle
250  ! local
251  integer(I4B) :: idiag
252  integer(I4B) :: inbr
253  integer(I4B) :: inface
254  integer(I4B) :: ipos
255 
256  inface = particle%iboundary(level_feature)
257  inbr = cell%defn%facenbr(inface)
258  idiag = this%fmi%dis%con%ia(cell%defn%icell)
259  ipos = idiag + inbr
260 
261  ! leaving old cell
262  this%flowja(ipos) = this%flowja(ipos) - done
263 
264  ! entering new cell
265  this%flowja(this%fmi%dis%con%isym(ipos)) &
266  = this%flowja(this%fmi%dis%con%isym(ipos)) + done