35 use mpi_f08,
only: mpi_wtime
64 type(
case_t),
intent(inout) :: c
66 character(len=LOG_SIZE) :: log_buf
68 character(len=:),
allocatable :: restart_file
71 call c%params%get(
'case.restart_file', restart_file, found)
72 if (found .and. len_trim(restart_file) .gt. 0)
then
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,
']'
82 if (.not. dt_controller%is_variable_dt)
then
83 write(log_buf,
'(A, E15.7)')
'dt : ', c%time%dt
85 write(log_buf,
'(A, E15.7)')
'CFL : ', dt_controller%cfl_trg
90 call neko_log%section(
'Preprocessing')
91 call c%output_controller%execute(c%time)
93 call c%user%initialize(c%time)
101 type(
case_t),
intent(inout) :: c
102 logical :: output_at_end
106 output_at_end, .true.)
107 call c%output_controller%execute(c%time, output_at_end)
109 if (.not. (output_at_end) .and. c%time%t .lt. c%time%end_time)
then
114 call c%user%finalize(c%time)
116 call neko_log%end_section(
'Normal end.')
122 type(
case_t),
intent(inout) :: c
124 real(kind=
dp),
optional,
intent(in) :: tstep_loop_start_time
125 real(kind=
dp) :: start_time, end_time, tstep_start_time
127 character(len=LOG_SIZE) :: log_buf
131 start_time = mpi_wtime()
132 tstep_start_time = start_time
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)
144 call neko_log%section(
'Preprocessing')
146 write(log_buf,
"(A,F8.4,2x,A,1P,E14.7)")
'CFL:', cfl,
'dt:', c%time%dt
150 call c%user%preprocess(c%time)
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
163 if (
allocated(c%scalars))
then
164 start_time = mpi_wtime()
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
173 call neko_log%section(
'Postprocessing')
179 call c%user%compute(c%time)
182 call c%output_controller%execute(c%time)
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
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
212 if (
allocated(ext_bdf))
then
214 time%tlag(i) = time%tlag(i-1)
215 time%dtlag(i) = time%dtlag(i-1)
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
225 call ext_bdf%set_coeffs(time%dtlag)
228 time%tstep = time%tstep + 1
229 time%t = time%t + time%dt
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
241 logical :: found, check_cont
244 call json_get(c%params,
'case.restart_file', restart_file)
246 restart_mesh_file,
"")
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)
253 c%chkp%mesh2mesh_tol, 1e-6_rp)
256 call neko_log%section(
'Restarting from checkpoint')
257 write(log_buf,
'(A,A)')
'File : ', trim(restart_file)
260 call chkpf%init(trim(restart_file))
261 call chkpf%read(c%chkp)
266 call c%chkp%previous_mesh%free()
268 write(log_buf,
'(A,E15.7)')
'Time : ', c%time%t
276 type(
case_t),
intent(inout) :: C
277 type(
chkp_t),
intent(inout) :: chkp
278 character(len=LOG_SIZE) :: log_buf
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)
288 call c%fluid%restart(chkp)
289 if (
allocated(c%scalars))
call c%scalars%restart(chkp)
292 call c%output_controller%set_counter(c%time)
301 type(
case_t),
intent(inout) :: C
302 real(kind=
rp),
intent(inout) :: t
304 character(len=:),
allocatable :: chkp_format
305 character(len=LOG_SIZE) :: log_buf
306 character(len=10) :: format_str
311 call c%chkp%sync_host()
312 if (chkp_format .eq.
'hdf5')
then
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 >>>'
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.
Module for file I/O operations.
logical function, public jobctrl_time_limit()
Check if the job's time limit has been reached.
Utilities for retrieving parameters from the case files.
type(log_t), public neko_log
Global log stream.
integer, parameter, public log_size
integer, parameter, public dp
integer, parameter, public rp
Global precision used in computations.
subroutine, public profiler_start
Start profiling.
subroutine, public profiler_start_region(name, region_id)
Started a named (name) profiler region.
subroutine, public profiler_end_region(name, region_id)
End the most recently started profiler region.
subroutine, public profiler_stop
Stop profiling.
Contains the simcomp_executor_t type.
type(simcomp_executor_t), target, public neko_simcomps
Global variable for the simulation component 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...
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.
Provides a tool to set time step dt.