Skip to content
Snippets Groups Projects
Commit b2cb0265 authored by Pierre Augier's avatar Pierre Augier
Browse files

Implement load_for_restart + add an example

parent c6bb9b9b
No related branches found
No related tags found
No related merge requests found
......@@ -280,7 +280,7 @@
autosummary_generate = True
autodoc_default_flags = ['show-inheritance']
autodoc_default_options = {"show-inheritance": None}
autodoc_member_order = 'bysource'
todo_include_todos = True
......@@ -6,6 +6,7 @@
examples/running-simul-onlineplot
examples/running-simul-cluster
examples/running-simul-restart
examples/init-fields-in-script
examples/init-fields-in-script-3d
examples/forcing-in-script
......
Restart a simulation
====================
Let's run a small simulation with this script:
.. literalinclude:: simul_ns2d.py
To restart this simulation and continue the time_stepping:
.. literalinclude:: simul_ns2d_restart.py
from fluidsim import load_for_restart, path_dir_results
path_dir_root = path_dir_results / "examples"
path_dir = sorted(path_dir_root.glob("NS2D_32x32_S2pix2pi*"))[-1]
params, Simul = load_for_restart(path_dir)
params.time_stepping.t_end += 2
sim = Simul(params)
sim.time_stepping.start()
......@@ -23,5 +23,7 @@
"""
from pathlib import Path
from ._version import __version__
......@@ -26,6 +28,9 @@
from ._version import __version__
from fluiddyn.io import FLUIDSIM_PATH as path_dir_results
from fluiddyn.io import FLUIDSIM_PATH
# has to be done before importing util
path_dir_results = Path(FLUIDSIM_PATH)
from .util.util import (
import_module_solver_from_key,
......@@ -34,6 +39,7 @@
load_state_phys_file,
modif_resolution_from_dir,
modif_resolution_all_dir,
load_for_restart,
)
from .base.params import load_params_simul
......@@ -54,4 +60,5 @@
"modif_resolution_from_dir",
"modif_resolution_all_dir",
"load_params_simul",
"load_for_restart",
]
......@@ -7,10 +7,7 @@
"""
from __future__ import division, print_function
from builtins import object
import os as _os
import glob as _glob
from copy import deepcopy as _deepcopy
import inspect
......@@ -13,7 +10,8 @@
import os as _os
import glob as _glob
from copy import deepcopy as _deepcopy
import inspect
from pathlib import Path
import numpy as _np
import h5py as _h5py
......@@ -151,7 +149,7 @@
return _os.getcwd()
if name_dir[0] != "/" and name_dir[0] != "~":
name_dir = path_dir_results + "/" + name_dir
name_dir = str(path_dir_results / name_dir)
return _os.path.expanduser(name_dir)
......@@ -278,4 +276,39 @@
"""
params, Simul = load_for_restart(name_dir, t_approx, merge_missing_params)
if modif_save_params:
params.output.HAS_TO_SAVE = False
params.output.ONLINE_PLOT_OK = False
sim = Simul(params)
return sim
def load_for_restart(name_dir=None, t_approx=None, merge_missing_params=False):
"""Load params and Simul for a restart.
Parameters
----------
name_dir : str (optional)
Name of the directory of the simulation. If nothing is given, we load the
data in the current directory.
t_approx : number (optional)
Approximate time of the file to be loaded.
merge_missing_params : bool (optional, default == False)
Can be used to load old simulations carried out with an old fluidsim
version.
"""
if isinstance(name_dir, Path):
name_dir = str(name_dir)
path_dir = pathdir_from_namedir(name_dir)
......@@ -281,8 +314,7 @@
path_dir = pathdir_from_namedir(name_dir)
solver = _import_solver_from_path(path_dir)
# choose the file with the time closer to t_approx
name_file = name_file_from_time_approx(path_dir, t_approx)
path_file = _os.path.join(path_dir, name_file)
......@@ -283,9 +315,12 @@
solver = _import_solver_from_path(path_dir)
# choose the file with the time closer to t_approx
name_file = name_file_from_time_approx(path_dir, t_approx)
path_file = _os.path.join(path_dir, name_file)
with _h5py.File(path_file, "r") as f:
params = Parameters(hdf5_object=f["info_simul"]["params"])
if mpi.rank > 0:
params = None
else:
with _h5py.File(path_file, "r") as f:
params = Parameters(hdf5_object=f["info_simul"]["params"])
......@@ -291,4 +326,4 @@
if merge_missing_params:
merge_params(params, solver.Simul.create_default_params())
if merge_missing_params:
merge_params(params, solver.Simul.create_default_params())
......@@ -294,14 +329,11 @@
params.path_run = path_dir
params.NEW_DIR_RESULTS = False
if modif_save_params:
params.output.HAS_TO_SAVE = False
params.output.ONLINE_PLOT_OK = False
params.init_fields.type = "from_file"
params.init_fields.from_file.path = path_file
params.init_fields.modif_after_init = False
try:
params.preprocess.enable = False
except AttributeError:
pass
params.path_run = path_dir
params.NEW_DIR_RESULTS = False
params.init_fields.type = "from_file"
params.init_fields.from_file.path = path_file
params.init_fields.modif_after_init = False
try:
params.preprocess.enable = False
except AttributeError:
pass
......@@ -307,3 +339,3 @@
fix_old_params(params)
fix_old_params(params)
......@@ -309,6 +341,8 @@
sim = solver.Simul(params)
return sim
if mpi.nb_proc > 1:
params = mpi.comm.bcast(params, root=0)
return params, solver.Simul
def modif_resolution_all_dir(t_approx=None, coef_modif_resol=2, dir_base=None):
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment