MODFLOW 6  version 6.7.0.dev0
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 201 of file MethodCellPollockQuad.f90.

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

◆ create_method_cell_quad()

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

Definition at line 30 of file MethodCellPollockQuad.f90.

31  ! dummy
32  type(MethodCellPollockQuadType), pointer :: method
33  ! local
34  type(CellRectQuadType), pointer :: cell
35  type(SubcellRectType), pointer :: subcell
36 
37  allocate (method)
38  call create_cell_rect_quad(cell)
39  method%cell => cell
40  method%name => method%cell%type
41  method%delegates = .true.
42  call create_subcell_rect(subcell)
43  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 47 of file MethodCellPollockQuad.f90.

48  class(MethodCellPollockQuadType), intent(inout) :: this
49  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 53 of file MethodCellPollockQuad.f90.

54  class(MethodCellPollockQuadType), intent(inout) :: this
55  type(ParticleType), pointer, intent(inout) :: particle
56  integer, intent(in) :: next_level
57  class(MethodType), pointer, intent(inout) :: submethod
58 
59  select type (subcell => this%subcell)
60  type is (subcellrecttype)
61  call this%load_subcell(particle, subcell)
62  end select
63  call method_subcell_plck%init( &
64  fmi=this%fmi, &
65  cell=this%cell, &
66  subcell=this%subcell, &
67  trackctl=this%trackctl, &
68  tracktimes=this%tracktimes)
69  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 236 of file MethodCellPollockQuad.f90.

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

◆ pass_mcpq()

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

Definition at line 73 of file MethodCellPollockQuad.f90.

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