Neko
0.8.1
A portable framework for high-order spectral element flow simulations
|
The user file is a fortran file where the user can implement their own functions to extend the capabilities of the default Neko executable. The user file can be used for setting advanced initial/boundary conditions, source terms, I/O operations, and interactions with the Neko framework.
The user file is a regular Fortran .f90
file that needs to be compiled with makeneko
, located in the bin
folder of your neko installation. To compile a user file user.f90
, run:
If everything goes well, you should observe the following output:
Compiling your user file with makeneko
will create a neko
executable, which you will need to execute with your case file as an argument. For example, if your case file is called user.case
:
Or in parallel using MPI:
The current high-level structure of the user file is shown below.
The user file implements the user
module. The user
modules contains a subroutine named user_setup
, which we use to interface the internal procedures defined in src/common/user_intf.f90
with the subroutines that you will implement in your user file. Each user subroutine should be implemented under the contains
statement, below user_setup
.
bin/neko
in your local neko installation folder.The following user functions, if defined in the user file, will always be executed, regardless of what is set in the case file:
rho
, mu
, cp
and lambda
.d=1
. For more information on the scalar, see the relevant section of the case file.The two subroutines user_init_modules
and user_finalize_modules
may be used to initialize/finalize any user defined variables, external objects, or processes. They are respectively executed right before/after the simulation time loop.
In the example above, the subroutines initialize
and finalize
contain the actual implementations. They must also be interfaced to the internal procedures user_init_modules
and user_finalize_modules
in user_setup
:
user_init_modules
and user_finalize_modules
are independent of each other. Using one does not require the use of the other.The subroutine user_check
is executed at the end of every time step. It can be used for computing and/or outputting your own variables/quantities at every time step.
In the example above, the subroutine usercheck
contains the actual implementation, and needs to be registered by adding:
to our user_setup
.
material_properties
allows for more complex computations and setting of various material properties, such as rho
, mu
for the fluid and cp
, lambda
for the scalar. The example below is taken from the rayleigh-benard-cylinder example.
And of course not forgetting to register our function in user_setup
by adding the following line:
This user function allows for the modification of the mesh at runtime, by acting on the element nodes of the mesh specified in the case file. This function is only called once before the simulation time loop. The example below is taken from the tgv example.
The registering of the above function in user_setup
should then be done as follows:
This user function can be used to specify the scalar boundary values, on all zones that are not already set to uniform Dirichlet or Neumann values e.g. d=1
or n=0
. For more information on the scalar, see the relevant section of the case file. The example below sets the scalar boundary condition values to be a linear function of the z
coordinate (taken from the rayleigh-benard example).
This function will be called on all the points on the relevant boundaries. The registering of the above function in user_setup
should be done as follows:
As explained in the case file page, certain components of the simulation can be set to be user defined. These components and their associated user functions are:
Description | User function | JSON Object in the case file |
---|---|---|
Fluid initial condition | fluid_user_ic | case.fluid.initial_condition |
Scalar initial condition | scalar_user_ic | case.scalar.initial_condition |
Fluid inflow boundary condition | fluid_user_if | case.fluid.inflow_condition |
Scalar boundary conditions | scalar_user_bc | (user function is always called) |
Fluid source term | fluid_user_f_vector or fluid_user_f | case.fluid.source_terms |
Scalar source term | scalar_user_f_vector or scalar_user_f | case.scalar.source_terms |
Fluid and Scalar boundary conditions | field_dirichlet_update | case.fluid.boundary_types and/or case.scalar.boundary_types |
Note that scalar_user_bc
is included for completeness but is technically not case-specific.
Enabling user defined initial conditions for the fluid and/or scalar is done by setting the initial_condition.type
to "user"
in the relevant sections of the case file, case.fluid
and/or case.scalar
.
See the relevant sections on the fluid and scalar initial conditions in the case file page for more details.
The associated user functions for the fluid and/or scalar initial conditions can then be added to the user file. An example for the fluid taken from the advecting cone example, is shown below.
NEKO_BCKND_DEVICE
flag, which will be set to 1 if running on GPUs, and the calls to device_memcpy
to transfer data between the host and the device. See Running on GPUs for more information on how this works.The same can be done for the scalar, with the example below also inspired from the advecting cone example:
NEKO_BCKND_DEVICE
flag, which will be set to 1 if running on GPUs, and the calls to device_memcpy
to transfer data between the host and the device. See Running on GPUs for more information on how this works.We should also add of the following lines in user_setup
, registering our user functions set_velocity
and set_s_ic
to be used as the fluid and scalar initial conditions:
Enabling user defined inflow condition for the fluid is done by setting the case.fluid.inflow_condition.type
to "user"
:
See the the relevant section in the case file page for more details. The associated user function for the fluid inflow condition can then be added to the user file. An example inspired from the lid-driven cavity example is shown below.
We should also add of the following line in user_setup
, registering our user function user_bc
to be used as the fluid inflow conditions:
Enabling user defined source terms for the fluid and/or scalar is done by adding JSON Objects to the case.fluid.source_terms
and/or case.scalar.source_terms
lists.
See the relevant sections on the fluid and scalar source terms in the case file page for more details.
_user_f
and _user_f_vector
. The former is called when setting "user_pointwise"
as the source term type, while the latter requires the use of the "user_vector"
keyword in the case file. The pointwise variant, fluid_user_f
is not supported on GPUs. In general, fluid_user_f_vector
is the prefered variant, and is the one which will be use in our examples below. The same applies for the scalar source term user functions.The associated user functions for the fluid and/or scalar source terms can then be added to the user file. An example for the fluid, taken from the rayleigh-benard-cylinder example, is shown below.
neko_field_registry
to retrieve the velocity and scalar fields. See Registries for more information about registries in neko. NEKO_BCKND_DEVICE
flag, which will be set to 1 if running on GPUs, and the use of device_
functions. See Running on GPUs for more information on how this works.The same can be done for the scalar, with the example below also taken from the scalar_mms example:
NEKO_BCKND_DEVICE
flag, which will be set to 1 if running on GPUs, and the call to device_memcpy
to transfer data between the host and the device. See Running on GPUs for more information on how this works.We should also add of the following lines in user_setup
, registering our user functions set_boussinesq_forcing_term
and set_source
to be used as the fluid and scalar source terms:
This user function can be used to specify dirichlet boundary values for velocity components u,v,w
, the pressure p
, and/or the scalar s
. This type of boundary condition allows for time-dependent velocity profiles (currently not possible with a standard user_inflow
) or non-uniform pressure profiles to e.g. impose an outlet pressure computed from another simulation.
The selection of such boundary condition is done in the case.fluid.boundary_types
array for the velocities and pressure, and in the case.scalar.boundary_types
array for the scalar. The case file outlines which keywords can be used for such purpose:
d_vel_u
for the u
component of the velocity fieldd_vel_v
for the v
component of the velocity fieldd_vel_w
for the w
component of the velocity fieldd_pres
for the pressure fieldd_s
for the scalar field (cannot be combined with the above)The separator "/"
can be used to combine the keywords related to u,v,w
and p
. For example, if one wants to only apply u,v
and p
values on a given boundary, one should use "d_vel_u/d_vel_v/d_pres"
. In this case, the w
component would be left untouched (not zeroed!). An example of case file from the cyl-boundary-layer example is shown below.
In this example, we indicate in case.fluid.boundary_types
that we would like to apply a velocity profile on all three components u,v,w
on the boundary number 1 (in this case, the inlet boundary). On boundary number 2 (the outlet boundary), we also indicate the three velocity components, with the addition of the pressure. In case.scalar.boundary_types
, we indicate the same for the scalar on boundaries 1 and 2 (inlet and outlet).
d_s
and d=x
boundary conditions for the scalar. The latter is to be used to specify a constant Dirichlet value x
along the relevant boundary.Once the appropriate boundaries have been identified and labeled, the user function field_dirichlet_update
should be used to compute and apply the desired values to our velocity/pressure/scalar field(s). The prefix "field" in field_dirichlet_update
refers to the fact that a list of entire fields is passed down for the user to edit.
The fields that are passed down are tied to the boundary_types
keywords passed in the case file. The function field_dirichlet_update
is then called internally, one time in the fluid
solver and one time in the scalar
solver (if enabled).
Finally, depending on which boundary labels were input, the fields given to the user are copied onto the solution field boundaries.
The header of the user function is given in the code snippet below.
The arguments and their purpose are as follows:
field_bc_list
is the list of the field that can be edited. It is a list of field_t
objects.i
contained in field_bc_list
is accessed using field_bc_listfields(i)f
and will refer to a field_t
object.which_solver = "fluid"
, it will contain the 4 fields u,v,w,p
. They are retrieved in that order in field_bc_list
, i.e. u
corresponds to field_bc_listfields(1)f
, etc.which_solver = "scalar"
, it will only contain the scalar field s
.bc_bc_list
contains a list of the bc_t
objects to help access the boundary indices through the boundary mask
.i
contained in bc_bc_list
is accessed with bc_bc_listbc(i)bcp
.i
-th bc_t
object contained in bc_bc_list
is accessed with bc_bc_listbc(i)bcpmsk
. It contains the linear indices of each GLL point on the i
-th boundary facets. msk(0)
contains the size of the array. The first boundary index is msk(1)
.which_solver = "fluid"
, it will contain the 4 bc_t
objects corresponding to d_vel_u
, d_vel_v
, d_vel_w
, and d_pres
. They can be retrieved in that order, in the same way as for field_bc_list
.which_solver = "scalar"
, it will only the 1 bc_t
object corresponding to d_s
.coef
is a coef_t
object containing various numerical parameters and variables, such as the polynomial order lx
, derivatives, facet normals...t
, tstep
are self-explanatory.which_solver
takes the value "fluid"
when the user function is called in the fluid solver. It takes the value "scalar"
when it is called in the scalar solver.Links to the documentation to learn more about what the types mentioned above contain and how to use them: field::field_t, bc::bc_t, coefs::coef_t.
The user function should be registered in user_setup
with the following line:
A very simple example illustrating the above is shown below, which is taken from the cyl_boundary_layer example
This example is applying constant dirichlet values at the selected boundaries for the velocity components and presure. The scalar is applied a function s(y,z) = sin(y)*sin(z)
to demonstrate the usage of boundary masks.
u = 1.0_rp
is only possible because of the overloading of the assignement operator =
in field_t
. In general, a field's array should be accessed and modified with u%x
.Note that we are only applying our boundary values at the first timestep, which is done simply with the line if (tstep .ne. 1) return
. This is a trick that can be used for time independent boundary profiles that require some kind of time consuming operation like interpolation or reading from a file, which would add overhead if executed at every time step.
Observe that we always check if the fields are allocated before manipulating them. This is to prevent accidental memory access if only part of the velocity components or pressure are given in case.fluid.boundary_types
. Fields in the lists are only allocated if they are present in the case file.For example, if we removed the d_pres
condition in the JSON case file code snippet above, the pressure field for our boundary condition would not be allocated ( in the example above, allocated(px)
would never be true
). "boundary_types": ["d_vel_u", "d_vel_v"]
will allocate the two first fields in field_bc_list
, which is the same behaviour as "boundary_types": ["d_vel_u/d_vel_v", ""]
.
device_memcpy
to make sure the device arrays are also updated.When running on GPUs, special care must be taken when using certain user functions. The short explanation is that the device (GPU) has its own memory and cannot directly access the memory on the host (CPU). This means that data and more specifically arrays must be copied manually from the host to the device (see device::device_memcpy).
device_memcpy
is avoidable. Neko has some device math functions implemented that operate directly on device arrays. If you can decompose whatever operations you are performing in a user function into a set of instructions from the math
module (e.g. cadd
, cfill
, sub2
, ...), you may use the corresponding device_math
functions to offload work to the GPU. See the fluid forcing code snippet for a simple example. For more advanced examples, see the rayleigh-benard example or the tgv example.To illustrate this, let us have a look at the fluid initial condition code snippet:
The code above is used to set the fluid initial condition, by specifying the values of fields u,v,w
(and p
) at all points in the domain. Notice that we have divided the above code into two parts.
In the first part, we set the velocity components u=-y*pi*
, v=x*pi*
, and w=0
, which updates the velocity field arrays u%x, v%x, w%x
allocated on the host (CPU). If we were to run on GPUs, these lines of code would only act on the velocity arrays on the host (CPU), leaving the device (GPU) arrays untouched.
We take care of this in the second part, for all three velocity arrays. To update the device (GPU) arrays, we use device_memcpy
to copy the data contained in a host (CPU) array to a device (GPU) array. Looking at the details of the device_memcpy
calls, we note the following:
_d
to the host array variable name (e.g. u%x
and u%x_d
). This is the standard in Neko.HOST_TO_DEVICE
. Other flags can also be used to move data from device to host (DEVICE_TO_HOST
) or device to device (DEVICE_TO_DEVICE
). See the accelerators page for more details on this.sync
argument is a non-optional argument which dictates wether or not to perform the data transfer synchronously.sync = .true.
as a starting point.Finally, observe that we use the flag NEKO_BCKND_DEVICE
to check if we are indeed running on GPUs. In that case, NEKO_BCKND_DEVICE
would be equal to 1.
Neko uses the concept of registry
as a practical way to retrieve fields and point zones anywhere in the user file.
The field registry neko_field_registry
is often used in user functions where certain fields are not directly accessible as arguments. One can retrieve any field in the registry by its name
with neko_field_registry%get_field(name)
. Default fields that are added to the registry are u,v,w,p
and s
if running with the scalar enabled. For a practical example of usage, see the rayleigh benard example
Other fields may be added to the registry by various simulation components. For example:
simulation_components.vorticity
enabled, the fields omega_x, omega_y, omega_z
will be accessible in the registry.simulation_components.lambda2
enabled, the field lambda2
will be accessible in the registry.neko_field_registry%add_field
(see field_registry::add_field).The point zone registry, neko_point_zone_registry
, can be used to retrieve pointers to point_zone_t
objects defined in the case file. See using point zones for detailed instructions.