MODFLOW 6  version 6.7.0.dev2
USGS Modular Hydrologic Model
MethodCellPassToBot.f90
Go to the documentation of this file.
2 
3  use kindmodule, only: dp, i4b
4  use methodmodule, only: level_feature
7  use prtfmimodule, only: prtfmitype
8  use basedismodule, only: disbasetype
9  use dismodule, only: distype
10  use disvmodule, only: disvtype
11  use particlemodule, only: particletype
12  use cellmodule, only: celltype
13  use subcellmodule, only: subcelltype
14  implicit none
15 
16  private
17  public :: methodcellpasstobottype
18  public :: create_method_cell_ptb
19 
21  private
22  contains
23  procedure, public :: apply => apply_ptb
24  procedure, public :: deallocate
26 
27 contains
28 
29  !> @brief Create a new pass-to-bottom tracking method
30  subroutine create_method_cell_ptb(method)
31  type(methodcellpasstobottype), pointer :: method
32  allocate (method)
33  allocate (method%name)
34  method%name = "passtobottom"
35  method%delegates = .false.
36  end subroutine create_method_cell_ptb
37 
38  !> @brief Deallocate the pass-to-bottom tracking method
39  subroutine deallocate (this)
40  class(methodcellpasstobottype), intent(inout) :: this
41  deallocate (this%name)
42  end subroutine deallocate
43 
44  !> @brief Pass particle vertically and instantaneously to the cell bottom
45  subroutine apply_ptb(this, particle, tmax)
46  use particlemodule, only: term_no_exits
47  ! dummy
48  class(methodcellpasstobottype), intent(inout) :: this
49  type(particletype), pointer, intent(inout) :: particle
50  real(DP), intent(in) :: tmax
51  ! local
52  integer(I4B) :: nlay
53 
54  call this%assess(particle, this%cell%defn, tmax)
55  if (.not. particle%advancing) return
56 
57  ! Pass to bottom face
58  particle%z = this%cell%defn%bot
59  particle%iboundary(level_feature) = this%cell%defn%npolyverts + 2
60 
61  ! Terminate if in bottom layer
62  select type (dis => this%fmi%dis)
63  type is (distype)
64  nlay = dis%nlay
65  type is (disvtype)
66  nlay = dis%nlay
67  end select
68  if (particle%ilay == nlay) then
69  call this%terminate(particle, status=term_no_exits)
70  else
71  call this%cellexit(particle)
72  end if
73  end subroutine apply_ptb
74 
subroutine, public create_defn(cellDefn)
Create a new cell definition object.
Definition: CellDefn.f90:42
Definition: Dis.f90:1
This module defines variable data types.
Definition: kind.f90:8
subroutine, public create_method_cell_ptb(method)
Create a new pass-to-bottom tracking method.
subroutine apply_ptb(this, particle, tmax)
Pass particle vertically and instantaneously to the cell bottom.
Particle tracking strategies.
Definition: Method.f90:2
@, public level_feature
Definition: Method.f90:37
@ term_no_exits
terminated in a cell with no exit face
Definition: Particle.f90:32
Base grid cell definition.
Definition: CellDefn.f90:10
Base type for grid cells of a concrete type. Contains a cell-definition which is information shared b...
Definition: Cell.f90:11
Structured grid discretization.
Definition: Dis.f90:23
Vertex grid discretization.
Definition: Disv.f90:24
Particle tracked by the PRT model.
Definition: Particle.f90:56
A subcell of a cell.
Definition: Subcell.f90:9