59 type(
case_t),
target,
intent(inout) :: c
60 real(kind=rp) :: t, cfl
61 real(kind=dp) :: start_time_org, start_time, end_time
62 character(len=LOG_SIZE) :: log_buf
64 character(len=:),
allocatable :: restart_file
65 logical :: output_at_end, found
67 real(kind=rp) :: cfl_avrg = 0.0_rp
72 call neko_log%section(
'Starting simulation')
73 write(log_buf,
'(A, E15.7,A,E15.7,A)')
'T : [', 0d0,
',',c%end_time,
')'
75 call dt_controller%init(c%params)
76 if (.not. dt_controller%if_variable_dt)
then
77 write(log_buf,
'(A, E15.7)')
'dt : ', c%dt
80 write(log_buf,
'(A, E15.7)')
'CFL : ', dt_controller%set_cfl
84 call c%params%get(
'case.restart_file', restart_file, found)
85 if (found .and. len_trim(restart_file) .gt. 0)
then
94 call neko_log%section(
'Postprocessing')
95 call c%q%eval(t, c%dt, tstep)
96 call c%s%sample(t, tstep)
98 call c%usr%user_init_modules(t, c%fluid%u, c%fluid%v, c%fluid%w,&
99 c%fluid%p, c%fluid%c_Xh, c%params)
104 cfl = c%fluid%compute_cfl(c%dt)
105 start_time_org = mpi_wtime()
110 start_time = mpi_wtime()
111 if (dt_controller%dt_last_change .eq. 0)
then
114 call dt_controller%set_dt(c%dt, cfl, cfl_avrg, tstep)
116 cfl = c%fluid%compute_cfl(c%dt)
119 write(log_buf,
'(A,I6)')
'Time-step: ', tstep
123 write(log_buf,
'(A,E15.7,1x,A,E15.7)')
'CFL:', cfl,
'dt:', c%dt
128 call c%fluid%step(t, tstep, c%dt, c%ext_bdf, dt_controller)
129 end_time = mpi_wtime()
130 write(log_buf,
'(A,E15.7,A,E15.7)') &
131 'Elapsed time (s):', end_time-start_time_org,
' Step time:', &
136 if (
allocated(c%scalar))
then
137 start_time = mpi_wtime()
139 call c%scalar%step(t, tstep, c%dt, c%ext_bdf, dt_controller)
140 end_time = mpi_wtime()
141 write(log_buf,
'(A,E15.7,A,E15.7)') &
142 'Elapsed time (s):', end_time-start_time_org,
' Step time:', &
147 call neko_log%section(
'Postprocessing')
151 call c%q%eval(t, c%dt, tstep)
152 call c%s%sample(t, tstep)
155 call c%usr%material_properties(t, tstep, c%material_properties%rho,&
156 c%material_properties%mu, &
157 c%material_properties%cp, &
158 c%material_properties%lambda, &
161 call c%usr%user_check(t, tstep,&
162 c%fluid%u, c%fluid%v, c%fluid%w, c%fluid%p, c%fluid%c_Xh, c%params)
172 output_at_end, .true.)
173 call c%s%sample(t, tstep, output_at_end)
175 if (.not. (output_at_end) .and. t .lt. c%end_time)
then
179 call c%usr%user_finalize_modules(t, c%params)
181 call neko_log%end_section(
'Normal end.')
186 real(kind=rp),
intent(inout) :: t
187 real(kind=rp),
intent(in) :: dt
189 real(kind=rp),
dimension(10) :: tlag
190 real(kind=rp),
dimension(10) :: dtlag
191 integer,
intent(in) :: step
197 dtlag(i) = dtlag(i-1)
202 if (ext_bdf%ndiff .eq. 0)
then
209 call ext_bdf%set_coeffs(dtlag)
216 type(
case_t),
intent(inout) :: C
217 real(kind=rp),
intent(inout) :: t
219 type(
file_t) :: chkpf, previous_meshf
220 character(len=LOG_SIZE) :: log_buf
221 character(len=:),
allocatable :: restart_file
222 character(len=:),
allocatable :: restart_mesh_file
226 call c%params%get(
'case.restart_file', restart_file, found)
227 call c%params%get(
'case.restart_mesh_file', restart_mesh_file,&
231 previous_meshf =
file_t(trim(restart_mesh_file))
232 call previous_meshf%read(c%fluid%chkp%previous_mesh)
235 call c%params%get(
'case.mesh2mesh_tolerance', tol,&
238 if (found) c%fluid%chkp%mesh2mesh_tol = tol
240 chkpf =
file_t(trim(restart_file))
241 call chkpf%read(c%fluid%chkp)
242 c%dtlag = c%fluid%chkp%dtlag
243 c%tlag = c%fluid%chkp%tlag
246 call c%fluid%chkp%previous_mesh%free()
247 do i = 1,
size(c%dtlag)
248 call c%ext_bdf%set_coeffs(c%dtlag)
251 call c%fluid%restart(c%dtlag, c%tlag)
252 if (
allocated(c%scalar))
call c%scalar%restart( c%dtlag, c%tlag)
254 t = c%fluid%chkp%restart_time()
255 call neko_log%section(
'Restarting from checkpoint')
256 write(log_buf,
'(A,A)')
'File : ', &
259 write(log_buf,
'(A,E15.7)')
'Time : ', t
263 call c%s%set_counter(t)
268 type(
case_t),
intent(inout) :: C
269 real(kind=rp),
intent(inout) :: t
271 character(len=LOG_SIZE) :: log_buf
273 call c%fluid%chkp%sync_host()
274 chkpf =
file_t(
'joblimit.chkp')
275 call chkpf%write(c%fluid%chkp, t)
276 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.
Device abstraction, common interface for various accelerators.
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.
subroutine, public col2(a, b, n)
Vector multiplication .
subroutine, public profiler_start
Start profiling.
subroutine, public profiler_end_region
End the most recently started profiler region.
subroutine, public profiler_start_region(name, region_id)
Started a named (name) profiler region.
subroutine, public profiler_stop
Stop profiling.
Contains the simcomp_executor_t type.
type(simcomp_executor_t), 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.