57 type(
case_t),
target,
intent(inout) :: c
58 real(kind=
rp) :: t, cfl
59 real(kind=
dp) :: start_time_org, start_time, end_time, tstep_start_time
60 character(len=LOG_SIZE) :: log_buf
62 character(len=:),
allocatable :: restart_file
63 logical :: output_at_end, found
65 real(kind=
rp) :: cfl_avrg = 0.0_rp
70 call neko_log%section(
'Starting simulation')
71 write(log_buf,
'(A, E15.7,A,E15.7,A)')
'T : [', 0d0,
',', c%end_time,
')'
73 call dt_controller%init(c%params)
74 if (.not. dt_controller%if_variable_dt)
then
75 write(log_buf,
'(A, E15.7)')
'dt : ', c%dt
78 write(log_buf,
'(A, E15.7)')
'CFL : ', dt_controller%set_cfl
82 call c%params%get(
'case.restart_file', restart_file, found)
83 if (found .and. len_trim(restart_file) .gt. 0)
then
92 call neko_log%section(
'Postprocessing')
93 call c%q%eval(t, c%dt, tstep)
94 call c%s%sample(t, tstep)
96 call c%usr%user_init_modules(t, c%fluid%u, c%fluid%v, c%fluid%w,&
97 c%fluid%p, c%fluid%c_Xh, c%params)
102 cfl = c%fluid%compute_cfl(c%dt)
103 start_time_org = mpi_wtime()
108 start_time = mpi_wtime()
109 tstep_start_time = start_time
110 if (dt_controller%dt_last_change .eq. 0)
then
113 call dt_controller%set_dt(c%dt, cfl, cfl_avrg, tstep)
115 cfl = c%fluid%compute_cfl(c%dt)
118 write(log_buf,
'(A,I6)')
'Time-step: ', tstep
122 write(log_buf,
'(A,E15.7,1x,A,E15.7)')
'CFL:', cfl,
'dt:', c%dt
127 call neko_log%section(
'Preprocessing')
132 call c%fluid%step(t, tstep, c%dt, c%ext_bdf, dt_controller)
133 end_time = mpi_wtime()
134 write(log_buf,
'(A,E15.7)') &
135 'Fluid step time (s): ', end_time-start_time
137 write(log_buf,
'(A,E15.7)') &
138 'Total elapsed time (s):', end_time-start_time_org
142 if (
allocated(c%scalar))
then
143 start_time = mpi_wtime()
145 call c%scalar%step(t, tstep, c%dt, c%ext_bdf, dt_controller)
146 end_time = mpi_wtime()
147 write(log_buf,
'(A,E15.7)') &
148 'Scalar step time: ', end_time-start_time
150 write(log_buf,
'(A,E15.7)') &
151 'Total elapsed time (s):', end_time-start_time_org
155 call neko_log%section(
'Postprocessing')
159 call c%q%eval(t, c%dt, tstep)
160 call c%s%sample(t, tstep)
163 call c%usr%material_properties(t, tstep, c%material_properties%rho,&
164 c%material_properties%mu, &
165 c%material_properties%cp, &
166 c%material_properties%lambda, &
169 call c%usr%user_check(t, tstep, c%fluid%u, c%fluid%v, c%fluid%w, &
170 c%fluid%p, c%fluid%c_Xh, c%params)
173 end_time = mpi_wtime()
174 call neko_log%section(
'Step summary')
175 write(log_buf,
'(A,I8,A,E15.7)') &
176 'Total time for step ', tstep,
' (s): ', end_time-tstep_start_time
178 write(log_buf,
'(A,E15.7)') &
179 'Total elapsed time (s): ', end_time-start_time_org
188 output_at_end, .true.)
189 call c%s%sample(t, tstep, output_at_end)
191 if (.not. (output_at_end) .and. t .lt. c%end_time)
then
195 call c%usr%user_finalize_modules(t, c%params)
197 call neko_log%end_section(
'Normal end.')
202 real(kind=
rp),
intent(inout) :: t
203 real(kind=
rp),
intent(in) :: dt
205 real(kind=
rp),
dimension(10) :: tlag
206 real(kind=
rp),
dimension(10) :: dtlag
207 integer,
intent(in) :: step
213 dtlag(i) = dtlag(i-1)
218 if (ext_bdf%ndiff .eq. 0)
then
225 call ext_bdf%set_coeffs(dtlag)
232 type(
case_t),
intent(inout) :: C
233 real(kind=
rp),
intent(inout) :: t
235 type(
file_t) :: chkpf, previous_meshf
236 character(len=LOG_SIZE) :: log_buf
237 character(len=:),
allocatable :: restart_file
238 character(len=:),
allocatable :: restart_mesh_file
242 call c%params%get(
'case.restart_file', restart_file, found)
243 call c%params%get(
'case.restart_mesh_file', restart_mesh_file,&
247 previous_meshf =
file_t(trim(restart_mesh_file))
248 call previous_meshf%read(c%fluid%chkp%previous_mesh)
251 call c%params%get(
'case.mesh2mesh_tolerance', tol,&
254 if (found) c%fluid%chkp%mesh2mesh_tol = tol
256 chkpf =
file_t(trim(restart_file))
257 call chkpf%read(c%fluid%chkp)
258 c%dtlag = c%fluid%chkp%dtlag
259 c%tlag = c%fluid%chkp%tlag
262 call c%fluid%chkp%previous_mesh%free()
263 do i = 1,
size(c%dtlag)
264 call c%ext_bdf%set_coeffs(c%dtlag)
267 call c%fluid%restart(c%dtlag, c%tlag)
268 if (
allocated(c%scalar))
call c%scalar%restart( c%dtlag, c%tlag)
270 t = c%fluid%chkp%restart_time()
271 call neko_log%section(
'Restarting from checkpoint')
272 write(log_buf,
'(A,A)')
'File : ', &
275 write(log_buf,
'(A,E15.7)')
'Time : ', t
279 call c%s%set_counter(t)
284 type(
case_t),
intent(inout) :: C
285 real(kind=
rp),
intent(inout) :: t
287 character(len=:),
allocatable :: chkp_format
288 character(len=LOG_SIZE) :: log_buf
289 character(len=10) :: format_str
292 call c%params%get(
'case.checkpoint_format', chkp_format, found)
293 call c%fluid%chkp%sync_host()
296 if (chkp_format .eq.
'hdf5')
then
300 chkpf =
file_t(
'joblimit'//trim(format_str))
301 call chkpf%write(c%fluid%chkp, t)
302 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...
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 simulation_settime(t, dt, ext_bdf, tlag, dtlag, step)
subroutine, public neko_solve(C)
Main driver to solve a case C.
subroutine simulation_joblimit_chkp(C, t)
Write a checkpoint at joblimit.
subroutine simulation_restart(C, t)
Restart a case C from a given checkpoint.
Compound scheme for the advection and diffusion operators in a transport equation.
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 ...
Provides a tool to set time step dt.