- Generated on Mon Sep 1 2025 03:49:15 for Neko by
1.9.8
|
Neko 1.99.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.
This section will provide a written explanation of the user file, its structure, and how to compile and run it. For a more hands-on approach, see the examples section, which contains several user file examples that demonstrate various aspects of programming with Neko.
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. This 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:
compute but executed at the beginning of the time step.rho, mu, cp and lambda.Three subroutines startup, initialize and finalize may be used to initialize/finalize any user defined variables, external objects, or processes. The startup routine is called immediately after the case file is read in, meaning that no important objects are initialized yet (e.g. the mesh, fluid, etc.). The initialize and finalize are respectively executed right before/after the simulation time loop.
In most cases, one can use initialize routine and not the startup. The latter is only necessary when you want to either manipulate the case file programmatically before it is used in the simulation, or define some variables that will be used in the constructors of some of the types. An example of the latter is defining some material property constants that can be used in the user material_properties routine, which is run by the constructor of the casefluid object.
In the example above, the subroutines startup, initialize and finalize contain the actual implementations. They must also be associated to the internal procedures inside, inside user_setup:
The subroutine compute 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 user_check contains the actual implementation, as a homage to the Nek5000 routine, which served a similar purpose. As usual, the routine needs to be registered by adding:
to our user_setup. The preprocess routine can be implemented similarly.
The material_properties routine 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:
Note the usage of the scheme_name argument. This is used in several user routines, to provide the information about what solver (scheme) called the routine. In this example, we make use of it to distinguish the fluid and scalar solvers. The default scheme names are fluid and scalar, but can be controlled by the name entry in the JSON configuration of the solver. Here, the name for the scalar was set to temperature.
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:
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 |
|---|---|---|
| Initial conditions | initial_conditions | case.fluid.initial_condition or case.scalar.initial_condition |
| Source terms | source_term | case.fluid.source_terms or case.scalar.source_terms |
| Dirichlet boundary conditions | dirichlet_conditions | case.fluid.boundary_types or case.scalar.boundary_types |
case.scalar in the JSON. For multiple scalars, the same things apply, but the configuration is inside each individual element of case.scalars. In the case of multiple scalars, users should use "checkpoint_format": "hdf5" when enabling checkpointing.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 function for initial conditions must then be added to the user file. An example for the fluid taken from the advecting cone example, is shown below.
Not again the usage of scheme_name to distinguish between the fluid and the scalar. Depending on the scheme, the contents of the field list fields changes, and we extract individual fields via field pointers accordingly. The incompressible fluid solver always generates solution fields, u, v and w. For the scalar, the name of the field coincides with scheme_name. For single scalar cases, this defaults to s. For multiple scalar cases, the field name is set to the scalar name specified in the JSON configuration (e.g., "s1", "s2", etc.).
We should also add of the following lines in user_setup, as usual.
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](Case File) for more details.
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 scalar_mms example, is shown below.
NEKO_BCKND_DEVICE flag, which will be set to 1 if running on GPUs, and the use of device_ functions. See [Running on GPUs](Running on GPUs) for more information on how this works.As usual, we point the user routine to our implementation in user_setup.
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 or non-uniform pressure profiles to e.g. impose an outlet pressure computed from another simulation.
The user routine is called by the user_velocity and user_pressure boundary conditions for the fluid, and the user boundary condtition for the scalar. Once the appropriate boundaries have been identified, the user function dirichlet_conditions should be used to compute and apply the desired values to our velocity/pressure/scalar field(s).
The structure of the interface is very similar to e.g. the initial conditions. One gets a list of solution fields, the contents of which depends on what scheme owns the boundary condition. For user_velocity, it will be a list of 3 fields with names u, v, w; For user_pressure, a list with a single field p; For the scalar, also a single field with the same name as the solution field for the scalar (s by default).
It is crucial to understand that all three boundary condition will call the same routine! So, if one has, for example, both user_velocity for the fluid and user for the scalar, it is necessary to have an if statement in the user routine to distinguish between the two cases. The convenient way to do that is to check the size of the passed field list and the names of the fields inside. For example, if there is one field and it is called s, one executes the code setting the boundary values for the scalar s.
Note that the fields that one manipulates in the user routine are not the actual velocity fields, etc. Instead, the code does a masked copy from the dummy fields manipulated in the user routine to the actual fields of unknowns, with the mask corresponding to the boundary faces. So, even if you somehow manipulate the fields elsewhere in the domain inside the user routine, that will not affect the solution.
In the following example, we indicate in case.fluid.boundary_conditions that we would like to apply a velocity profile on the boundary number 1 (in this case, the inlet boundary). On boundary number 2 (the outlet boundary), we change the pressure. In case.scalar.boundary_conditions, we indicate the same for the scalar on boundaries 1 and 2 (inlet and outlet).
The header of the user function is given in the code snippet below.
The arguments and their purpose are as follows:
fields is the list of the fields that can be edited. It is a list of field_t objects.i contained in fields is accessed using fields%items(i)%ptr and will refer to a field_t object. Alternatively, one can use the get method to retrieve a field by name or index, as done in the examples above for other routines.bc contains a field_dirichlet_t object to help access the boundary indices through the boundary mask, msk.bcobject is accessed via bc%msk. It contains the linear indices of each GLL point on the boundary facets. msk(0) contains the size of the array. The first boundary index is msk(1).time, is a simple structure that contains various info on time stepping, notably, the current time iteration and time value.Links to the documentation to learn more about what the types mentioned above contain and how to use them: src/field/field.f90, src/bc/bc.f90.
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 (timetstep .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.
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. Additionally, when you operate on field_t objects (which is very common), you can use routines in the field_math module, that will automatically run on the right backend. This is highly recommended and avoids if statements.An example where manual transfer to the GPU was necessary has already been given for the source term. It is reiterated below.
Note that f%x is the field's data on the CPU. Therefore, to populate the correspoding array on GPU, we need to call to device_memcpy. Looking at the detailof the device_memcpy call, we note the following:
_d to the host array variable name (e.g. f%x and f%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.
When running on GPUs it is possible to call an own custom kernel. This could be more performant for more complex source terms or boundary conditions.
Assume we have a custom kernel called mydevice_kernel in a CUDA file called mykernel.cu, to set the initial velicity field. To call the custom kernel, the user file must define a C interface to the routine (inside the .cu file) that launches the kernel
Furthermore, the CUDA/HIP file must allow for C linkage, hence the routine mydevice_kernel must be inside an extern "C" block.
The user defined routine for initial conditions calls the kernel in the following way.
Finally, compile using makeneko and provide both the user.f90 and mykernel.cu file.
makeneko does currently only support custom kernels written in CUDA or HIP.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::field_add).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.
Neko proives a special variable called neko_user_access that can be used in the user file to access the internal state of the entire case and the underlying solvers.
In particular the case object (case_t) is reached via neko_user_access%case. There-in the fluid and scalars components refer to the asscociated schemes. You are encouraged to look a bit at case.f90 to see the overall structure.
A very common use case is to get access to various SEM-related objects, as well as the mesh object.
neko_user_access%case%fluid%msh – the mesh.neko_user_access%case%fluid%dm_Xh – the function space used in the SEM.neko_user_access%case%fluid%dm_Xh – the map of degrees of freedom, contains the locations of the GLL nodes.neko_user_access%case%fluid%c_Xh – the coefficients of the SEM,neko_user_access%case%fluid%gs_Xh – the gather-scatter kernels used for direct stiffness summation. 1.9.8