MODFLOW 6  version 6.7.0.dev1
USGS Modular Hydrologic Model
methodcellpollockquadmodule Module Reference

Data Types

type  methodcellpollockquadtype
 

Functions/Subroutines

procedure subroutine, public create_method_cell_quad (method)
 Create a new Pollock quad-refined cell method. More...
 
subroutine deallocate (this)
 Deallocate the Pollock quad-refined cell method. More...
 
subroutine load_mcpq (this, particle, next_level, submethod)
 Load subcell into tracking method. More...
 
subroutine pass_mcpq (this, particle)
 Pass particle to next subcell if there is one, or to the cell face. More...
 
subroutine apply_mcpq (this, particle, tmax)
 Solve the quad-rectangular cell via Pollock's method. More...
 
subroutine load_subcell (this, particle, subcell)
 Load the rectangular subcell from the rectangular cell. More...
 

Function/Subroutine Documentation

◆ apply_mcpq()

subroutine methodcellpollockquadmodule::apply_mcpq ( class(methodcellpollockquadtype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle,
real(dp), intent(in)  tmax 
)
private

Definition at line 202 of file MethodCellPollockQuad.f90.

203  ! dummy
204  class(MethodCellPollockQuadType), intent(inout) :: this
205  type(ParticleType), pointer, intent(inout) :: particle
206  real(DP), intent(in) :: tmax
207  ! local
208  double precision :: xOrigin, yOrigin, zOrigin, sinrot, cosrot
209 
210  select type (cell => this%cell)
211  type is (cellrectquadtype)
212  ! Check termination/reporting conditions
213  call this%check(particle, cell%defn, tmax)
214  if (.not. particle%advancing) return
215 
216  ! Transform model coordinates to local cell coordinates
217  ! (translated/rotated but not scaled relative to model)
218  xorigin = cell%xOrigin
219  yorigin = cell%yOrigin
220  zorigin = cell%zOrigin
221  sinrot = cell%sinrot
222  cosrot = cell%cosrot
223  call particle%transform(xorigin, yorigin, zorigin, &
224  sinrot, cosrot)
225 
226  ! Track the particle over the cell
227  call this%track(particle, 2, tmax)
228 
229  ! Transform cell coordinates back to model coordinates
230  call particle%transform(xorigin, yorigin, zorigin, &
231  sinrot, cosrot, invert=.true.)
232  call particle%reset_transform()
233  end select

◆ create_method_cell_quad()

procedure subroutine, public methodcellpollockquadmodule::create_method_cell_quad ( type(methodcellpollockquadtype), pointer  method)

Definition at line 31 of file MethodCellPollockQuad.f90.

32  ! dummy
33  type(MethodCellPollockQuadType), pointer :: method
34  ! local
35  type(CellRectQuadType), pointer :: cell
36  type(SubcellRectType), pointer :: subcell
37 
38  allocate (method)
39  call create_cell_rect_quad(cell)
40  method%cell => cell
41  method%name => method%cell%type
42  method%delegates = .true.
43  call create_subcell_rect(subcell)
44  method%subcell => subcell
Here is the call graph for this function:
Here is the caller graph for this function:

◆ deallocate()

subroutine methodcellpollockquadmodule::deallocate ( class(methodcellpollockquadtype), intent(inout)  this)
private

Definition at line 48 of file MethodCellPollockQuad.f90.

49  class(MethodCellPollockQuadType), intent(inout) :: this
50  deallocate (this%name)

◆ load_mcpq()

subroutine methodcellpollockquadmodule::load_mcpq ( class(methodcellpollockquadtype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle,
integer, intent(in)  next_level,
class(methodtype), intent(inout), pointer  submethod 
)
private

Definition at line 54 of file MethodCellPollockQuad.f90.

55  class(MethodCellPollockQuadType), intent(inout) :: this
56  type(ParticleType), pointer, intent(inout) :: particle
57  integer, intent(in) :: next_level
58  class(MethodType), pointer, intent(inout) :: submethod
59 
60  select type (subcell => this%subcell)
61  type is (subcellrecttype)
62  call this%load_subcell(particle, subcell)
63  end select
64  call method_subcell_plck%init( &
65  fmi=this%fmi, &
66  cell=this%cell, &
67  subcell=this%subcell, &
68  trackctl=this%trackctl, &
69  tracktimes=this%tracktimes)
70  submethod => method_subcell_plck

◆ load_subcell()

subroutine methodcellpollockquadmodule::load_subcell ( class(methodcellpollockquadtype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle,
class(subcellrecttype), intent(inout)  subcell 
)
private

Definition at line 237 of file MethodCellPollockQuad.f90.

238  ! dummy
239  class(MethodCellPollockQuadType), intent(inout) :: this
240  type(ParticleType), pointer, intent(inout) :: particle
241  class(SubcellRectType), intent(inout) :: subcell
242  ! local
243  real(DP) :: dx, dy, dz, areax, areay, areaz
244  real(DP) :: dxprel, dyprel
245  integer(I4B) :: isc, npolyverts, m1, m2
246  real(DP) :: qextl1, qextl2, qintl1, qintl2
247  real(DP) :: factor, term
248 
249  select type (cell => this%cell)
250  type is (cellrectquadtype)
251  factor = done / cell%defn%retfactor
252  factor = factor / cell%defn%porosity
253  npolyverts = cell%defn%npolyverts
254 
255  isc = particle%idomain(3)
256  ! Subcells 1, 2, 3, and 4 are Pollock's subcells A, B, C, and D,
257  ! respectively
258 
259  dx = cell%dx
260  dy = cell%dy
261  ! If not already known, determine subcell number
262  if (isc .le. 0) then
263  dxprel = particle%x / dx
264  dyprel = particle%y / dy
265 
266  if (dyprel .ge. 5d-1) then
267  if (dxprel .le. 5d-1) then
268  isc = 4
269  else
270  isc = 1
271  end if
272  else
273  if (dxprel .le. 5d-1) then
274  isc = 3
275  else
276  isc = 2
277  end if
278  end if
279 
280  subcell%isubcell = isc
281  particle%idomain(3) = isc
282  end if
283  dx = 5d-1 * dx
284  dy = 5d-1 * dy
285  dz = cell%defn%top - &
286  cell%defn%bot
287  areax = dy * dz
288  areay = dx * dz
289  areaz = dx * dy
290  qintl1 = cell%qintl(isc)
291  ! qintl list wraps around, so isc+1=5 is ok
292  qintl2 = cell%qintl(isc + 1)
293  qextl1 = cell%qextl1(isc)
294  qextl2 = cell%qextl2(isc)
295 
296  subcell%dx = dx
297  subcell%dy = dy
298  subcell%dz = dz
299  subcell%sinrot = dzero
300  subcell%cosrot = done
301  subcell%zOrigin = dzero
302  select case (isc)
303  case (1)
304  subcell%xOrigin = dx
305  subcell%yOrigin = dy
306  term = factor / areax
307  subcell%vx1 = qintl1 * term
308  subcell%vx2 = -qextl2 * term
309  term = factor / areay
310  subcell%vy1 = -qintl2 * term
311  subcell%vy2 = -qextl1 * term
312  case (2)
313  subcell%xOrigin = dx
314  subcell%yOrigin = dzero
315  term = factor / areax
316  subcell%vx1 = -qintl2 * term
317  subcell%vx2 = -qextl1 * term
318  term = factor / areay
319  subcell%vy1 = qextl2 * term
320  subcell%vy2 = -qintl1 * term
321  case (3)
322  subcell%xOrigin = dzero
323  subcell%yOrigin = dzero
324  term = factor / areax
325  subcell%vx1 = qextl2 * term
326  subcell%vx2 = -qintl1 * term
327  term = factor / areay
328  subcell%vy1 = qextl1 * term
329  subcell%vy2 = qintl2 * term
330  case (4)
331  subcell%xOrigin = dzero
332  subcell%yOrigin = dy
333  term = factor / areax
334  subcell%vx1 = qextl1 * term
335  subcell%vx2 = qintl2 * term
336  term = factor / areay
337  subcell%vy1 = qintl1 * term
338  subcell%vy2 = -qextl2 * term
339  end select
340  m1 = npolyverts + 2
341  m2 = m1 + 1
342  term = factor / areaz
343  subcell%vz1 = 2.5d-1 * cell%defn%faceflow(m1) * term
344  subcell%vz2 = -2.5d-1 * cell%defn%faceflow(m2) * term
345  end select

◆ pass_mcpq()

subroutine methodcellpollockquadmodule::pass_mcpq ( class(methodcellpollockquadtype), intent(inout)  this,
type(particletype), intent(inout), pointer  particle 
)
private

Definition at line 74 of file MethodCellPollockQuad.f90.

75  ! dummy
76  class(MethodCellPollockQuadType), intent(inout) :: this
77  type(ParticleType), pointer, intent(inout) :: particle
78  ! local
79  integer(I4B) :: isc, exitFace, npolyverts, inface, infaceoff
80 
81  select type (cell => this%cell)
82  type is (cellrectquadtype)
83  exitface = particle%iboundary(3)
84  isc = particle%idomain(3)
85  npolyverts = cell%defn%npolyverts
86 
87  ! exitFace uses MODPATH 7 iface convention here
88  select case (exitface)
89  case (0)
90  ! Subcell interior (cell interior)
91  inface = -1
92  case (1)
93  select case (isc)
94  case (1)
95  ! W face, subcell 1 --> E face, subcell 4 (cell interior)
96  particle%idomain(3) = 4
97  particle%iboundary(3) = 2
98  inface = 0 ! want Domain(2) unchanged; Boundary(2) = 0
99  case (2)
100  ! W face, subcell 2 --> E face, subcell 3 (cell interior)
101  particle%idomain(3) = 3
102  particle%iboundary(3) = 2
103  inface = 0 ! want Domain(2) unchanged; Boundary(2) = 0
104  case (3)
105  ! W face, subcell 3 (cell face)
106  inface = 1 ! want Domain(2) = -Domain(2); Boundary(2) = inface
107  infaceoff = 0
108  case (4)
109  ! W face, subcell 4 (cell face)
110  inface = 2 ! want Domain(2) = -Domain(2); Boundary(2) = inface
111  infaceoff = -1
112  end select
113  case (2)
114  select case (isc)
115  case (1)
116  ! E face, subcell 1 (cell face)
117  inface = 3 ! want Domain(2) = -Domain(2); Boundary(2) = inface
118  infaceoff = 0
119  case (2)
120  ! E face, subcell 2 (cell face)
121  inface = 4 ! want Domain(2) = -Domain(2); Boundary(2) = inface
122  infaceoff = -1
123  case (3)
124  ! E face, subcell 3 --> W face, subcell 2 (cell interior)
125  particle%idomain(3) = 2
126  particle%iboundary(3) = 1
127  inface = 0 ! want Domain(2) unchanged; Boundary(2) = 0
128  case (4)
129  ! E face, subcell 4 --> W face subcell 1 (cell interior)
130  particle%idomain(3) = 1
131  particle%iboundary(3) = 1
132  inface = 0 ! want Domain(2) unchanged; Boundary(2) = 0
133  end select
134  case (3)
135  select case (isc)
136  case (1)
137  ! S face, subcell 1 --> N face, subcell 2 (cell interior)
138  particle%idomain(3) = 2
139  particle%iboundary(3) = 4
140  inface = 0 ! want Domain(2) unchanged; Boundary(2) = 0
141  case (2)
142  ! S face, subcell 2 (cell face)
143  inface = 4 ! want Domain(2) = -Domain(2); Boundary(2) = inface
144  infaceoff = 0
145  case (3)
146  ! S face, subcell 3 (cell face)
147  inface = 1 ! want Domain(2) = -Domain(2); Boundary(2) = inface
148  infaceoff = -1
149  case (4)
150  ! S face, subcell 4 --> N face, subcell 3 (cell interior)
151  particle%idomain(3) = 3
152  particle%iboundary(3) = 4
153  inface = 0 ! want Domain(2) unchanged; Boundary(2) = 0
154  end select
155  case (4)
156  select case (isc)
157  case (1)
158  ! N face, subcell 1 (cell face)
159  inface = 3 ! want Domain(2) = -Domain(2); Boundary(2) = inface
160  infaceoff = -1
161  case (2)
162  ! N face, subcell 2 --> S face, subcell 1 (cell interior)
163  particle%idomain(3) = 1
164  particle%iboundary(3) = 3
165  inface = 0 ! want Domain(2) unchanged; Boundary(2) = 0
166  case (3)
167  ! N face, subcell 3 --> S face, subcell 4 (cell interior)
168  particle%idomain(3) = 4
169  particle%iboundary(3) = 3
170  inface = 0 ! want Domain(2) unchanged; Boundary(2) = 0
171  case (4)
172  ! N face, subcell 4 (cell face)
173  inface = 2 ! want Domain(2) = -Domain(2); Boundary(2) = inface
174  infaceoff = 0
175  end select
176  case (5)
177  ! Subcell bottom (cell bottom)
178  inface = npolyverts + 2 ! want Domain(2) = -Domain(2); Boundary(2) = inface
179  case (6)
180  ! Subcell top (cell top)
181  inface = npolyverts + 3 ! want Domain(2) = -Domain(2); Boundary(2) = inface
182  end select
183 
184  if (inface .eq. -1) then
185  particle%iboundary(2) = 0
186  else if (inface .eq. 0) then
187  particle%iboundary(2) = 0
188  else
189  if ((inface .ge. 1) .and. (inface .le. 4)) then
190  ! Account for local cell rotation
191  inface = inface + cell%irvOrigin - 1
192  if (inface .gt. 4) inface = inface - 4
193  inface = cell%irectvert(inface) + infaceoff
194  if (inface .lt. 1) inface = inface + npolyverts
195  end if
196  particle%iboundary(2) = inface
197  end if
198  end select