MODFLOW 6  version 6.7.0.dev3
USGS Modular Hydrologic Model
MethodCellPollock.f90
Go to the documentation of this file.
2 
3  use kindmodule, only: dp, i4b
4  use constantsmodule, only: done, dzero
11  use particlemodule, only: particletype
12  implicit none
13 
14  private
15  public :: methodcellpollocktype
17 
19  contains
20  procedure, public :: apply => apply_mcp
21  procedure, public :: deallocate => destroy_mcp
22  procedure, public :: load => load_mcp
23  procedure, public :: load_subcell
24  procedure, public :: pass => pass_mcp
25  end type methodcellpollocktype
26 
27 contains
28 
29  !> @brief Create a tracking method
30  subroutine create_method_cell_pollock(method)
31  ! dummy
32  type(methodcellpollocktype), pointer :: method
33  ! local
34  type(cellrecttype), pointer :: cell
35  type(subcellrecttype), pointer :: subcell
36 
37  allocate (method)
38  call create_cell_rect(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
44  end subroutine create_method_cell_pollock
45 
46  !> @brief Destroy the tracking method
47  subroutine destroy_mcp(this)
48  class(methodcellpollocktype), intent(inout) :: this
49  deallocate (this%name)
50  end subroutine destroy_mcp
51 
52  !> @brief Load subcell tracking method
53  subroutine load_mcp(this, particle, next_level, submethod)
54  ! modules
55  use subcellmodule, only: subcelltype
56  ! dummy
57  class(methodcellpollocktype), intent(inout) :: this
58  type(particletype), pointer, intent(inout) :: particle
59  integer, intent(in) :: next_level
60  class(methodtype), pointer, intent(inout) :: submethod
61 
62  select type (subcell => this%subcell)
63  type is (subcellrecttype)
64  call this%load_subcell(particle, subcell)
65  end select
66  call method_subcell_plck%init( &
67  fmi=this%fmi, &
68  cell=this%cell, &
69  subcell=this%subcell, &
70  events=this%events, &
71  tracktimes=this%tracktimes)
72  submethod => method_subcell_plck
73  particle%itrdomain(next_level) = 1
74  end subroutine load_mcp
75 
76  !> @brief Having exited the lone subcell, pass the particle to the cell face
77  !! In this case the lone subcell is the cell.
78  subroutine pass_mcp(this, particle)
79  ! dummy
80  class(methodcellpollocktype), intent(inout) :: this
81  type(particletype), pointer, intent(inout) :: particle
82  ! local
83  integer(I4B) :: exitface
84  integer(I4B) :: entryface
85 
86  exitface = particle%iboundary(level_subfeature)
87  ! Map subcell exit face to cell face
88  select case (exitface) ! note: exitFace uses Dave's iface convention
89  case (0)
90  entryface = -1
91  case (1)
92  entryface = 1
93  case (2)
94  entryface = 3
95  case (3)
96  entryface = 4
97  case (4)
98  entryface = 2
99  case (5)
100  entryface = 6 ! note: inface=5 same as inface=1 due to wraparound
101  case (6)
102  entryface = 7
103  end select
104  if (entryface .eq. -1) then
105  particle%iboundary(level_feature) = 0
106  else
107  if ((entryface .ge. 1) .and. (entryface .le. 4)) then
108  ! Account for local cell rotation
109  select type (cell => this%cell)
110  type is (cellrecttype)
111  entryface = entryface + cell%ipvOrigin - 1
112  end select
113  if (entryface .gt. 4) entryface = entryface - 4
114  end if
115  particle%iboundary(level_feature) = entryface
116  end if
117  end subroutine pass_mcp
118 
119  !> @brief Apply Pollock's method to a rectangular cell
120  subroutine apply_mcp(this, particle, tmax)
121  ! dummy
122  class(methodcellpollocktype), intent(inout) :: this
123  type(particletype), pointer, intent(inout) :: particle
124  real(DP), intent(in) :: tmax
125  ! local
126  real(DP) :: xOrigin
127  real(DP) :: yOrigin
128  real(DP) :: zOrigin
129  real(DP) :: sinrot
130  real(DP) :: cosrot
131 
132  select type (cell => this%cell)
133  type is (cellrecttype)
134  call this%assess(particle, this%cell%defn, tmax)
135  if (.not. particle%advancing) return
136 
137  ! Transform model coordinates to local cell coordinates
138  ! (translated/rotated but not scaled relative to model)
139  xorigin = cell%xOrigin
140  yorigin = cell%yOrigin
141  zorigin = cell%zOrigin
142  sinrot = cell%sinrot
143  cosrot = cell%cosrot
144  call particle%transform(xorigin, yorigin, zorigin, &
145  sinrot, cosrot)
146 
147  ! Track the particle over the cell
148  call this%track(particle, 2, tmax)
149 
150  ! Transform cell coordinates back to model coordinates
151  call particle%transform(xorigin, yorigin, zorigin, &
152  sinrot, cosrot, invert=.true.)
153  call particle%reset_transform()
154  end select
155  end subroutine apply_mcp
156 
157  !> @brief Loads the lone rectangular subcell from the rectangular cell
158  subroutine load_subcell(this, particle, subcell) !
159  ! dummy
160  class(methodcellpollocktype), intent(inout) :: this
161  type(particletype), pointer, intent(inout) :: particle
162  type(subcellrecttype), intent(inout) :: subcell
163 
164  select type (cell => this%cell)
165  type is (cellrecttype)
166  ! Set cell/subcell numbers
167  subcell%icell = cell%defn%icell
168  subcell%isubcell = 1
169 
170  ! Subcell calculations will be done in local subcell coordinates
171  subcell%dx = cell%dx
172  subcell%dy = cell%dy
173  subcell%dz = cell%dz
174  subcell%sinrot = dzero
175  subcell%cosrot = done
176  subcell%xOrigin = dzero
177  subcell%yOrigin = dzero
178  subcell%zOrigin = dzero
179 
180  ! Set subcell edge velocities
181  subcell%vx1 = cell%vx1 ! cell velocities already account for retfactor and porosity
182  subcell%vx2 = cell%vx2
183  subcell%vy1 = cell%vy1
184  subcell%vy2 = cell%vy2
185  subcell%vz1 = cell%vz1
186  subcell%vz2 = cell%vz2
187  end select
188  end subroutine load_subcell
189 
190 end module methodcellpollockmodule
subroutine, public create_cell_rect(cell)
Create a new rectangular cell.
Definition: CellRect.f90:40
This module contains simulation constants.
Definition: Constants.f90:9
real(dp), parameter dzero
real constant zero
Definition: Constants.f90:65
real(dp), parameter done
real constant 1
Definition: Constants.f90:76
This module defines variable data types.
Definition: kind.f90:8
subroutine pass_mcp(this, particle)
Having exited the lone subcell, pass the particle to the cell face In this case the lone subcell is t...
procedure subroutine, public create_method_cell_pollock(method)
Create a tracking method.
subroutine load_subcell(this, particle, subcell)
Loads the lone rectangular subcell from the rectangular cell.
subroutine destroy_mcp(this)
Destroy the tracking method.
subroutine load_mcp(this, particle, next_level, submethod)
Load subcell tracking method.
subroutine apply_mcp(this, particle, tmax)
Apply Pollock's method to a rectangular cell.
Particle tracking strategies.
Definition: Method.f90:2
@, public level_feature
Definition: Method.f90:41
@, public level_subfeature
Definition: Method.f90:42
Subcell-level tracking methods.
type(methodsubcellpollocktype), pointer, public method_subcell_plck
type(methodsubcellternarytype), pointer, public method_subcell_tern
subroutine, public create_subcell_rect(subcell)
Create a new rectangular subcell.
Definition: SubcellRect.f90:28
Base type for particle tracking methods.
Definition: Method.f90:58
Particle tracked by the PRT model.
Definition: Particle.f90:56
A subcell of a cell.
Definition: Subcell.f90:10