Chi-Tech

Modules

 lbs.LBSGroupset
 

Functions

void chi_lua::chiLBSCreateGroupset (int SolverIndex)
 
void chi_lua::chiLBSCreateGroup (int SolverIndex, int GroupId)
 
void chi_lua::chiLBSGroupsetAddGroups (int SolverIndex, int GroupsetIndex, int FromIndex, int ToIndex)
 
void chi_lua::chiLBSGroupsetSetQuadrature (int SolverIndex, int GroupsetIndex, int QuadratureIndex)
 
void chi_lua::chiLBSGroupsetSetAngleAggregationType (int SolverIndex, int GroupsetIndex, int AggregationType)
 
void chi_lua::chiLBSGroupsetSetAngleAggDiv (int SolverIndex, int GroupsetIndex, int NumDiv)
 
void chi_lua::chiLBSGroupsetSetGroupSubsets (int SolverIndex, int GroupsetIndex, int NumDiv)
 
void chi_lua::chiLBSGroupsetSetIterativeMethod (int SolverIndex, int GroupsetIndex, int IterativeMethod)
 
void chi_lua::chiLBSGroupsetSetResidualTolerance (int SolverIndex, int GroupsetIndex, float ResidualTol)
 
void chi_lua::chiLBSGroupsetSetMaxIterations (int SolverIndex, int GroupsetIndex, int Numiter)
 
void chi_lua::chiLBSGroupsetSetGMRESRestartIntvl (int SolverIndex, int GroupsetIndex, int Intvl)
 
void chi_lua::chiLBSGroupsetSetEnableSweepLog (int SolverIndex, int GroupsetIndex, bool flag)
 
void chi_lua::chiLBSGroupsetSetWGDSA (int SolverIndex, int GroupsetIndex, int MaxIters, float ResTol, bool Verbose, char PETSCString)
 
void chi_lua::chiLBSGroupsetSetTGDSA (int SolverIndex, int GroupsetIndex, int MaxIters, float ResTol, bool Verbose, char PETSCString)
 

Detailed Description

The code below is an example of a complete specification of a groupset.

--===================================== Setup physics
phys1 = chiLBSCreateSolver()
chiSolverAddRegion(phys1,region1)
chiLBSSetProperty(phys1,DISCRETIZATION_METHOD,PWLD3D)
chiLBSSetProperty(phys1,SCATTERING_ORDER,1)
--========== Groups
grp = {}
for g=1,num_groups do
grp[g] = chiLBSCreateGroup(phys1)
end
--========== ProdQuad
pquad0 = chiCreateProductQuadrature(GAUSS_LEGENDRE_CHEBYSHEV,2, 2)
pquad1 = chiCreateProductQuadrature(GAUSS_LEGENDRE_CHEBYSHEV,8, 8)
--========== Groupset def
cur_gs = gs0
chiLBSGroupsetAddGroups(phys1,cur_gs,0,15)
chiLBSGroupsetSetQuadrature(phys1,cur_gs,pquad0)
chiLBSGroupsetSetIterativeMethod(phys1,cur_gs,NPT_GMRES)
chiLBSGroupsetSetWGDSA(phys1,cur_gs,30,1.0e-4,false," ")
chiLBSGroupsetSetTGDSA(phys1,cur_gs,30,1.0e-4,false," ")
void chiLBSSetProperty(int SolverIndex, int PropertyIndex)
void chiLBSGroupsetAddGroups(int SolverIndex, int GroupsetIndex, int FromIndex, int ToIndex)
void chiLBSGroupsetSetQuadrature(int SolverIndex, int GroupsetIndex, int QuadratureIndex)
void chiLBSGroupsetSetTGDSA(int SolverIndex, int GroupsetIndex, int MaxIters, float ResTol, bool Verbose, char PETSCString)
void chiLBSCreateGroup(int SolverIndex, int GroupId)
void chiLBSGroupsetSetMaxIterations(int SolverIndex, int GroupsetIndex, int Numiter)
void chiLBSGroupsetSetAngleAggDiv(int SolverIndex, int GroupsetIndex, int NumDiv)
void chiLBSGroupsetSetGroupSubsets(int SolverIndex, int GroupsetIndex, int NumDiv)
void chiLBSGroupsetSetResidualTolerance(int SolverIndex, int GroupsetIndex, float ResidualTol)
void chiLBSGroupsetSetWGDSA(int SolverIndex, int GroupsetIndex, int MaxIters, float ResTol, bool Verbose, char PETSCString)
void chiLBSGroupsetSetIterativeMethod(int SolverIndex, int GroupsetIndex, int IterativeMethod)
void chiLBSCreateGroupset(int SolverIndex)
void chiLBSGroupsetSetGMRESRestartIntvl(int SolverIndex, int GroupsetIndex, int Intvl)
Returns chiCreateProductQuadrature(int QuadratureType, varying values)

Groupsets segregate the code into pieces arranged by the number of groups it contains. A great deal of care must be taken with intergroupset transfer since the order in which the groupsets are executed determine what information will be available to them.

Function Documentation

◆ chiLBSCreateGroup()

◆ chiLBSCreateGroupset()

◆ chiLBSGroupsetAddGroups()

void chi_lua::chiLBSGroupsetAddGroups ( int  SolverIndex,
int  GroupsetIndex,
int  FromIndex,
int  ToIndex 
)

◆ chiLBSGroupsetSetAngleAggDiv()

void chi_lua::chiLBSGroupsetSetAngleAggDiv ( int  SolverIndex,
int  GroupsetIndex,
int  NumDiv 
)

Sets the angle aggregation divisions

Parameters
SolverIndexint Handle to the solver for which the group is to be created.
GroupsetIndexint Handle to the groupset to which the group is to be added.
NumDivint Number of divisions to use for the angle aggregation.

Note: by default polar aggregation will combine all polar angles in a hemisphere for a given azimuthal angleset. Therefore if there are 24 polar angles and 4 azimuthal angles the default polar aggregation will create 8 anglesets (2 per quadrant to allow top and bottom hemisphere) and each angleset will have the 12 polar angles associated with a hemisphere. When the number of divisions is greater than 1 then the polar angles will be split into divisions. For example if the number of divisions is 2 then more angleset will be created, this time having 6 polar angles per angleset.

_

Example:

Usage Examples:

test/framework/tutorials/tutorial_93_raytracing.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_1.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_2.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_2.lua

Definition at line 3777 of file lua_functions.c.

◆ chiLBSGroupsetSetAngleAggregationType()

void chi_lua::chiLBSGroupsetSetAngleAggregationType ( int  SolverIndex,
int  GroupsetIndex,
int  AggregationType 
)

Sets the the type of angle aggregation to use for this groupset.

Parameters
SolverIndexint Handle to the solver for which the group is to be created.
GroupsetIndexint Index to the groupset to which this function should apply
AggregationTypeint See AggregationType.

_

AggregationType

LBSGroupset.ANGLE_AGG_POLAR
Use Polar angle aggregation. This is the default.

LBSGroupset.ANGLE_AGG_SINGLE
Use Single angle aggregation.

LBSGroupset.ANGLE_AGG_AZIMUTHAL
Use Azimuthal angle aggregation.

Example:

chiLBSGroupsetSetAngleAggregationType(phys1,cur_gs,LBSGroupset.ANGLE_AGG_POLAR)
void chiLBSGroupsetSetAngleAggregationType(int SolverIndex, int GroupsetIndex, int AggregationType)

Usage Examples:

test/modules/LinearBoltzmannSolvers/Transport_Steady_Cyl/Transport2DCyl_1Monoenergetic.lua

Definition at line 3744 of file lua_functions.c.

◆ chiLBSGroupsetSetEnableSweepLog()

void chi_lua::chiLBSGroupsetSetEnableSweepLog ( int  SolverIndex,
int  GroupsetIndex,
bool  flag 
)

Enables or disables the printing of a sweep log.

Parameters
SolverIndexint Handle to the solver for which the group is to be created.
GroupsetIndexint Index to the groupset to which this function should apply
flagbool Flag indicating whether to print sweep log. Default false.

_

Example:

void chiLBSGroupsetSetEnableSweepLog(int SolverIndex, int GroupsetIndex, bool flag)

Definition at line 4006 of file lua_functions.c.

◆ chiLBSGroupsetSetGMRESRestartIntvl()

void chi_lua::chiLBSGroupsetSetGMRESRestartIntvl ( int  SolverIndex,
int  GroupsetIndex,
int  Intvl 
)

◆ chiLBSGroupsetSetGroupSubsets()

void chi_lua::chiLBSGroupsetSetGroupSubsets ( int  SolverIndex,
int  GroupsetIndex,
int  NumDiv 
)

Sets the number of group-subsets to use for groupset. Default 1.

Parameters
SolverIndexint Handle to the solver for which the group is to be created.
GroupsetIndexint Index to the groupset to which this function should apply
NumDivint Number of divisions of the groupset to use.

_

Example:

Usage Examples:

test/framework/tutorials/tutorial_93_raytracing.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_1.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_2.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_2.lua

Definition at line 3802 of file lua_functions.c.

◆ chiLBSGroupsetSetIterativeMethod()

void chi_lua::chiLBSGroupsetSetIterativeMethod ( int  SolverIndex,
int  GroupsetIndex,
int  IterativeMethod 
)

Sets the number of group-subsets to use for groupset. Default 1.

Parameters
SolverIndexint Handle to the solver for which the group is to be created.
GroupsetIndexint Index to the groupset to which this function should apply
IterativeMethodint Iteritve method identifier.

_

IterativeMethod

NPT_CLASSICRICHARDSON
Standard source iteration, without using PETSc.

NPT_CLASSICRICHARDSON_CYCLES
Standard source iteration, without using PETSc, with cyclic dependency convergence.

NPT_GMRES
Legacy Generalized Minimal Residual formulation for iterations.

NPT_GMRES_CYCLES
Legacy Generalized Minimal Residual formulation for iterations with cyclic dependency convergence.

KRYLOV_RICHARDSON
Richardson iteration.

KRYLOV_RICHARDSON_CYCLES
Richardson iteration with cyclic dependency convergence.

KRYLOV_GMRES
Generalized Minimal Residual method.

KRYLOV_GMRES_CYCLES
Generalized Minimal Residual method with cyclic dependency convergence.

KRYLOV_BICGSTAB
Biconjugate Gradient Stabilized method.

KRYLOV_BICGSTAB_CYCLES
Biconjugate Gradient Stabilized method with cyclic dependency convergence.

_

Notes on selecting iterative methods

The iterative methods NPT_CLASSICRICHARDSON, NPT_CLASSICRICHARDSON_CYCLES, NPT_GMRES and NPT_GMRES_CYCLES are considered legacy. The NPT_GMRES and NPT_GMRES_CYCLES are now considered deprecated with the inclusion of the generalized Krylov iteration method-class (which supports all the options prepended with KRYLOV_).

RICHARDSON is probably the least memory consuming but has the poorest convergence rate.

GMRES generally has the best convergence rate but it builds a basis comprising multiple solutions vectors, the amount of which is controlled via the gmres-restart parameter, which can dramatically increase memory consumption. GMRES restarts, i.e. the amount of iterations before the basis is destroyed and restarted, influences both memory consumptions and convergence behavior, e.g., lower restart numbers generally lowers memory consumption but increases the amount of iteration required to convergence.

The required memory and the computational time for one iteration with BiCGStab is constant, i.e., the time and memory requirements do not increase with the number of iterations as they do for restarted GMRES. BiCGStab uses approximately the same amount of memory as GMRES uses for two iterations. Therefore, BiCGStab typically uses less memory than GMRES. The convergence behavior of BiCGStab is often more irregular than that of GMRES. Intermediate residuals can even be orders of magnitude larger than the initial residual, which can affect the numerical accuracy as well as the rate of convergence. If the algorithm detects poor accuracy in the residual or the risk of stagnation, it restarts the iterations with the current solution as the initial guess. In contrast to GMRES, BiCGStab uses two matrix-vector multiplications each iteration (requiring two transport sweeps). Also, when using the left-preconditioned BiCGStab, an additional preconditioning step is required each iteration. That is, left-preconditioned BiCGStab requires a total of three preconditioning steps in each iteration. We generally apply Within-group Diffusion Synthetic Acceleration (WGDSA) and Two-Grid Acceleration (TGDSA) as left-preconditioners and therefore the total cost of these pre-conditioners will increase when using BiCGStab. Use BiCGStab when you are running problem with a high scattering order (i.e., L is large) because this will dramatically increase the GMRES basis.

Example:

chiLBSGroupsetSetIterativeMethod(phys1,cur_gs,NPT_CLASSICRICHARDSON)
chiLBSGroupsetSetIterativeMethod(phys1,cur_gs,NPT_GMRES)

Usage Examples:

test/framework/tutorials/tutorial_93_raytracing.lua
test/modules/LinearBoltzmannSolvers/Transport_Steady_Cyl/Transport2DCyl_1Monoenergetic.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_1.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_2.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_2.lua
test/modules/LinearBoltzmannSolvers/Transport_Keigen/c5g7/c5g7.lua

Definition at line 3905 of file lua_functions.c.

◆ chiLBSGroupsetSetMaxIterations()

void chi_lua::chiLBSGroupsetSetMaxIterations ( int  SolverIndex,
int  GroupsetIndex,
int  Numiter 
)

◆ chiLBSGroupsetSetQuadrature()

void chi_lua::chiLBSGroupsetSetQuadrature ( int  SolverIndex,
int  GroupsetIndex,
int  QuadratureIndex 
)

◆ chiLBSGroupsetSetResidualTolerance()

void chi_lua::chiLBSGroupsetSetResidualTolerance ( int  SolverIndex,
int  GroupsetIndex,
float  ResidualTol 
)

Sets the residual tolerance for the iterative method of the groupset.

Parameters
SolverIndexint Handle to the solver for which the group is to be created.
GroupsetIndexint Index to the groupset to which this function should apply
ResidualTolfloat residual tolerance (default 1.0e-6)

Note this tolerance also gets used for classic-richardson pointwise convergence tolerance.

_

Example:

Usage Examples:

test/framework/tutorials/tutorial_93_raytracing.lua
test/modules/LinearBoltzmannSolvers/Transport_Steady_Cyl/Transport2DCyl_1Monoenergetic.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_1.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_2.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_2.lua

Definition at line 3935 of file lua_functions.c.

◆ chiLBSGroupsetSetTGDSA()

void chi_lua::chiLBSGroupsetSetTGDSA ( int  SolverIndex,
int  GroupsetIndex,
int  MaxIters,
float  ResTol,
bool  Verbose,
char  PETSCString 
)

Sets the Two-Grid Diffusion Synthetic Acceleration parameters for this groupset. If this call is being made then it is assumed TGDSA is being applied.

Parameters
SolverIndexint Handle to the solver for which the group is to be created.
GroupsetIndexint Index to the groupset to which this function should apply
MaxItersint Maximum amount of iterations to use for TGDSA solvers. Default 30.
ResTolfloat Residual tolerance to use for the TGDSA solve.
Verbosebool Optional flag indicating verbose output of TGDSA. Default false.
PETSCStringchar Optional. Options string to be inserted during initialization.

_

Example:

petsc_options = " -pc_hypre_boomeramg_strong_threshold 0.8"
petsc_options = petsc_options .. " -pc_hypre_boomeramg_max_levels 25"
chiLBSGroupsetSetTGDSA(phys1,cur_gs,30,1.0e-4,false,petsc_options)

Usage Examples:

test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_1.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_2.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_2.lua

Definition at line 4082 of file lua_functions.c.

◆ chiLBSGroupsetSetWGDSA()

void chi_lua::chiLBSGroupsetSetWGDSA ( int  SolverIndex,
int  GroupsetIndex,
int  MaxIters,
float  ResTol,
bool  Verbose,
char  PETSCString 
)

Sets the Within-Group Diffusion Synthetic Acceleration parameters for this groupset. If this call is being made then it is assumed WGDSA is being applied.

Parameters
SolverIndexint Handle to the solver for which the group is to be created.
GroupsetIndexint Index to the groupset to which this function should apply
MaxItersint Maximum amount of iterations to use for WGDSA solvers. Default 30.
ResTolfloat Residual tolerance to use for the WGDSA solve.
Verbosebool Optional flag indicating verbose output of WGDSA. Default false.
PETSCStringchar Optional. Options string to be inserted during initialization.

_

Example:

petsc_options = " -pc_hypre_boomeramg_strong_threshold 0.8"
petsc_options = petsc_options .. " -pc_hypre_boomeramg_max_levels 25"
chiLBSGroupsetSetWGDSA(phys1,cur_gs,30,1.0e-4,false,petsc_options)

Usage Examples:

test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_1.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_3.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport2D_2.lua
test/modules/LinearBoltzmannSolvers/Transport_Transient/TransientTransport1D_2.lua

Definition at line 4044 of file lua_functions.c.