Parameter file
Summary
Many functionalities of an AMDiS
program can be controlled using a parameter file. This file is passed as first argument to the AMDiS
executable. If we want to run the ellipt.2d
example using the parameter file init/ellipt.dat.2d
we can do so using the following console command in the examples directory:
./ellipt.2d init/ellipt.dat.2d [other arguments]
We will give a list of possible parameters at the end of this document.
Structure of a parameter file
Parameters in the parameter file have keys in the form of strings separated by ->
. This allows grouping of parameters that are related to the same things. High-level classes like ProblemStat
have a string constructor argument where a name can be provided that relates to an initfile key.
examples/ellipt.cc:
using Param = LagrangeBasis<Grid, 2>;
ProblemStat<Param> prob("ellipt"); // make a ProblemStat with the name 'ellipt'
examples/init/ellipt.dat.2d:
ellipt->mesh: elliptMesh % name of the mesh created by 'ellipt'
elliptMesh->num cells: 8 8 % size of the mesh 'elliptMesh'
ellipt->solver: pcg % solver type used by 'ellipt'
ellipt->solver->max iteration: 10000 % solver parameter related to 'ellipt->solver'
As we saw above, parameters are separated from the keys with the :
character and can be strings and numbers. Anything following a %
character is considered a comment.
Reading from a parameter file
In the section above we learned that AMDiS
classes can be given names that can then be used to specify parameters in the parameter file. In addition to the options AMDiS
classes already provide we can add user-defined parameters. For this we simply have to add a key: value
pair to the parameter file as we would do with the built-in options and read the value somewhere in our program using the function Parameters::get<T>
:
// (1) Overwrites 'value' with the value related to 'key' if found
template <class T>
static void get(std::string const& key, T& value);
// (2) Returns a std::optional containing the value related to 'key' if found
template <class T>
static std::optional<T> get(std::string const& key);
We will now give an example of reading an int
parameter called myParameter
from an initfile that defaults to 3
if no parameter file containing the key is provided.
// example parameter file entry
myParameter: 42
// using (1)
int i = 3;
Parameters::get("myParameter", i);
// using (2)
int i = Parameters::get<int>("myParameter").value_or(3);
Writing to a parameter file
It is possible to change to value of a parameter given in a parameter file using the function
template <class T>
static void set(std::string const& key, T const& value);
This can be useful if we want to override certain values within our program.
List of parameters and values
In the following section <prefix>
denotes a placeholder for the name given to the class instance. <treepath>
is a (possibly empty) comma-separated list of integers specifying a TreePath
.
ProblemStat
key | type | values | default value |
---|---|---|---|
<prefix>->mesh |
string | Prefix used by the mesh | mesh |
<prefix>->restore->grid |
string | Filename of the grid file | Not optional1 |
<prefix>->restore->solution |
string | Filename of the solution file | Not optional1 |
<prefix>->symmetry |
string | Symmetry of the stiffness matrix, one of:spd ,symmetric ,hermitian ,structurally_symmetric ,unknown |
unknown |
<prefix>->solver |
string | Solver type, see Solvers and Preconditioners | default |
<prefix>->marker[<treepath>] |
string | Prefix for marker parameters, see Markers | None2 |
<prefix>->output[<treepath>] |
string | Prefix for filewriter parameters, see Filewriters | None3 |
<prefix>->output |
string | as <prefix>->output[] |
None3 |
1Only required when ProblemStat::restore
is called
2No marker will be set up if the ->strategy
subkey is not specified
3No filewriter will be set up if the ->format
subkey is not specified
Markers
key | type | values | default value |
---|---|---|---|
<prefix>->strategy |
string | Marking strategy, one of:0 (no marking),1 (global refinement),2 (maximum strategy),3 (equidistribution strategy),4 (guaranteed error reduction strategy) |
Not optional |
<prefix>->info |
int | Verbosity | 0 |
<prefix>->max refinement level |
int | Maximum refinement level | No maximum |
<prefix>->min refinement level |
int | Minimum refinement level | 0 |
<prefix>->p 123 |
double | Power in estimator norm | 2.0 |
<prefix>->MSGamma 1 |
double | Marking parameter | 0.5 |
<prefix>->MSGammaC 1 |
double | Marking parameter | 0.1 |
<prefix>->ESTheta 2 |
double | Marking parameter | 0.9 |
<prefix>->EsThetaC 2 |
double | Marking parameter | 0.2 |
<prefix>->GERSThetaStar 3 |
double | Marking parameter | 0.6 |
<prefix>->GERSNu 3 |
double | Marking parameter | 0.1 |
<prefix>->GERSThetaC 3 |
double | Marking parameter | 0.1 |
1Only for maximum strategy
2Only for equidistribution strategy
3Only for guaranteed error reduction strategy
Solvers and Preconditioners
Parameters used by the SolverInfo class
key | type | values | default value |
---|---|---|---|
<prefix>->info |
int | Solver verbosity level | 0 |
<prefix>->break if tolerance not reached |
bool | If true throw an error if tolerance could not be reached |
false |
PETSc
key | type | values | default value |
---|---|---|---|
<prefix>->info |
int | Verbosity, one of:0 (no convergence information),1 (minimal output, print residual every 10th iteration),2 (full convergence output) |
0 |
<prefix>->ksp |
string | Name of a PETSc Krylov methoda | default |
<prefix> |
string | alternative to <prefix>->ksp when not set |
as above |
<prefix>->max it |
PetscInt | Verbosity | PETSC_DEFAULT |
<prefix>->rtol |
PetscReal | Relative convergence toleranceb | PETSC_DEFAULT |
<prefix>->atol |
PetscReal | Absolute convergence toleranceb | PETSC_DEFAULT |
<prefix>->divtol |
PetscReal | Divergence toleranceb | PETSC_DEFAULT |
<prefix>->prefix |
string | The options prefix for the solverc | No prefix set |
<prefix>->mat solver |
string | Sets the software that is used to perform the factorization d | umfpack (sequential matrix) or mumps (otherwise) |
<prefix>->pc |
string | The preconditioner typee | default |
<prefix>->pc->mat solver 1 |
string | Solver used by the LU preconditionerd | None |
<prefix>->pc->prefix |
string | The options prefix for the preconditionerf | No prefix set |
<prefix>->scale 2 |
PetscReal | Damping factor used by richardson solver | 1.0 |
<prefix>->restart 3 |
PetscInt | Number of iterations until restart | 30 |
<prefix>->gramschmidt 3 |
string | Orthogonalization routine used by gmres solversg: classical or modified |
Use PETSc default |
<prefix>->restart 3 |
string | Type of iterative refinement to use in the classical Gram Schmidt orthogonalizationh, one of:never (never refine),ifneeded (refine if needed),always (always refine) |
Use PETSc default |
<prefix>->ell 4 |
PetscReal | Number of search directions in BiCGStab_ell | Use PETSc default |
<prefix>->restart 5 |
PetscInt | Number of iterations until restart | Use PETSc default |
1Only for lu
preconditioner
2Only for richardson solver
3Only for gmres solver types
4Only for BiCGStab_ell solver
5Only for GCR solver
aSee PETSc documentation for Krylov methods
bSee PETSc documentation for Krylov method tolerances
cSee PETSc documentation for Krylov method options
dSee PETSc documentation for preconditioners
eSee PETSc documentation for preconditioners
fSee PETSc documentation for preconditioners
gSee PETSc documentation for GMRes solver
hSee PETSc documentation for GMRes solver
EIGEN
key | type | values | default value |
---|---|---|---|
<prefix> |
string | Solver type, one of:cg (conjugate gradient solver),bicgstab (bi-conjugate gradient stabilized solver),bcgs (same as bicgstab ),minres (minimal residual algorithm for symmetric matrices),gmres (generalized minimal residual algorithm),dgmres (restarted gmres with deflation),default (same as gmres ),umfpack 1 (sparse LU factorization based on UmfPack),superlu 2 (sparse LU factorization based on SuperLU library),direct (same as umfpack or superlu ) |
default |
<prefix>->reuse pattern |
bool | Reuse matrix pattern if true |
false |
<prefix>->relative tolerance |
RealScalar | Tolerance threshold used by the stopping criteria | floating-point epsilon |
<prefix>->max iteration |
int | Maximum number of iterations | 500 |
<prefix>->restart 34 |
int | Number of iterations until restart | 30 |
<prefix>->num eigenvalues 4 |
int | Number of eigenvalues to deflate at each restart | 0 |
<prefix>->max num eigenvalues 4 |
int | Maximum size of the deflation subspace | 5 |
<prefix>->precon->name 5 |
string | Preconditioner, one of:no (no preconditioner),diag (diagonal preconditioner),jacobi (diagonal preconditioner),ilu (incomplete LU factorization),ic (incomplete Choleski factorization) |
no |
<prefix>->precon->ordering 6 |
string | Ordering used in Choleski factorization: amd or natural |
amd |
<prefix>->precon->drop tolerance 7 |
double | Drop any element whose magnitude is less than this tolerance | 1e-15l |
<prefix>->precon->fill factor 7 |
int | This is used to compute the number of largest elements to keep on each row | 10 |
1Requires UmfPack
2Requires SuperLU
3Only for gmres
solver
4Only for dgmres
solver
5Does not apply to umfpack
or superlu
solvers
6Only for incomplete Choleski factorization
7Only for ilu
solver
MTL
key | type | values | default value |
---|---|---|---|
<prefix> |
string | Solver type, one of:cg (conjugate gradient method),cgs (stabilized conjugate gradient method),bicg (BiCG method),bcgs (stabilized BiCG(1) method),bicgstab (as bcgs ),bcgs2 (stabilized BiCG(2) method),bicgstab2 (as bcgs2 ),qmr (quasi-minimal residual method),tfqmr (transposed-free quasi-minimal residual method),bcgsl (stabilized BiCG(ell) method),bicgstab_ell (as bcgsl ),gmres (generalized minimal residual method),fgmres (flexible gmres),minres (minimal residual method),idrs (induced dimension reduction method),gcr (generalized conjugate residual method),preonly (just apply a preconditioner),default (as gmres ),umfpack 1 (external UMFPACK solver),direct (as umfpack ) |
default |
<prefix>->absolute tolerance 2 |
type of matrix values | Absolute solver tolerance: ∥b-Ax∥ < absolute tolerance |
0 |
<prefix>->relative tolerance 2 |
type of matrix values | Relative solver tolerance: ∥b-Ax∥ / ∥b-Ax0∥ < relative tolerance |
1.e-6 |
<prefix>->max iteration 2 |
int | Maximal number of iterations | 1000 |
<prefix>->print cycle 2 |
int | Print solver residuals every print cycle iterations |
100 |
<prefix>->store symbolic 3 |
bool | false |
|
<prefix>->symmetric strategy 3 |
int | 0 |
|
<prefix>->alloc init 3 |
double | 0.7 |
|
<prefix>->ell 4 |
int | Parameter of stabilized BiCG(ell ) |
3 |
<prefix>->restart 5 |
int | Maximal number of orthogonalized vectors | 30 |
<prefix>->orthogonalization 6 |
int | 1 (gram schmidt) or 2 (householder) |
1 |
<prefix>->s 7 |
int | Parameter of IDR(s ) |
30 |
<prefix>->precon 2 |
string | Left preconditioner type, one of:diag (diagonal preconditioner),jacobi (as diag ),masslumping (mass lumping preconditioner),ilu (incomplete LU preconditioner),ilu0 (as ilu ),no (no preconditioner),identity (as no ),solver (use a solver as preconditioner),default (as no ),hypre 8 (Hypre preconditioner) |
default |
<prefix>->precon->solver 9 |
int | Solver type used for the preconditioner, same options as for <prefix> |
default |
1Requires UmfPack
2Does not apply to umfpack
solver
3Only umfpack
solver
4Only bcgsl
solver
5Only gmres
, fgmres
and gcr
solvers
6Only gmres
and fgmres
solvers
7Only idrs
solver
8Requires Hypre, see amdis/linearalgebra/mtl/HyprePrecon.hpp
for subkeys and possible values
9Only solver
preconditioner, subkeys use the same values as above e.g. <prefix>->precon->solver->max iteration: 500
ISTL
key | type | values | default value |
---|---|---|---|
<prefix> |
string | Solver type, one of:cg (conjugate gradient method),pcg (generalized preconditioned conjugate gradient solver),fcg 1 (accelerated flexible conjugate gradient method),cfcg 1 (complete flexible conjugate gradient method),bcgs (stabilized bi-conjugate gradient method),bicgstab (as bcgs ),default (as bcgs ),minres (minimal residul method),gmres (generalized minimal residual method),fgmres 1 (flexible generalized minimal residual (FGMRes) method),umfpack 2 (external UMFPACK solver),ldl 2 (external LDL solver),spqr 2 (external SQPR solver),cholmod 2 (external Cholmod solver),superlu 2 (external SuperLU solver),direct (as umfpack or superlu ) |
default |
<prefix>->info |
int | Verbosity | 0 |
<prefix>->maxit 3 |
int | Maximal number of solver iterations | 500 |
<prefix>->reduction 3 |
type of matrix values | Relative break tolerance | 1.e-6 |
<prefix>->restart 4 |
int | Restart parameter | 30 |
<prefix>->reuse vector 5 |
bool | Reuse vectors in subsequent calls to apply if true |
true |
<prefix>->category |
string | sequential (sequential solver),overlapping (overlapping parallel solver),nonoverlapping (nonoverlapping parallel solver),default (chooses depending on grid and process count) |
None9 |
<prefix>->precon 3 |
string | Preconditioner type, one of:diag (diagonal preconditioner),jacobi (as diag ),gs (Gauss-Seidel preconditioner),gauss_seidel (as gs ),sor (Successive Overrelaxation methods),ssor (Symmetric Successive Overrelaxation methods),pssor (A parallel SSOR preconditioner (requires overlap)),richardson (Richardson methods),default (as richardson ),solver (Turns a solver into a preconditioner),bjacobi (Block-Jacobi methods),ilu (Incomplete LU factorization),ilu0 (as ilu ),ildl (Incomplete LDL factorization),amg 6 (Algebraic multigrid method),fastamg 6 (Algebraic multigrid method),kamg 6 (Algebraic multigrid method) |
default |
<prefix>->precon->relaxation 3 |
double | Dumping/relaxation factor | 1.0 |
<prefix>->precon->iterations 3 |
int | Number of iterations | 1 |
<prefix>->precon->solver 7 |
string | Solver type used for the preconditioner, same options as for <prefix> |
default |
<prefix>->precon->sub precon 8 |
string | Preconditioner used in each block, same options as for <prefix>->precon |
default |
solver category |
string | as <prefix>->category |
default |
1Requires dune-istl
2.7
2Requires external solver package
3Only for cg
, pcg
, fcg
, cfcg
, bcgs
, minres
, gmres
and fgmres
solvers
4Only for pcg
and gmres
solvers
5Only for superlu
solver
6See amdis/linearalgebra/istl/AMGPrecon.hpp
for subkeys and possible values
7Only solver
preconditioner, subkeys use the same values as above e.g. <prefix>->precon->solver->max iteration: 500
8Only for bjacobi
preconditioner, subkeys use the same values as above e.g. <prefix>->precon->sub precon->iterations: 1
9Checks global parameter solver category
if not set
Grids
key | type | values | default value |
---|---|---|---|
<prefix>->macro file name |
string | Filename of the mesh file | None1 |
<prefix>->structured 2 |
string | Type of the structured grid, one of:cube (rectangular or cube grid),simplex (triangular or tetrahedral grid) |
Grid dependant |
<prefix>->global refinements |
int | Number of initial global refinement steps | 0 |
<prefix>->load balance |
bool | Perform initial load balance if true |
No load balance |
<prefix>->min corner 2 |
dim-sized array of global coordinates | lower left corner of the domain | (0,...,0) |
<prefix>->max corner 2 |
dim-sized array of global coordinates | upper right corner of the domain | (1,...,1) |
<prefix>->num cells 2 |
dim-sized array of int | number of blocks in each coordinate direction | (1,...,1) |
<prefix>->overlap 3 |
int | Number of overlap elements on a distributed grid | 0 |
<prefix>->periodic 3 |
string | String of ones and zeros indicating if the grid is periodic in the xyz-direction, e.g. 010 periodic in y-direction |
0...0 (not periodic) |
1If no file name is given a structured grid will be constructed
2Only for structured grid generation
3YaspGrid
only
Adaptivity
Parameters used by the AdaptInfo
class
key | type | values | default value |
---|---|---|---|
<prefix>->tolerance |
double | Tolerance for the (absolute or relative) error | 0.0 |
<prefix>->time tolerance |
double | Time tolerance | 0.0 |
<prefix>->time relative tolerance |
double | Relative time tolerance | 0.0 |
<prefix>->coarsen allowed |
int | 0 (coarsening not allowed),1 (coarsening allowed) |
0 |
<prefix>->refinement allowed |
int | 0 (refinement not allowed),1 (refinement allowed) |
0 |
<prefix>->sum factor |
double | Factor to combine max and integral time estimate | 1.0 |
<prefix>->max factor |
double | Factor to combine max and integral time estimate | 0.0 |
<prefix>->start time |
double | Initial time | 0.0 |
<prefix>->timestep |
double | Time step size to be used | 0.0 |
<prefix>->end time |
double | Final time | 0.0 |
<prefix>->max iteration |
int | Maximal allowed number of iterations of the adaptive procedure | No maximum |
<prefix>->max timestep iteration |
int | Maximal number of iterations for choosing a timestep | 30 |
<prefix>->max time iteration |
int | Maximal number of time iterations | 30 |
<prefix>->min timestep |
double | Minimal step size | Square root of floating-point epsilon |
<prefix>->max timestep |
double | Maximal step size | Square root of floating-point maximal value |
<prefix>->number of timesteps |
int | Number of fixed timestep iterations | None1 |
<prefix>->time tolerance |
double | Tolerance for the overall time error | 1.0 |
1If this is set to a nonzero value the computation is done number of timesteps
times with a fixed time step size
Parameters used by the AdaptInstationary
class
key | type | values | default value |
---|---|---|---|
<prefix>->strategy |
int | 0 (explicit time strategy),1 (implicit time strategy),2 (simple adaptive time strategy) |
0 |
<prefix>->time delta 1 |
double | Parameter Δ1 used in time step reduction | 1/√2 |
<prefix>->time delta 2 |
double | Parameter Δ2 used in time step enlargement | √2 |
<prefix>->break when stable |
bool | Stops the iteration when the instationary problem is stable when set to true |
false |
Filewriters
key | type | values | default value |
---|---|---|---|
<prefix>->format |
vector<string> | List of format strings from the following list1:vtk (VTK format using the dune-grid writer),dune-vtk (VTK format using the dune-vtk writer)2,gmsh (GMsh file format),backup (backup writer) |
Not optional |
<prefix>->filename |
string | Filename of the output file excluding the file extension | solution |
<prefix>->output directory |
string | Directory where to put the files | . |
<prefix>->name |
string | Name of the data vector in the output file | solution |
<prefix>->write every i-th timestep |
int | Timestep number interval | 0 |
<prefix>->write after timestep |
double | Time interval | 0.0 |
<prefix>->animation 345 |
bool | write .pvd file if true |
false |
<prefix>->mode 3 |
int | 0 (output in ASCII),1 (output appended base64 binary) |
0 |
<prefix>->subsampling 3 |
int | Number of refinement intervals for subsampling | No subsampling |
<prefix>->mode 4 |
int | 0 (ASCII),1 (binary),2 (compressed) |
0 |
<prefix>->precision 45 |
int | 0 (float),1 (double) |
0 |
<prefix>->animation 6 |
bool | Append timestepnumber if true , otherwise always override |
false |
1One filewriter will be set up for each format string in the list
2Requires the dune-vtk
module
3Only for vtk
format
4Only for dune-vtk
format
5Only for gmsh
format
6Only for backup
format