Path: blob/main/src/callbacks_step/save_restart.jl
2802 views
# By default, Julia/LLVM does not use fused multiply-add operations (FMAs).1# Since these FMAs can increase the performance of many numerical algorithms,2# we need to opt-in explicitly.3# See https://ranocha.de/blog/Optimizing_EC_Trixi for further details.4@muladd begin5#! format: noindent67"""8SaveRestartCallback(; interval=0,9save_final_restart=true,10output_directory="out")1112Save the current numerical solution in a restart file every `interval` time steps.13"""14struct SaveRestartCallback15interval::Int16save_final_restart::Bool17output_directory::String18end1920function Base.show(io::IO, cb::DiscreteCallback{<:Any, <:SaveRestartCallback})21@nospecialize cb # reduce precompilation time2223restart_callback = cb.affect!24print(io, "SaveRestartCallback(interval=", restart_callback.interval, ")")25return nothing26end2728function Base.show(io::IO, ::MIME"text/plain",29cb::DiscreteCallback{<:Any, <:SaveRestartCallback})30@nospecialize cb # reduce precompilation time3132if get(io, :compact, false)33show(io, cb)34else35save_restart_callback = cb.affect!3637setup = [38"interval" => save_restart_callback.interval,39"save final solution" => save_restart_callback.save_final_restart ? "yes" :40"no",41"output directory" => abspath(normpath(save_restart_callback.output_directory))42]43summary_box(io, "SaveRestartCallback", setup)44end45end4647function SaveRestartCallback(; interval = 0,48save_final_restart = true,49output_directory = "out")50restart_callback = SaveRestartCallback(interval, save_final_restart,51output_directory)5253return DiscreteCallback(restart_callback, restart_callback, # the first one is the condition, the second the affect!54save_positions = (false, false),55initialize = initialize!)56end5758function initialize!(cb::DiscreteCallback{Condition, Affect!}, u, t,59integrator) where {Condition, Affect! <: SaveRestartCallback}60restart_callback = cb.affect!6162mpi_isroot() && mkpath(restart_callback.output_directory)6364semi = integrator.p6566@trixi_timeit timer() "I/O" begin67save_mesh(semi, restart_callback.output_directory)68end6970return nothing71end7273# this method is called to determine whether the callback should be activated74function (restart_callback::SaveRestartCallback)(u, t, integrator)75@unpack interval, save_final_restart = restart_callback7677# With error-based step size control, some steps can be rejected. Thus,78# `integrator.iter >= integrator.stats.naccept`79# (total #steps) (#accepted steps)80# We need to check the number of accepted steps since callbacks are not81# activated after a rejected step.82return interval > 0 && (integrator.stats.naccept % interval == 0 ||83(save_final_restart && isfinished(integrator)))84end8586# this method is called when the callback is activated87function (restart_callback::SaveRestartCallback)(integrator)88u_ode = integrator.u89@unpack t, dt = integrator90iter = integrator.stats.naccept91semi = integrator.p9293@trixi_timeit timer() "I/O" begin94save_mesh(semi, restart_callback.output_directory, iter)95save_restart_file(u_ode, t, dt, iter, semi, restart_callback)96# If using an adaptive time stepping scheme, store controller values for restart97if integrator.opts.adaptive98save_adaptive_time_integrator(integrator, integrator.opts.controller,99restart_callback)100end101end102103# avoid re-evaluating possible FSAL stages104u_modified!(integrator, false)105return nothing106end107108@inline function save_restart_file(u_ode, t, dt, iter,109semi::AbstractSemidiscretization, restart_callback)110mesh, equations, solver, cache = mesh_equations_solver_cache(semi)111u = wrap_array_native(u_ode, mesh, equations, solver, cache)112return save_restart_file(u, t, dt, iter, mesh, equations, solver, cache,113restart_callback)114end115116"""117load_time(restart_file::AbstractString)118119Load the time saved in a `restart_file`.120"""121function load_time(restart_file::AbstractString)122h5open(restart_file, "r") do file123return read(attributes(file)["time"])124end125end126127"""128load_timestep(restart_file::AbstractString)129130Load the time step number (`iter` in OrdinaryDiffEq.jl) saved in a `restart_file`.131"""132function load_timestep(restart_file::AbstractString)133h5open(restart_file, "r") do file134return read(attributes(file)["timestep"])135end136end137138"""139load_timestep!(integrator, restart_file::AbstractString)140141Load the time step number saved in a `restart_file` and assign it to both the time step142number and and the number of accepted steps143(`iter` and `stats.naccept` in OrdinaryDiffEq.jl, respectively) in `integrator`.144"""145function load_timestep!(integrator, restart_file::AbstractString)146integrator.iter = load_timestep(restart_file)147return integrator.stats.naccept = integrator.iter148end149150"""151load_dt(restart_file::AbstractString)152153Load the time step size (`dt` in OrdinaryDiffEq.jl) saved in a `restart_file`.154"""155function load_dt(restart_file::AbstractString)156h5open(restart_file, "r") do file157return read(attributes(file)["dt"])158end159end160161function load_restart_file(semi::AbstractSemidiscretization, restart_file)162return load_restart_file(mesh_equations_solver_cache(semi)..., restart_file)163end164165"""166load_adaptive_time_integrator!(integrator, restart_file::AbstractString)167168Load the context information for time integrators with error-based step size control169saved in a `restart_file`.170"""171function load_adaptive_time_integrator!(integrator, restart_file::AbstractString)172controller = integrator.opts.controller173# Read context information for controller174h5open(restart_file, "r") do file175# Ensure that the necessary information was saved176if !("time_integrator_qold" in keys(attributes(file))) ||177!("time_integrator_dtpropose" in keys(attributes(file))) ||178(hasproperty(controller, :err) &&179!("time_integrator_controller_err" in keys(attributes(file))))180error("Missing data in restart file: check the consistency of adaptive time controller with initial setup!")181end182# Load data that is required both for PIController and PIDController183integrator.qold = read(attributes(file)["time_integrator_qold"])184integrator.dtpropose = read(attributes(file)["time_integrator_dtpropose"])185# Accept step to use dtpropose already in the first step186integrator.accept_step = true187# Reevaluate integrator.fsal_first on the first step188integrator.reeval_fsal = true189# Load additional parameters for PIDController190if hasproperty(controller, :err) # Distinguish PIDController from PIController191controller.err[:] = read(attributes(file)["time_integrator_controller_err"])192end193end194end195196include("save_restart_dg.jl")197end # @muladd198199200