Neko  0.9.0
A portable framework for high-order spectral element flow simulations
time_based_controller.f90
Go to the documentation of this file.
1 ! Copyright (c) 2023, The Neko Authors
2 ! All rights reserved.
3 !
4 ! Redistribution and use in source and binary forms, with or without
5 ! modification, are permitted provided that the following conditions
6 ! are met:
7 !
8 ! * Redistributions of source code must retain the above copyright
9 ! notice, this list of conditions and the following disclaimer.
10 !
11 ! * Redistributions in binary form must reproduce the above
12 ! copyright notice, this list of conditions and the following
13 ! disclaimer in the documentation and/or other materials provided
14 ! with the distribution.
15 !
16 ! * Neither the name of the authors nor the names of its
17 ! contributors may be used to endorse or promote products derived
18 ! from this software without specific prior written permission.
19 !
20 ! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21 ! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22 ! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23 ! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24 ! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25 ! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26 ! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27 ! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28 ! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29 ! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30 ! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31 ! POSSIBILITY OF SUCH DAMAGE.
32 !
35  use num_types, only : rp
36  use utils, only : neko_error
37  implicit none
38  private
39 
47  type, public :: time_based_controller_t
49  real(kind=rp) :: frequency = 0
51  real(kind=rp) :: time_interval = 0
53  integer :: nsteps = 0
55  real(kind=rp) :: end_time = 0
57  integer :: nexecutions = 0
59  logical :: never = .false.
62  character(len=:), allocatable :: control_mode
64  real(kind=rp) :: control_value
65 
66  contains
68  procedure, pass(this) :: init => time_based_controller_init
70  procedure, pass(this) :: check => time_based_controller_check
72  procedure, pass(this) :: register_execution => &
75  procedure, pass(this) :: set_counter => &
77 
79 
80  interface assignment(=)
81  module procedure time_based_controller_assignment
82  end interface assignment(=)
83 
84 contains
85 
90  subroutine time_based_controller_init(this, end_time, control_mode, &
91  control_value)
92  class(time_based_controller_t), intent(inout) :: this
93  real(kind=rp), intent(in) :: end_time
94  character(len=*), intent(in) :: control_mode
95  real(kind=rp), intent(in) :: control_value
96 
97  this%end_time = end_time
98  this%control_mode = control_mode
99  this%control_value = control_value
100 
101  if (trim(control_mode) .eq. 'simulationtime') then
102  this%time_interval = control_value
103  this%frequency = 1/this%time_interval
104  this%nsteps = 0
105  else if (trim(control_mode) .eq. 'nsamples') then
106  if (control_value .le. 0) then
107  call neko_error("nsamples must be positive")
108  end if
109 
110  this%frequency = control_value / end_time
111  this%time_interval = 1.0_rp / this%frequency
112  this%nsteps = 0
113  else if (trim(control_mode) .eq. 'tsteps') then
114  this%nsteps = control_value
115  ! if the timestep will be variable, we cannot compute these.
116  this%frequency = 0
117  this%time_interval = 0
118  else if (trim(control_mode) .eq. 'never') then
119  this%never = .true.
120  else
121  call neko_error("The control parameter must be simulationtime, nsamples&
122  & tsteps, or never, but received "//trim(control_mode))
123  end if
124  end subroutine time_based_controller_init
125 
133  function time_based_controller_check(this, t, tstep, force) result(check)
134  class(time_based_controller_t), intent(inout) :: this
135  real(kind=rp), intent(in) :: t
136  integer, intent(in) :: tstep
137  logical, intent(in), optional :: force
138  logical :: check
139  logical :: ifforce
140 
141  if (present(force)) then
142  ifforce = force
143  else
144  ifforce = .false.
145  end if
146 
147  check = .false.
148  if (ifforce) then
149  check = .true.
150  else if (this%never) then
151  check = .false.
152  else if ( (this%nsteps .eq. 0) .and. &
153  (t .ge. this%nexecutions * this%time_interval) ) then
154  check = .true.
155  else if (this%nsteps .gt. 0) then
156  if (mod(tstep, this%nsteps) .eq. 0) then
157  check = .true.
158  end if
159  end if
160  end function time_based_controller_check
161 
165  subroutine time_based_controller_assignment(ctrl1, ctrl2)
166  type(time_based_controller_t), intent(inout) :: ctrl1
167  type(time_based_controller_t), intent(in) :: ctrl2
168 
169  ctrl1%end_time = ctrl2%end_time
170  ctrl1%frequency = ctrl2%frequency
171  ctrl1%nsteps = ctrl2%nsteps
172  ctrl1%time_interval = ctrl2%time_interval
173  ctrl1%nexecutions = ctrl2%nexecutions
174 
175  end subroutine time_based_controller_assignment
176 
179  class(time_based_controller_t), intent(inout) :: this
180 
181  this%nexecutions = this%nexecutions + 1
182 
184 
188  class(time_based_controller_t), intent(inout) :: this
189  real(kind=rp), intent(in) :: t
190 
191  if (this%nsteps .eq. 0) then
192  this%nexecutions = int(t / this%time_interval) + 1
193  end if
194 
195  end subroutine time_based_controller_set_counter
196 
197 
198 end module time_based_controller
integer, parameter, public rp
Global precision used in computations.
Definition: num_types.f90:12
Contains the time_based_controller_t type.
subroutine time_based_controller_init(this, end_time, control_mode, control_value)
Constructor.
subroutine time_based_controller_assignment(ctrl1, ctrl2)
Assignment operator. Simply copies attribute values.
subroutine time_based_controller_set_counter(this, t)
Set the counter based on a time (for restarts)
subroutine time_based_controller_register_execution(this)
Increment nexectutions.
logical function time_based_controller_check(this, t, tstep, force)
Check if the execution should be performed.
Utilities.
Definition: utils.f90:35
A utility type for determening whether an action should be executed based on the current time value....