Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
simulation.f90
Go to the documentation of this file.
1! Copyright (c) 2020-2021, 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 mpi_f08, only: mpi_wtime
36 use case, only : case_t
37 use checkpoint, only : chkp_t
38 use num_types, only : rp, dp
40 use file, only : file_t
41 use logger, only : log_size, neko_log
42 use jobctrl, only : jobctrl_time_limit
47 use time_state, only : time_state_t
49 implicit none
50 private
51
55 end interface simulation_restart
56
59
60contains
61
63 subroutine simulation_init(C, dt_controller)
64 type(case_t), intent(inout) :: c
65 type(time_step_controller_t), intent(inout) :: dt_controller
66 character(len=LOG_SIZE) :: log_buf
67 logical :: found
68 character(len=:), allocatable :: restart_file
69
70 ! Restart the case if needed
71 call c%params%get('case.restart_file', restart_file, found)
72 if (found .and. len_trim(restart_file) .gt. 0) then
73 ! Restart the case
74 call simulation_restart(c)
75 end if
76
77 ! Write the initial logging message
78 call neko_log%section('Starting simulation')
79 write(log_buf, '(A, E15.7,A,E15.7,A)') &
80 'T : [', c%time%t, ',', c%time%end_time, ']'
81 call neko_log%message(log_buf)
82 if (.not. dt_controller%is_variable_dt) then
83 write(log_buf, '(A, E15.7)') 'dt : ', c%time%dt
84 else
85 write(log_buf, '(A, E15.7)') 'CFL : ', dt_controller%cfl_trg
86 end if
87 call neko_log%message(log_buf)
88
89 ! Execute outputs and user-init before time loop
90 call neko_log%section('Preprocessing')
91 call c%output_controller%execute(c%time)
92
93 call c%user%initialize(c%time)
94 call neko_log%end_section()
95 call neko_log%newline()
96
97 end subroutine simulation_init
98
101 type(case_t), intent(inout) :: c
102 logical :: output_at_end
103
104 ! Run a final output if specified in the json
105 call json_get_or_default(c%params, 'case.output_at_end', &
106 output_at_end, .true.)
107 call c%output_controller%execute(c%time, output_at_end)
108
109 if (.not. (output_at_end) .and. c%time%t .lt. c%time%end_time) then
110 call simulation_joblimit_chkp(c, c%time%t)
111 end if
112
113 ! Finalize the user modules
114 call c%user%finalize(c%time)
115
116 call neko_log%end_section('Normal end.')
117
118 end subroutine simulation_finalize
119
121 subroutine simulation_step(C, dt_controller, tstep_loop_start_time)
122 type(case_t), intent(inout) :: c
123 type(time_step_controller_t), intent(inout) :: dt_controller
124 real(kind=dp), optional, intent(in) :: tstep_loop_start_time
125 real(kind=dp) :: start_time, end_time, tstep_start_time
126 real(kind=rp) :: cfl
127 character(len=LOG_SIZE) :: log_buf
128
129 ! Setup the time step, and start time
130 call profiler_start_region('Time-Step')
131 start_time = mpi_wtime()
132 tstep_start_time = start_time
133
134 ! Compute the next time step size
135 cfl = c%fluid%compute_cfl(c%time%dt)
136 call dt_controller%set_dt(c%time, cfl)
137 if (dt_controller%is_variable_dt) cfl = c%fluid%compute_cfl(c%time%dt)
138
139 ! Advance time step from t to t+dt and print the status
140 call simulation_settime(c%time, c%fluid%ext_bdf)
141 call c%time%status()
142 call neko_log%begin()
143
144 call neko_log%section('Preprocessing')
145
146 write(log_buf, "(A,F8.4,2x,A,1P,E14.7)") 'CFL:', cfl, 'dt:', c%time%dt
147 call neko_log%message(log_buf)
148
149 ! Run the preprocessing
150 call c%user%preprocess(c%time)
151 call neko_simcomps%preprocess(c%time)
152 call neko_log%end_section()
153
154 ! Fluid step
155 call neko_log%section('Fluid')
156 start_time = mpi_wtime()
157 call c%fluid%step(c%time, dt_controller)
158 end_time = mpi_wtime()
159 write(log_buf, '(A,3X,E15.7)') 'Fluid step time (s):', end_time - start_time
160 call neko_log%end_section(log_buf)
161
162 ! Scalar step
163 if (allocated(c%scalars)) then
164 start_time = mpi_wtime()
165 call neko_log%section('Scalar')
166 call c%scalars%step(c%time, c%fluid%ext_bdf, dt_controller)
167 end_time = mpi_wtime()
168 write(log_buf, '(A,6X,E15.7)') 'Scalar step time:', end_time - start_time
169 call neko_log%end_section(log_buf)
170 end if
171
172 ! Postprocessing
173 call neko_log%section('Postprocessing')
174
175 ! Execute all simulation components
176 call neko_simcomps%compute(c%time)
177
178 ! Run the user compute routine
179 call c%user%compute(c%time)
180
181 ! Run any IO needed.
182 call c%output_controller%execute(c%time)
183
184 call neko_log%end_section()
185
186 ! End the step and print summary
187 end_time = mpi_wtime()
188 call neko_log%section('Step summary')
189 write(log_buf, '(A20,I8,A6,E15.7)') &
190 'Total time for step ', c%time%tstep, ' (s): ', &
191 end_time-tstep_start_time
192 call neko_log%message(log_buf)
193
194 if (present(tstep_loop_start_time)) then
195 write(log_buf, '(A34,E15.7)') &
196 'Total elapsed time (s): ', &
197 end_time - tstep_loop_start_time
198 call neko_log%message(log_buf)
199 end if
200
201 call neko_log%end_section()
202 call neko_log%end()
204
205 end subroutine simulation_step
206
207 subroutine simulation_settime(time, ext_bdf)
208 type(time_state_t), intent(inout) :: time
209 type(time_scheme_controller_t), intent(inout), allocatable :: ext_bdf
210 integer :: i
211
212 if (allocated(ext_bdf)) then
213 do i = 10, 2, -1
214 time%tlag(i) = time%tlag(i-1)
215 time%dtlag(i) = time%dtlag(i-1)
216 end do
217
218 time%dtlag(1) = time%dt
219 time%tlag(1) = time%t
220 if (ext_bdf%ndiff .eq. 0) then
221 time%dtlag(2) = time%dt
222 time%tlag(2) = time%t
223 end if
224
225 call ext_bdf%set_coeffs(time%dtlag)
226 end if
227
228 time%tstep = time%tstep + 1
229 time%t = time%t + time%dt
230
231 end subroutine simulation_settime
232
235 type(case_t), intent(inout) :: C
236 type(file_t) :: chkpf, previous_meshf
237 character(len=LOG_SIZE) :: log_buf
238 character(len=:), allocatable :: restart_file
239 character(len=:), allocatable :: restart_mesh_file
240 real(kind=rp) :: tol
241 logical :: found, check_cont
242 integer :: i
243
244 call json_get(c%params, 'case.restart_file', restart_file)
245 call json_get_or_default(c%params, 'case.restart_mesh_file', &
246 restart_mesh_file, "")
247
248 if (restart_mesh_file .ne. "") then
249 call previous_meshf%init(trim(restart_mesh_file))
250 call previous_meshf%read(c%chkp%previous_mesh)
251
252 call json_get_or_default(c%params, 'case.mesh2mesh_tolerance', &
253 c%chkp%mesh2mesh_tol, 1e-6_rp)
254 end if
255
256 call neko_log%section('Restarting from checkpoint')
257 write(log_buf, '(A,A)') 'File : ', trim(restart_file)
258 call neko_log%message(log_buf)
259
260 call chkpf%init(trim(restart_file))
261 call chkpf%read(c%chkp)
262
263 call case_restart_from_checkpoint(c, c%chkp)
264
265 ! Free the previous mesh
266 call c%chkp%previous_mesh%free()
267
268 write(log_buf, '(A,E15.7)') 'Time : ', c%time%t
269 call neko_log%message(log_buf)
270 call neko_log%end_section()
271
272 end subroutine case_restart_from_parameters
273
276 type(case_t), intent(inout) :: C
277 type(chkp_t), intent(inout) :: chkp
278 character(len=LOG_SIZE) :: log_buf
279 integer :: i
280
281 ! Restart the time state and BDF coefficients
282 call c%time%restart(chkp)
283 do i = 1, size(c%time%dtlag)
284 call c%fluid%ext_bdf%set_coeffs(c%time%dtlag)
285 end do
286
287 ! Restart the fluid and scalars
288 call c%fluid%restart(chkp)
289 if (allocated(c%scalars)) call c%scalars%restart(chkp)
290
291 ! Restart the output controller
292 call c%output_controller%set_counter(c%time)
293
294 ! Restart the simulation components
295 call neko_simcomps%restart(c%time)
296
297 end subroutine case_restart_from_checkpoint
298
301 type(case_t), intent(inout) :: C
302 real(kind=rp), intent(inout) :: t
303 type(file_t) :: chkpf
304 character(len=:), allocatable :: chkp_format
305 character(len=LOG_SIZE) :: log_buf
306 character(len=10) :: format_str
307 logical :: found
308
309 call json_get_or_default(c%params, 'case.checkpoint_format', chkp_format, &
310 'default')
311 call c%chkp%sync_host()
312 if (chkp_format .eq. 'hdf5') then
313 format_str = '.h5'
314 else
315 format_str = '.chkp'
316 end if
317 call chkpf%init(c%output_directory // 'joblimit'//trim(format_str))
318 call chkpf%write(c%chkp, t)
319 write(log_buf, '(A)') '! saving checkpoint >>>'
320 call neko_log%message(log_buf)
321
322 end subroutine simulation_joblimit_chkp
323
324end module simulation
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
Retrieves a parameter by name or throws an error.
Defines a simulation case.
Definition case.f90:34
Defines a checkpoint.
Module for file I/O operations.
Definition file.f90:34
Job control.
Definition jobctrl.f90:34
logical function, public jobctrl_time_limit()
Check if the job's time limit has been reached.
Definition jobctrl.f90:108
Utilities for retrieving parameters from the case files.
Logging routines.
Definition log.f90:34
type(log_t), public neko_log
Global log stream.
Definition log.f90:70
integer, parameter, public log_size
Definition log.f90:46
integer, parameter, public dp
Definition num_types.f90:9
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Profiling interface.
Definition profiler.F90:34
subroutine, public profiler_start
Start profiling.
Definition profiler.F90:52
subroutine, public profiler_start_region(name, region_id)
Started a named (name) profiler region.
Definition profiler.F90:78
subroutine, public profiler_end_region(name, region_id)
End the most recently started profiler region.
Definition profiler.F90:109
subroutine, public profiler_stop
Stop profiling.
Definition profiler.F90:65
Contains the simcomp_executor_t type.
type(simcomp_executor_t), target, public neko_simcomps
Global variable for the simulation component driver.
Simulation driver.
subroutine, public simulation_step(c, dt_controller, tstep_loop_start_time)
Compute a single time-step of a case.
subroutine simulation_settime(time, ext_bdf)
subroutine, public simulation_finalize(c)
Finalize a simulation of a case.
subroutine case_restart_from_checkpoint(c, chkp)
Restart a case C from a given checkpoint.
subroutine, public simulation_init(c, dt_controller)
Initialise a simulation of a case.
subroutine simulation_joblimit_chkp(c, t)
Write a checkpoint at joblimit.
subroutine case_restart_from_parameters(c)
Restart a case C from a given checkpoint.
Compound scheme for the advection and diffusion operators in a transport equation.
Module with things related to the simulation time.
Implements type time_step_controller.
A wrapper around a polymorphic generic_file_t that handles its init. This is essentially a factory fo...
Definition file.f90:55
Implements the logic to compute the time coefficients for the advection and diffusion operators in a ...
A struct that contains all info about the time, expand as needed.