First step with snek5000-phill
: setting up a simulation#
Installation of snek5000-phill
#
The
Nek5000 phill example
(phill
means periodic hill) has been adapted for a workflow using Snek5000. We created
a small Python project called snek5000-phill
using the Snek5000 API to define the phill
case. In this tutorial and in the tutorials
Second step with snek5000-phill: Snakemake rules and Demo periodic hill (snek5000-phill solver): running a simulation from a script, we will show how this solver
can be used. Let’s first recall that the installation procedure looks like this:
export NEK_SOURCE_ROOT="/path/to/Nek5000"
python -m venv venv
source venv/bin/activate
pip install snek5000-phill
Note
The repository snek5000-phill
provides a Python package phill
- the same
name as the usr
files. More on packaging Nek5000 files can be found in the
packaging tutorial.
params
object for simulation parameters#
Creating the params
object containing default parameters#
from phill.solver import Simul
params = Simul.create_default_params()
The params
object gives you a consolidated view of the parameters which are spread out
in a typical Nek5000 case into .par
, .box
and SIZE
file. We will see that the
parameters are more verbose and easier to understand. As a bonus, some Nek5000
parameters which depend on others are automatically set. For example, see
snek5000.operators
.
Studying the parameters (print and get help)#
Now let us take a look at the parameters defining a simulation. In a IPython or Jupyter
console, the params
object would also output as follows:
print(params)
<snek5000.params.Parameters object at 0x7f98743b6070>
<params NEW_DIR_RESULTS="True" short_name_type_run="run">
<nek>
<general dealiasing="True" dt="-0.0002" end_time="nan"
extrapolation="standard" filter_cutoff_ratio="0.67"
filter_modes="2" filter_weight="0.02" filtering="explicit"
log_level="2" num_steps="20" opt_level="2" start_from=""
stop_at="num_steps" target_cfl="0.3" time_stepper="bdf3"
variable_dt="False" write_control="timestep"
write_double_precision="True" write_interval="100"/>
<problemtype equation="incompNS" stress_formulation="False"
variable_properties="False"/>
<velocity advection="True" density="1.0" residual_proj="1e-08"
residual_tol="1e-08" viscosity="-700" write_to_field_file="True"/>
<pressure preconditioner="semg_xxt" residual_proj="1e-08"
residual_tol="1e-08" write_to_field_file="True"/>
<mesh write_to_field_file="True"/>
<temperature absolute_tol="nan" advection="True" conductivity="nan"
conjugate_heat_transfer="False" residual_proj="False"
residual_tol="nan" rho_cp="nan" solver="helm"
write_to_field_file="True"/>
<scalar01 absolute_tol="nan" advection="True" density="nan"
diffusivity="nan" residual_proj="False" residual_tol="nan"
solver="helm" write_to_field_file="True"/>
<cvode/>
<runpar parf_name="outparfile" parf_write="False"/>
<monitor log_level="4" wall_time="23:45"/>
<chkpoint chkp_fnumber="1" chkp_interval="250" read_chkpt="False"/>
<stat av_step="10" io_step="50"/>
</nek>
<oper Lx="1.0" Ly="1.0" Lz="1.0" boundary="['P', 'P', 'W', 'W', 'P', 'P']"
boundary_scalars="[]" dim="3" nproc_max="32" nproc_min="8" nx="22"
ny="16" nz="19" origin_x="0.0" origin_y="0.0" origin_z="0.0"
ratio_x="1.0" ratio_y="1.0" ratio_z="1.0" scalars="1">
<max dim_krylov="30" dim_proj="20" hist="1" obj="1" perturb="1"
scalars_cons="1" scalars_proj="1" sessions="1"/>
<elem coef_dealiasing="0.6666666666666666" order="6" order_out="6"
staggered="True"/>
<misc fast_diag="False"/>
</oper>
<output HAS_TO_SAVE="True" session_id="0" sub_directory="">
<phys_fields available_readers="['pymech', 'pymech_stats']"
reader="pymech"/>
<history_points coords="None" write_interval="100"/>
<remaining_clock_time period_save_in_seconds="5.0"/>
</output>
</params>
One can also print only one section:
print(params.nek.general)
<snek5000.params.Parameters object at 0x7f9848e26cd0>
<general dealiasing="True" dt="-0.0002" end_time="nan" extrapolation="standard"
filter_cutoff_ratio="0.67" filter_modes="2" filter_weight="0.02"
filtering="explicit" log_level="2" num_steps="20" opt_level="2"
start_from="" stop_at="num_steps" target_cfl="0.3" time_stepper="bdf3"
variable_dt="False" write_control="timestep"
write_double_precision="True" write_interval="100"/>
or alternatively, with the method _print_as_code
(useful for copy-pasting):
params.nek.general._print_as_code()
general.dealiasing = True
general.dt = -0.0002
general.end_time = nan
general.extrapolation = "standard"
general.filter_cutoff_ratio = 0.67
general.filter_modes = 2
general.filter_weight = 0.02
general.filtering = "explicit"
general.log_level = 2
general.num_steps = 20
general.opt_level = 2
general.start_from = ""
general.stop_at = "num_steps"
general.target_cfl = 0.3
general.time_stepper = "bdf3"
general.variable_dt = False
general.write_control = "timestep"
general.write_double_precision = True
general.write_interval = 100
Note that it is easy to print some help about some parameters:
params.oper._print_docs()
Documentation for params.oper
-----------------------------
Parameters for mesh description:
- ``nx``, ``ny``, ``nz``: int
Number of elements in each directions
- ``origin_x``, ``origin_y``, ``origin_z``: float
Starting coordinate of the mesh (default: 0.0)
- ``ratio_x``, ``ratio_y``, ``ratio_z``: float
Mesh stretching ratio (default: 1.0)
- ``Lx``, ``Ly``, ``Lz``: float
Length of the domain
Parameters for boundary conditions:
- ``boundary``: list[str]
:ref:`Velocity boundary conditions <nek:tab:ubcf>`
- ``boundary_scalars``: list[str]
:ref:`Temperature and passive scalar boundary conditions <nek:tab:userbct>`
The following table matches counterpart of mandatory :ref:`SIZE
<nek:case_files_size>` variables.
========== =================== =============================================
SIZE params.oper Comment
========== =================== =============================================
``ldim`` ``dim`` Domain dimensions (2 or 3)
``lpmin`` ``nproc_min`` Min MPI ranks
``lpmax`` ``nproc_max`` Max MPI ranks
``ldimt`` ``scalars`` Number of auxilary fields (minimum 1).
========== =================== =============================================
Documentation for params.oper.max
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
The following table matches counterpart of optional :ref:`SIZE
<nek:case_files_size>` variables. These refer to upper bound number of
`something`. The parameters are considered "optional" and would be ignored with
the default values.
============== =================== =========================================
SIZE params.oper.max Comment
============== =================== =========================================
``mxprev`` ``dim_proj`` Max. dimension of projection space
``lgmres`` ``dim_krylov`` Max. dimension of Krylov space for GMRES
``lhis`` ``hist`` Max. number of history (i.e. monitoring)
points.
``maxobj`` ``obj`` Max. number of objects?
``lpert`` ``perturb`` Max. number of perturbations
``toteq`` ``scalars_cons`` Max. conserved scalars
``ldimt_proj`` ``scalars_proj`` Max. scalars for residual projection
``nsessmax`` ``sessions`` Max. sessions to NEKNEK
============== =================== =========================================
Documentation for params.oper.elem
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
The following table relate to element configuration for certain operations.
The parameters are considered "optional" (except for ``lx1`` / ``order`` and
``lxd`` / ``coef_dealiasing`` which are mandatory) and would be ignored with
the default values.
========== =================== =========================================
SIZE params.oper.elem Comment
========== =================== =========================================
``lxd`` ``coef_dealiasing`` | p-order [#f1]_ for over-integration.
**Automatically computed** from float
| ``coef_dealiasing``. See
:any:`order_dealiasing`
``lx1`` ``order`` p-order (avoid uneven and values <6).
``lxo`` ``order_out`` Max. p-order on output (should be ``>=
order``. See :any:`order_out`)
``lx2`` ``staggered`` | p-order for pressure. **Automatically
computed**. See :any:`order_pressure`
========== =================== =========================================
.. [#f1] Polynomial order which means the number of GLL / GL points per element.
Documentation for params.oper.misc
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Other miscellaneous parameters:
========== =================== =========================================
SIZE params.oper.misc Comment
========== =================== =========================================
``lfdm`` ``fast_diag`` | Equals to True for global tensor
product solver (that uses fast
| diagonalization method). ``False``
otherwise.
========== =================== =========================================
Warning
We will see later that the same object params
can also be obtained from a simulation
object (with sim.params
), but the help on the parameters can only be printed from
params
objects obtained from Simul.create_default_params()
.
Modifying the parameters#
Of course, the parameters can be modified. For instance, let us tweak the number of elements, time-stepping and I/O parameters
# This affects both the .box and SIZE files
params.oper.nx = 12
params.oper.ny = 10
params.oper.nz = 8
# This affects the .par file
params.nek.velocity.residual_tol = 1e-8
params.nek.pressure.residual_tol = 1e-6
params.nek.general.num_steps = 10
params.nek.general.time_stepper = "bdf2"
params.nek.general.write_interval = 10
params.nek.stat.av_step = 3
params.nek.stat.io_step = 10
Creation of the simulation object and directory#
Now let’s create the simulation object (We usually use the name sim
.):
sim = Simul(params)
path_run: /home/docs/Sim_data/phill_run_22x16x19_V1.x1.x1._2023-08-24_07-55-59
INFO: session_id: 0
INFO: Writing params files... /home/docs/Sim_data/phill_run_22x16x19_V1.x1.x1._2023-08-24_07-55-59/phill.par, params_simul.xml, info_solver.xml
INFO: sim: <class 'phill.solver.SimulPhill'>
sim.oper: <class 'snek5000.operators.Operators'>
sim.output.print_stdout: <class 'snek5000.output.print_stdout.PrintStdOut'>
sim.output.phys_fields: <class 'snek5000.output.phys_fields.PhysFields'>
sim.output.history_points: <class 'snek5000.output.history_points.HistoryPoints'>
sim.output.remaining_clock_time: <class 'snek5000.output.remaining_clock_time.RemainingClockTime'>
sim.output: <class 'phill.output.OutputPhill'>
sim.make: <class 'snek5000.make.Make'>
INFO: Writing box file... /home/docs/Sim_data/phill_run_22x16x19_V1.x1.x1._2023-08-24_07-55-59/phill.box
INFO: Writing SIZE file... /home/docs/Sim_data/phill_run_22x16x19_V1.x1.x1._2023-08-24_07-55-59/SIZE
Information is given about the structure of the sim
object (which has attributes
sim.oper
, sim.output
and sim.make
corresponding to different aspects of the
simulation) and about the creation of some files.
We see that a directory has been created with all the files (.usr
, .par
, .box
,
SIZE
, etc.) necessary to compile and run the simulation. Note that these files has
been created based on templates files in the snek5000
and phill
packages. The path
towards the simulation directory is:
sim.path_run
PosixPath('/home/docs/Sim_data/phill_run_22x16x19_V1.x1.x1._2023-08-24_07-55-59')
Let’s check that the .par
, .box
and SIZE
files are present:
!ls {sim.path_run}
SIZE etc map_user_params.json phill.log toolbox
Snakefile info_solver.xml params_simul.xml phill.par
config_simul.yml makefile_usr.inc phill.box phill.usr
Note
The path of the directory where the simulation directories are saved can be
changed through the environment variable FLUIDSIM_PATH
.