Skip to content

Understanding the Parameter file (.prm)

Luckabarbeau edited this page Nov 24, 2021 · 48 revisions

Launching of an application requires an executable of the required solver, and a parameters file (with extension .prm if it uses the parameter file format, or .json if it uses the JSON file format). The parameter file controls all aspect of Lethe and drives the simulation. This section aims at describing the various parameters available within Lethe.

The parameter file is composed of different subsections that have meaningful titles. In the following, the principal subsections of a parameter template file are explained. The complete options for the parameter file can be obtained by doing the following:

  • Go to /lethe/build/applications/navier_stokes_parameter_template/
  • Run the navier_stokes_parameter_template using the command ./navier_stokes_parameter_template. Two parameter files are obtained for 2D and 3D respectively. These files contain all options and information regarding the parameters that can be used to run a simulation.

In the parameter file format, the parameters are established one by one using the following syntax: 'set time end = 10' would fix the end time to 10. The arguments can be either doubles, integers, strings or a choice between a predefined number of variable. In the parameter files, comments are preceded by the sharp symbol (e.g. '#comment').

  1. Simulation Control
  2. Physical Properties
  3. Multiphysics
  4. Initial Condition
  5. Mesh
  6. Nitsche
  7. FEM
  8. Boundary conditions - CFD
  9. Boundary conditions - Multiphysics
  10. Source term
  11. Analytical Solution
  12. Timer
  13. Mesh Adaptation Control
  14. Restart
  15. Non-Linear Solver Control
  16. Linear Solver Control
  17. Post-Processing
  18. Force and Torques calculation
  19. Dynamic flow control
  20. Void Fraction
  21. CFD-DEM

Simulation Control

This subsection contains the general information of the simulation, including the time integration method and the general output names. It is the most commonly modified section for a simulation.

subsection simulation control
  # Adaptative time-stepping 
  set adapt                        = false

  # Adaptative time step scaling
  set adaptative time step scaling = 1.1

  # The kind of method use to startup high order bdf methods Choices are
  # .
  set bdf startup method           = multiple step bdf

  # Maximal number of vtu output files
  set group files                  = 1

  # log frequency
  set log frequency                = 1

  # Display precision when writing to log
  set log precision                = 6

  # Maximum CFL value
  set max cfl                      = 1

  # The kind of solver for the linear system. Choices are
  # .
  set method                       = steady

  # Number of mesh adaptation (for steady simulations)
  set number mesh adapt            = 0

  # Output the boundaries of the domain along with their ID
  set output boundaries            = false

  # The control for the output of the simulation resultsResults can be either
  # outputted at constant iteration frequency or at constant time
  set output control               = iteration

  # Output frequency
  set output frequency             = 1

  # File output prefix
  set output name                  = out

  # File output prefix
  set output path                  = ./

  # Output time
  set output time                  = 1

  # Scaling factor used in the iterations necessary to start-up the BDF
  # schemes.
  set startup time scaling         = 0.4

  # Tolerance at which the simulation is stopped
  set stop tolerance               = 1e-10

  # Subdivision of mesh cell in postprocessing
  set subdivision                  = 1

  # Time step value
  set time end                     = 1

  # Time step value
  set time step                    = 1.
end
  • The group files parameter controls the number of vtu files generated in a parallel simulation. This parameters allows to reduce the number of files generated when the simulation is ran with a large number of processors. Setting group files to 1 ensures that there is a single vtu file generated. In this case, the file is written using MPI IO functionalities. The value for this parameter should always be a compromise between keeping a low number of files but preventing excessive MPI communications. We have found that the default value of 1 does not have a significant impact on performance on Compute Canada clusters.

  • The method parameter controls the time-stepping method use. The available options are:

    • steady (steady-state simulation)
    • steady-bdf (steady-state simulation using adjoint time stepping with a bdf1 scheme)
    • bdf1 (1st order backward differentiation)
    • bdf2 (2nd order backward differentiation)
    • bdf3 (3rd order backward differentiation)
    • sdirk2 (2nd order singly diagonally implicit Runge Kutta)
    • sdirk3 (3rd order singly diagonally implicit Runge Kutta)
  • The bdf startup method parameter controls the scheme use to start high order bdf scheme.

  • The output frequency controls after which number of iterations the .pvd/.vtu results are written (forces/torques values are written at every time step).

  • The output name and output path control the name of the output files and their directory.

  • The output boundaries control if the boundaries of the domain are written to a file. This will write additional VTU files made of the contour of the domain. This is particularly useful for the visualisation of 3D flows with obstacles or objects.

  • The subdivision parameter enables sub-division of the mesh cells to enable visualisation of high-order elements with Paraview.

  • The time end parameter indicates at which value of the time to end the simulations if they are transient.

  • The time step parameter controls the value of the time step.

  • The adapt parameter controls if adaptive time-stepping is enabled. If adaptive time-stepping is enabled, the time-step will evolve to ensure that the 'max cfl' value is reached.

  • The max cfl parameter controls the maximal value of the CFL condition that can be reached during the simulation. This parameter is only used when adapt is set to true.

  • The adaptative time step scaling parameter controls the rate of increase of the time step value. The new time step value is fixed by adaptative time step scaling * previous value of the time step

  • The stop tolerance parameter controls the tolerance at which the adjoint time stepping steady state simulations (steady_bdf) stops. The adjoint time stepping will stop when the L2 norm of the initial residual is lower than stop tolerance at the start of a non-linear solution step.

  • The log frequency parameter controls the frequency at which information is written to the log (the terminal).

  • The log precision parameter controls the number of significant digit used when writing to the log (the terminal).

Physical Properties

In this subsection, the physical properties of the materials present in the simulation (for instance, viscosity), are defined. For the CFD solver, the only physical property defined is the kinematic viscosity (ref).

Monophasic simulations

NB: Default behavior.

subsection physical properties
  # Kinematic viscosity
  set kinematic viscosity = 1

  # Thermal conductivity
  set thermal conductivity = 1

  # Specific Heat
  set specific heat = 1

  # Thermal expansion coefficient
  set thermal expansion = 1

  # Density
  set density = 1

  # Tracer diffusivity
  set tracer diffusivity = 0
end
  • The kinematic viscosity parameter is the kinematic viscosity of the fluid in units of Length^{2} * Time^{1}. In SI this is m^{2} * s^{-1}

  • The density parameter is the density of the fluid in units of Mass * Length^{-3}

  • The specific heat parameter is the specific heat of the fluid in units of Energy * Temperature^{-1} * Mass^{-1}.

  • The thermal expansion parameter is the thermal expansion coefficient of the fluid with dimension of Temperature^{-1}.

  • The tracer diffusivity parameter is the diffusivity coefficient of the tracer in units of Length^{2} * Time^{1}. In SI this is m^{2} * s^{-1}

Multiphasic simulations

NB: Requires to define set free surface = true in Multiphysics subsection.

For multiphasic simulations, the properties are defined for each fluid. Default values are:

subsection physical properties
set number of fluids		= 2
    subsection fluid 0
       set density 		= 1
       set kinematic viscosity 	= 1
       set specific heat 	= 1
       set thermal conductivity = 1
       set tracer diffusivity   = 0
    end
    subsection fluid 1
       set density 		= 1
       set kinematic viscosity 	= 1
       set specific heat 	= 1
       set thermal conductivity = 1
       set tracer diffusivity   = 0
    end
end
  • number of fluids = 2 is required for a free surface simulations, otherwise an error will be thrown in the terminal.
  • subsection fluid 0 indicates the properties of fluid where the phase indicator = 0 (Volume of Fluid method), as defined when initializing the free surface (see Initial Condition), and correspondingly fluid 1 is located where phase indicator = 1.

Multiphysics

This subsection defines the multiphysics interface of Lethe and enables the solution of Auxiliary Physics in addition to traditional fluid dynamics simulations.

subsection multiphysics
  set heat transfer = false
  set buoyancy force = false
  set tracer = false
  set free surface = false
end
  • The heat transfer parameter adds the heat transfer auxiliary physics. This is an advection-diffusion equation with viscous dissipation.

  • The buoyancy force parameter adds the buoyancy force auxiliary physics. The buoyancy force is calculated using the Boussinesq approximation. Note that heat transfer is a prerequisite for buoyancy force, and it must be set to true when we set buoyancy force to true.

  • The tracer parameter adds a passive tracer auxiliary physics. This is an advection-diffusion equation.

  • The free-surface parameter enable defining a multiphasic simulation, with 2 fluids separated by a free surface. Volume-of-Fluid method is used, with the phase parameter following an advection equation. See Initial Condition for the free surface initial definition and Physical Properties for the property definition of both fluids.

    • Current scope (03/05/2021): 2D simulations, with 2 fluids.

Initial Condition

It is often necessary to set-up complex initial condition when simulating transient problems. Furthermore, setting a coherent initial condition can greatly help the convergence of steady-state problem. Lethe offers three different strategies to set initial conditions : Nodal (through a function), L2 projection (through a function) and viscous. The viscous strategy solves the same CFD problem but with a user-defined viscosity. This can allow the user to initialize the velocity field with a Stokes-like flow for which convergence can be obtained more easily.

subsection initial conditions
 # Type of initial conditions. Choices are L2projection, viscous or nodal
  set type      = nodal

  # viscosity for viscous initial conditions
  set viscosity = 1

  subsection uvwp
    set Function expression = 0; 0; 0; 0 
  end

  subsection free surface
    set Function expression = if (x<0.5 & y<1, 1, 0)
  end
end
  • The type parameter indicates which strategy is used to impose the initial conditions. The choices are : L2projection, viscous or nodal

  • The viscosity parameter controls the viscosity that is used to solver the Navier-Stokes equations when the viscous initial condition is chosen.

  • The subsection uvwp allows the use to select a function (velocity-pressure) to set a nodal or L2 projection initial condition.

  • The subsection free surface defines the areas where both fluids lay at the initial state (see section Multiphysics).

    • Function expression is used with a condition: if (condition, value if true, value if false)
    • in the above example, if (x<0.5 & y<1, 1, 0) means that fluid 1 fills initially the surface where x<0.5 and y<1, the rest being filled with fluid 0.
    • equation in the Cartesian coordinate system can be used for more complex free surface initial geometry. Example for a circle with its center being the center of the mesh: if ( (x^2+y^2)<=(r)^2 ,1,0)

Mesh

This subsection provides information of the simulation geometry and meshing. The simulation geometry shape and size, number of refinements and other required information can be provided here. It should be mentioned that meshes from gmsh can also be defined in this section by setting the type = gmsh.

subsection mesh
  # GMSH file name
  set file name          = none
  set grid arguments     = -1 : 1 : false
  set grid type          = hyper_cube

  # Initial refinement of the mesh
  set initial refinement = 0

  # Type of mesh Choices are gmsh or dealii
  set type               = dealii

  # Indicates that the mesh is a simplex mesh
  set simplex     = false
end
  • The following choices for the mesh type are available:

    • gmsh. If the gmsh type is chosen, a .msh file generated from GMSH can be used. In this case, the grid file name must be specified in the file name variable.

    • dealii. If the dealii type is chosen, dealii grid generator files can be used. For additional documentation on these grids, you can consult the deal.II documentation. In this case, the arguments for the grid generation are put in the grid arguments parameter.

  • The initial refinement number determines the number of refinement that the grid will undergo in the simulation before the simulation is ran. This allows one to automatically refine a coarse grid.

  • simplex. If simplex is set to true, it indicates that the mesh being read is made of only simplex elements. This is an experimental feature of deal.II.

Nitsche

Lethe can also run a simulation using the Nitsche immersed boundary method.

subsection nitsche
    set beta = 10
    set verbosity = verbose
    set calculate torques on solid = true
    set number of solids = 1

    subsection nitsche solid 0
	subsection mesh
		set type = dealii
		set grid type = hyper_ball
		set grid arguments = 0, 0 : 0.25 : true
		set initial refinement = 5
		set simplex = false
	end
	subsection solid velocity
		set Function expression = -y ; x
	end
    end
end
  • beta is a parameter needed to apply the conditions of the immersed mesh on the fluid solution, its value is normally between 1 and 1000. beta = 0 (the solid has no influence on the flow) can be used for debugging purposes. A classical value is beta = 10.
  • subsection nitsche solid 0 define a solid object on which Nitsche immersed boundary is applied. Multiple solids can be added in the same fashion (add subsection nitsche solid 1 etc.)
  • subsection mesh defines the solid mesh used to apply Nitsche immersed boundary. In the mesh subsection, see Mesh for more details.
    • NB: if set type = gmsh and a simplex mesh is given, pass argument set simplex = true (false by default)
  • subsection solid velocity defines the solid mesh velocity with a Function expression.
  • enable particles motion must be set to true if the shape moves. If the mesh movement is static, for example a rotating cylinder, the shape does not have to move within the fluid and this option can be set to false.
  • calculate forces on solid enables forces calculation, written in the file of name solid force name.
  • calculate torques on solid enables torques calculation, written in the file of name solid torque name.
  • particles sub iterations splits the particle time-stepping into n sub time steps. This enables the particles to displace less per iteration and ensures optimal working of the sort_particles_into_cells_and_subdomain() which is used to locate the cells in which the particles reside.
  • verbosity enables printing Nitsche intermediary results to the terminal
  • subsection test enables particle position printing in the terminal.

FEM

In this subsection, the general information of the finite element method utilized in the simulation are specified. The interpolation orders of velocity and pressure, as well as other options for using high order mapping all over the simulation domain, are defined in this section.

subsection FEM
  # interpolation order pressure
  set pressure order    = 1

  # Apply high order mapping everywhere
  set qmapping all      = false

  # interpolation order velocity
  set velocity order    = 1

  # interpolation order temperature
  set temperature order    = 1

  # interpolation order tracer
  set tracer order    = 1
end
  • The qmapping all option enables the use of isoparametric elements everywhere in the domain. By default these are only used on the elements which have a boundary face

Boundary conditions - CFD

This subsection's purpose is defining the boundary conditions associated to the problem.

subsection boundary conditions
  set number                  = 2
  set time dependent          = false
  subsection bc 0
        set type              = function
       
        subsection u
            set Function expression = -y
        end
        subsection v
            set Function expression = x
        end
        subsection w
            set Function expression = 0
        end

        # Center of rotation used for torque calculation
        subsection center of rotation
           set x = 0
           set y = 0
           set z = 0
        end
  end
  subsection bc 1
        set type              = noslip
  end
end
  • number: This is the number of boundary conditions of the problem. Periodicity between 2 boundaries counts as 1 condition.
  • time dependent: This is a bool to define if the boundary condition is time dependent.
  • type is the type of the boundary condition. At the moment, only slip, no-slip, periodic and function are enabled. The function boundary condition enables one to set the individual value of u,v,w using a user-defined function. At the moment, these function can only depend on position (x,y,z) and not on time.

Boundary conditions - Multiphysics

This subsection's purpose is defining the boundary conditions associated to multiphysic problems.

Heat Transfer

subsection boundary conditions heat transfer
  set number                  = 2
    subsection bc 0
	set type	      = temperature
        set value	      = 0
    end
    subsection bc 1
      set type		    = convection
      set h 		    = 10
      set Tinf	   	    = -13
    end
end
  • number: This is the number of boundary conditions of the problem.

  • type : This is the type of boundary condition been imposed. At the moment, choices are

    • temperature (Dirichlet BC), to impose a given temperature value at the boundary
    • convection (Robin BC) for cooling/heating, depending on the environment temperature at the boundary Tinf, with a given heat transfer coefficient h
  • NB : by default, the boundary condition type is temperature with value = 0

Tracer

subsection boundary conditions tracer
  set number                  = 1
  subsection bc 0
        set type              = dirichlet
        subsection dirichlet
            set Function expression = x
        end
  end
end
  • number: This is the number of boundary conditions of the problem.

  • type : This is the type of boundary condition been imposed. At the moment, only dirichlet boundary conditions can be imposed for tracer.

Source term

Possible present source terms can be added to the simulation using this subsection. With enabling the source term by setting the set enable = true, the specified function expression is added as a source term to the simulation.

subsection source term
  # Enable the calculation of the analytical solution and L2 error
  set enable = true


  subsection xyz
    # Sometimes it is convenient to use symbolic constants in the expression
    # that describes the function, rather than having to use its numeric value
    # everywhere the constant appears. These values can be defined using this
    # parameter, in the form `var1=value1, var2=value2, ...'.
    # 
    # A typical example would be to set this runtime parameter to
    # `pi=3.1415926536' and then use `pi' in the expression of the actual
    # formula. (That said, for convenience this class actually defines both
    # `pi' and `Pi' by default, but you get the idea.)
    set Function constants  = 

    # The formula that denotes the function you want to evaluate for
    # particular values of the independent variables. This expression may
    # contain any of the usual operations such as addition or multiplication,
    # as well as all of the common functions such as `sin' or `cos'. In
    # addition, it may contain expressions like `if(x>0, 1, -1)' where the
    # expression evaluates to the second argument if the first argument is
    # true, and to the third argument otherwise. For a full overview of
    # possible expressions accepted see the documentation of the muparser
    # library at http://muparser.beltoforion.de/.
    # 
    # If the function you are describing represents a vector-valued function
    # with multiple components, then separate the expressions for individual
    # components by a semicolon.
    set Function expression = 0; 0; 0 

    # The names of the variables as they will be used in the function,
    # separated by commas. By default, the names of variables at which the
    # function will be evaluated are `x' (in 1d), `x,y' (in 2d) or `x,y,z' (in
    # 3d) for spatial coordinates and `t' for time. You can then use these
    # variable names in your function expression and they will be replaced by
    # the values of these variables at which the function is currently
    # evaluated. 
    set Variable names      = x,y,t
  end
end

Analytical Solution

If the problem being simulated has a known analytical solution, it can be added in this section by first setting the set enable = true. This allows for calculation of the error in the simulation by comparing the analytical solution with the simulated results.

subsection analytical solution
  set enable                 = true
    subsection uvwp
            set Function constants = eta=0.25, ri=0.25
            set Function expression = eta*x; ri*y ; 0
    end
end
  • subsection uvwp allows to define the expected expression for the velocity field (u in x-direction, v in y-direction and w in z-direction) and the pressure field p.
    • Function constants can be defined and use for a more readable Function expression
    • NB: Function expression requires three entry for 2d-simulations (uvp-solution) and four for 3d-simulations (uvwp-solution)

Timer

The timer section controls the frequency at which the timing of the simulation is output. Setting its value to iteration outputs the timing of all elements of the simulation at the end of each iteration. Setting its value to end outputs the timing all elements only at the end of the simulation. The timing can be disabled completely, which is necessary for the unit and functional tests.

subsection timer
  # Clock monitoring methods Choices are none, iteration or end
  set type = none
end

Mesh Adaptation Control

This subsection controls the mesh adaptation method.

subsection mesh adaptation
  # Fraction of coarsened elements
  set fraction coarsening  = 0.05

  # Fraction of refined elements
  set fraction refinement  = 0.1

  # How the fraction of refinement/coarsening are interpreted
  # Choices are number or fraction 
  set fraction type        = number

  # Frequency of the mesh refinement
  set frequency            = 1

  # Maximum number of elements
  set max number elements  = 100000000

  # Maximum refinement level
  set max refinement level = 10

  # Minimum refinement level
  set min refinement level = 0

  # Type of mesh adaptation. Choices are  none, uniform or kelly.
  set type                 = none

  # Variable for kelly estimation. Choices are velocity or pressure.
  set variable             = velocity
end
  • Two type of mesh adaptation are available. The uniform mesh adaptation refines the mesh at every cell, whereas the kelly uses a kelly error estimator to decides which cell are refined.
  • The fraction of cell that are refined and coarsened are controlled with the fraction refinement and fraction coarsening parameters. Note that a cell cannot be coarsened more than it's initial level. Consequently, adaptively refined simulations should start with a mesh as coarse as possible.
  • The maximal and minimal refinement level reachable for a cell are controlled with the max refinement and min refinement parameters.
    • for deal.ii meshes, min refinement level should be equal to the initial refinement
    • for gmsh imported meshes, min refinement level = 0
    • for a good compromise between speed and precision, max refinement level should be set to 2 or 3 more than the min refinement level

If the option kelly is chosen, variables for use for the kelly estimation should be specified:

  • velocity
  • pressure

Restart

This subsection controls the restart and checkpointing of the simulation.

subsection restart
  # Enable checkpointing
  set checkpoint = true

  # Enable restart from previous checkpoint
  set restart    = false

  # Prefix for the filename of checkpoints
  set filename   = restart

  # Frequency for checkpointing
  set frequency  = 1

end
  • The checkpoint option enables (true) or disables (false) checkpointing.
  • The restart option controls if the simulation will be restarted from a previous checkpoint. Turn to false for the first checkpoint creation. The simulation will start from the last checkpoint and will finish at the provided time end (see Simulation Control).
  • The filename option fixes the prefix of all the checkpointing files
  • The frequency controls the checkpointing frequency. Checkpointing is very intensive from an I/O point of view, consequently it should not be done at every iteration.

Non-Linear Solver Control

The Navier-Stokes equations (and other) are non-linear equation. This section controls the non-linear solver used within lethe.

subsection non-linear solver
  # Maximum number of Newton Iterations
  set max iterations     = 10

  # Non-linear solver that will be used. Choices is newton.
  set solver             = newton

  # Newton solver tolerance
  set tolerance          = 1e-6

  # Tolerance for the acceptation of the step
  set step tolerance     = 0.9

  # State whether from the non-linear solver should be printed Choices are
  # .
  set verbosity          = verbose
end
  • Lethe currently has one solver for non-linear equations. The newton solver is a Newton-Raphson solver which recalculates the jacobian matrix at every iteration.
  • The verbosity options enables the display of the residual at each non-linear iteration. This allow ones to monitor how the non-linear iterations are progressing.
  • The step tolerance option controls how much the L2 norm of the residual must decrease for a non-linear step to be accepted. If the new_residual < old_residual * step tolerance, than a newton iteration is accepted. If this condition is not reached, then a theta weighting of the step is applied until this condition is reached.

Linear Solver Control

In these subsections, the control options of linear solvers are specified. These control options may include solution method, maximum number of iterations, tolerance, residual precision and preconditioner information.

Lethe currently supports three core solution strategy or method :

  • GMRES solver + ILU preconditioning. This is set by method = gmres
  • GMRES solver + AMG preconditioning. This is set by method = amg
  • BiCGStab solver + ILU preconditioning. This is set by method = bicgstab

For all methods, a number of parameters to control the preconditioner can be set.

On coarse meshes, the gmres/bicgstab solver with ILU preconditioning is more efficient. As the number of mesh elements increase, the amg solver is the most efficient. Generally, at 1M elements, the amg solver always outperforms the gmres or bicgstab for 2D steady-state problems.

subsection linear solver
  # amg aggregation threshold
  set amg aggregation threshold                 = 1e-14

  # amg number of cycles
  set amg n cycles                              = 1

  # amg preconditioner ilu smoother/coarsener absolute tolerance
  set amg preconditioner ilu absolute tolerance = 1e-12

  # amg preconditioner ilu smoother/coarsener fill
  set amg preconditioner ilu fill               = 1

  # amg preconditioner ilu smoother/coarsener relative tolerance
  set amg preconditioner ilu relative tolerance = 1.00

  # amg smoother overlap
  set amg smoother overlap                      = 1

  # amg smoother sweeps
  set amg smoother sweeps                       = 2

  # amg w cycling. If this is set to true, W cycling is used. Otherwise, V
  # cycling is used.
  set amg w cycles                              = false

  # Ilu preconditioner tolerance
  set ilu preconditioner absolute tolerance     = 1e-6

  # Ilu preconditioner fill
  set ilu preconditioner fill                   = 1

  # Ilu relative tolerance
  set ilu preconditioner relative tolerance     = 1.00

  # Maximum solver iterations
  set max iters                                 = 1000

  # Maximum number of krylov vectors for GMRES solver
  set max krylov vectors                        = 30

  # The iterative solver for the linear system of equations. Choices are
  # . gmres is a GMRES iterative solver with ILU
  # preconditioning. bicgstab is a BICGSTAB iterative solver with ILU
  # preconditioning. amg is GMRES + AMG preconditioning with an ILU coarsener
  # and smoother. On coarse meshes, the gmres/bicgstab solver with ILU
  # preconditioning is more efficient. As the number of mesh elements
  # increase, the amg solver is the most efficient. Generally, at 1M elements,
  # the amg solver always outperforms the gmres or bicgstab
  set method                                    = gmres

  # Linear solver minimum residual
  set minimum residual                          = 1e-8

  # Linear solver residual
  set relative residual                         = 1e-3

  # State whether output from solver runs should be printed. Choices are
  # .
  set verbosity                                 = verbose
end

Post-Processing

This subsection controls the post-processing options. For instance, the user can define to obtain enstrophy, kinetic energy or pressure drop in the post-processing. When calculate enstrophy, calculate kinetic energy and/or pressure drop is enabled, text files that are written for the various variables (e.g. enstrophy, kinetic energy, etc.) at a frequency controlled by output frequency. Those post-processing variables are not expensive to calculate and are thus calculated at every iteration. It is also possible to obtain time-averaged velocities and pressure and the Reynolds stresses when calculate average velocities is enable. However, the averaged solutions and Reynolds stresses are not written in text files but are displayable using visualization software as other solutions. initial time has to be a time where the simulation is stable, this is the time calculating time-averaged solutions begins.

subsection post-processing
  # Enable calculation of total enstrophy
  set calculate enstrophy         = false

  # Enable calculation of total kinetic energy
  set calculate kinetic energy    = false

  # Enable calculation of pressure drop
  set calculate pressure drop     = false

  # Boundary id of the inlet for the pressure drop calculation
  set inlet boundary id           = 0

  # Boundary id of the outlet for the pressure drop calculation
  set outlet boundary id          = 1

  # Calculation frequency
  set calculation frequency       = 1

  # File output enstrophy
  set enstrophy name              = enstrophy

  # File output kinetic energy
  set kinetic energy name         = kinetic_energy

  # File output pressure drop
  set pressure drop name          = pressure_drop

  # Output frequency
  set output frequency            = 1

  # State whether from the post-processing values should be printed Choices
  # are .
  set verbosity                   = quiet

  # Enable time-averaged velocities and pressure
  set calculate average velocities = true

  # Initial time of the averaging
  set initial time                 = 10
end

Force and Torques calculation

This subsection controls the post-processing of the forces and the torque on the boundary conditions. For instance, the user can obtain the force acting on a wall or the torque acting on an object­. This is enabled by setting the parameters calculate forces and calculate torques to true. The frequency at which this is written to a text file is controlled by the output frequency parameter. The force and the torque on all boundary condition is calculated and written to it's own text file. This text file will bear the name of the parameter force name with a suffix ".boundary_id" where boundary id is the number associated with the boundary condition.

subsection forces
  # Enable calculation of forces
  set calculate forces      = true

  # Enable calculation of torques
  set calculate torques     = true

  # Calculation frequency
  set calculation frequency = 1

  # File output force prefix
  set force name            = force

  # Output frequency
  set output frequency      = 1

  # Calculation frequency
  set output precision      = 10

  # File output torque prefix
  set torque name           = torque

  # State whether from the non-linear solver should be printed (quiet or verbose).
  set verbosity             = quiet
end

Dynamic flow control

This subsection's purpose is to enable a dynamic flow control. When a periodic boundary condition is set, the flow must be controlled at inlet. It calculates a Beta coefficient at each time step and may be added to a source term. If the chosen wall boundary to force the flow is at inlet, the volumetric flow rate targeted must be negative since outward normal vector is used. Therefore, if the chosen wall boundary is at outlet, volumetric flow rate must be positive.

subsection flow control
 # Enable the dynamic flow control
 set enable = true

 # Volumetric flow rate targeted (negative if at inlet)
 set volumetric flow rate = -91.75
 
 # Wall boundary where the flow is adjusted
 set boundary id = 0

 # Flow direction
 set flow direction = 0

 # Initial Beta coefficient
 set initial beta = 100
 
end
  • The volumetric flow rate parameter set a flow rate where the algorithm will target at each time step. Important : when there is a periodic boundary, the wall boundary should be at inlet only and if the chosen wall boundary is at inlet, flow rate must be negative and vice versa.
  • The boundary id parameter is the wall boundary where the flow rate intended is controlled. (inlet or outlet).
  • The flow direction parameter set the absolute direction. 0 = x, 1 = y and 2 = z.
  • The initial beta parameter set a initial beta coefficient that may speed up the convergence of the flow rate targeted. Some tests should be done to find a good initial beta coefficient.

Visualization

Once the solver executable and parameter file are provided, the simulation can be started. After finishing the simulation, user can open the output files (with extension .vtu in the test folder) using visualization software, such as ParaView.

Void Fraction

This subsection is required for the gas_vans solver as well as for the cfd-dem coupling. It allows the user to specify the method for void fraction calculation when there exists more than one phase in the simulation domain. Based on the method chosen, further parameters may be required as shown the the section below.

subsection void fraction
   set mode = dem
   set read dem = true
   set dem file name = dem
   set l2 smoothing factor = 0.0
   set l2 lower bound = 0.35
   set l2 upper bound = 1
end
  • The mode parameter allows the user to choose the method for void fraction calculation. Currently, there are two methods implemented. The first one is to calculate the void fraction from function. In this case, an additional subsection is required to insert the function such as:

     subsection function
     set Function expression = 0.5 + 0.25*sin(pi*x)*sin(pi*y)
     end

If the mode chosen is dem, then the void fraction is calculated using the Particle Centered Method. In this method, the remaining parameters are required.

  • The read dem allows us to read an already existing dem simulation result which can be obtained from checkpointing the Lethe-DEM simulation. This is important as the gls_vans solver requires reading an initial dem triangulation and particle information to simulate flows in the presence of particles.
  • The dem_file_name parameter specifies the prefix of the dem files that must be read.
  • The l2 smoothing factor is a smoothing length used for smoothing the L2 projection of the void fraction to avoid sharp discontinuities which can lead to instabilities in the simulation.
  • The l2 lower bound and l2 upper bound are the minimum and maximum values around which the void fraction is bounded. This is important especially for upper bounds as the void fraction can sometimes slightly exceed a value of 1 when projected.

CFD-DEM

This subsection includes parameters related to multiphase flow simulations using the both the gas_vans solver and the cfd-dem solver within Lethe.


subsection cfd-dem
   set grad div = true
   set drag model = rong
   set post processing = true
   set coupling frequency = 100
end
  • The grad div parameter allows the enabling of the grad div stabilization for the Volume Averaged Navier Stokes equations. This allows for a much better mass conservation of the system.
  • The drag model parameter allows to choose the type of drag model to be implemented for the calculation of the drag force between the particles and the fluids.
  • The post processing parameter, when enabled allows the calculation of the pressure drop, void fraction in the packed region, and the mass conservation in a packed bed.
  • The coupling frequency parameter is only applicable for the cfd-dem solver and it determines the number of DEM iterations per 1 CFD iteration.
Clone this wiki locally