MODFLOW 6  version 6.8.0.dev0
USGS Modular Hydrologic Model
ParticleEvents.f90
Go to the documentation of this file.
2  use kindmodule, only: dp, i4b, lgp
3  use listmodule, only: listtype
6  implicit none
7 
8  private
9  public :: handle_event
10 
11  !> @brief Subscription to particle events: a procedure to
12  !! handle the event with an unlimited pointer for storing
13  !! arbitrary context the handling procedure may reference.
15  procedure(handle_event), pointer, nopass :: handler
16  class(*), pointer :: context
18 
19  !> @brief Dispatcher for particle events. Consumers subscribe
20  !! handlers to the dispatcher. Events may be dispatched, with
21  !! the first handler to handle the event stopping propagation,
22  !! or broadcast, with all handlers receiving the event.
24  type(listtype) :: subscriptions
25  contains
26  procedure, public :: subscribe
27  procedure, public :: dispatch
28  procedure, public :: broadcast
29  procedure :: prep_event
30  procedure :: destroy
32 
33  abstract interface
34  !> @brief Event handler interface. Handlers may signal
35  !! to the dispatching caller whether they have handled
36  !! the event, but the signal is ignored for broadcasts.
37  logical function handle_event(context, particle, event)
39  class(*), pointer :: context
40  type(particletype), pointer, intent(inout) :: particle
41  class(particleeventtype), pointer, intent(in) :: event
42  end function
43  end interface
44 
45 contains
46  !> @brief Add a subscription to the dispatcher.
47  subroutine subscribe(this, handler, context)
48  class(particleeventdispatchertype), intent(inout) :: this
49  procedure(handle_event) :: handler
50  class(*), pointer :: context
51  ! local
52  type(particleeventsubscriptiontype), pointer :: subscription
53  class(*), pointer :: p
54 
55  allocate (subscription)
56  subscription%handler => handler
57  subscription%context => context
58  p => subscription
59  call this%subscriptions%Add(p)
60  end subroutine subscribe
61 
62  !> @brief Prepare an event for dispatching, loading it with
63  !! the current state of the particle. For internal use only.
64  subroutine prep_event(this, particle, event)
65  use tdismodule, only: kper, kstp
66  ! dummy
67  class(particleeventdispatchertype), intent(inout) :: this
68  type(particletype), pointer, intent(inout) :: particle
69  class(particleeventtype), pointer, intent(in) :: event
70  ! local
71  real(DP) :: x, y, z
72 
73  x = particle%x
74  y = particle%y
75  z = particle%z
76  ! switch to model coordinates if needed
77  call particle%get_model_coords(x, y, z)
78 
79  event%kper = kper
80  event%kstp = kstp
81  event%imdl = particle%imdl
82  event%iprp = particle%iprp
83  event%irpt = particle%irpt
84  event%ilay = particle%ilay
85  event%icu = particle%icu
86  event%izone = particle%izone
87  event%trelease = particle%trelease
88  event%ttrack = particle%ttrack
89  event%x = x
90  event%y = y
91  event%z = z
92  event%istatus = particle%istatus
93  end subroutine prep_event
94 
95  !> @brief Dispatch an event for handling. The first
96  !! subscriber to handle the event stops propagation.
97  subroutine dispatch(this, particle, event)
98  class(particleeventdispatchertype), intent(inout) :: this
99  type(particletype), pointer, intent(inout) :: particle
100  class(particleeventtype), pointer, intent(in) :: event
101  ! local
102  logical(LGP) :: handled
103  integer(I4B) :: i
104  class(*), pointer :: p
105 
106  call this%prep_event(particle, event)
107 
108  do i = 1, this%subscriptions%Count()
109  p => this%subscriptions%GetItem(i)
110  select type (subscription => p)
112  handled = subscription%handler( &
113  subscription%context, &
114  particle, &
115  event)
116  if (handled) exit
117  end select
118  end do
119  end subroutine dispatch
120 
121  !> @brief Broadcast an event to all subscribers so
122  !! all receive the event and a chance to handle it.
123  subroutine broadcast(this, particle, event)
124  class(particleeventdispatchertype), intent(inout) :: this
125  type(particletype), pointer, intent(inout) :: particle
126  class(particleeventtype), pointer, intent(in) :: event
127  ! local
128  logical(LGP) :: handled
129  integer(I4B) :: i
130  class(*), pointer :: p
131 
132  call this%prep_event(particle, event)
133 
134  do i = 1, this%subscriptions%Count()
135  p => this%subscriptions%GetItem(i)
136  select type (subscription => p)
138  handled = subscription%handler(subscription%context, particle, event)
139  end select
140  end do
141  end subroutine broadcast
142 
143  !> @brief Destroy the dispatcher.
144  subroutine destroy(this)
145  class(particleeventdispatchertype), intent(inout) :: this
146  call this%subscriptions%Clear(destroy=.true.)
147  end subroutine destroy
148 
149 end module particleeventsmodule
Event handler interface. Handlers may signal to the dispatching caller whether they have handled the ...
This module defines variable data types.
Definition: kind.f90:8
subroutine subscribe(this, handler, context)
Add a subscription to the dispatcher.
subroutine prep_event(this, particle, event)
Prepare an event for dispatching, loading it with the current state of the particle....
subroutine broadcast(this, particle, event)
Broadcast an event to all subscribers so all receive the event and a chance to handle it.
subroutine dispatch(this, particle, event)
Dispatch an event for handling. The first subscriber to handle the event stops propagation.
subroutine destroy(this)
Destroy the dispatcher.
integer(i4b), pointer, public kstp
current time step number
Definition: tdis.f90:24
integer(i4b), pointer, public kper
current stress period number
Definition: tdis.f90:23
A generic heterogeneous doubly-linked list.
Definition: List.f90:14
Base type for particle events.
Dispatcher for particle events. Consumers subscribe handlers to the dispatcher. Events may be dispatc...
Subscription to particle events: a procedure to handle the event with an unlimited pointer for storin...
Particle tracked by the PRT model.
Definition: Particle.f90:62