|
| BaseLevelHeatSolver (const Vector< DisjointBoxLayout > &a_grids, const Vector< int > &a_refRat, const ProblemDomain &a_level0Domain, RefCountedPtr< AMRLevelOpFactory< LevelDataType > > &a_opFact, const RefCountedPtr< AMRMultiGrid< LevelDataType > > &a_solver) |
|
virtual | ~BaseLevelHeatSolver () |
| Destructor, called after destructors of BaseLevelHeatSolver subclasses. More...
|
|
virtual void | updateSoln (LevelDataType &a_phiNew, LevelDataType &a_phiOld, LevelDataType &a_src, LevelData< FluxDataType > &a_flux, FluxRegisterType *a_fineFluxRegPtr, FluxRegisterType *a_crseFluxRegPtr, const LevelDataType *a_crsePhiOldPtr, const LevelDataType *a_crsePhiNewPtr, Real a_oldTime, Real a_crseOldTime, Real a_crseNewTime, Real a_dt, int a_level, bool a_zeroPhi=true, bool a_rhsAlreadyKappaWeighted=false, int a_fluxStartComponent=0)=0 |
|
virtual void | updateSoln (LevelDataType &a_phiNew, LevelDataType &a_phiOld, LevelDataType &a_src, FluxRegisterType *a_fineFluxRegPtr, FluxRegisterType *a_crseFluxRegPtr, const LevelDataType *a_crsePhiOldPtr, const LevelDataType *a_crsePhiNewPtr, Real a_oldTime, Real a_crseOldTime, Real a_crseNewTime, Real a_dt, int a_level, bool a_zeroPhi=true, bool a_rhsAlreadyKappaWeighted=false, int a_fluxStartComponent=0) |
|
virtual void | computeDiffusion (LevelDataType &a_diffusiveTerm, LevelDataType &a_phiOld, LevelDataType &a_src, LevelData< FluxDataType > &a_flux, FluxRegisterType *a_fineFluxRegPtr, FluxRegisterType *a_crseFluxRegPtr, const LevelDataType *a_crsePhiOldPtr, const LevelDataType *a_crsePhiNewPtr, Real a_oldTime, Real a_crseOldTime, Real a_crseNewTime, Real a_dt, int a_level, bool a_zeroPhi=true, bool a_rhsAlreadyKappaWeighted=false) |
|
virtual void | applyOperator (LevelDataType &a_ans, const LevelDataType &a_phi, const LevelDataType *a_phiC, int a_level, Real a_alpha, Real a_beta, bool a_applyBC) |
|
void | applyHelm (LevelDataType &a_ans, const LevelDataType &a_phi, const LevelDataType *a_phiC, int a_level, Real a_mu, Real a_dt, bool a_homogeneousBC) |
|
void | incrementFlux (LevelData< FluxDataType > &a_diffusiveFlux, LevelDataType &a_phi, int a_level, Real a_mu, Real a_dt, Real a_sign, bool a_setToZero) |
|
void | solveHelm (LevelDataType &a_phi, LevelDataType &a_phiC, LevelDataType &a_rhs, int a_level, Real a_mu, Real a_dt, bool a_zeroPhi=true) |
|
void | resetSolverAlphaAndBeta (const Real &a_alpha, const Real &a_beta) |
|
virtual LevelTGAHelmOp< LevelDataType, FluxDataType > * | newOp (const ProblemDomain &a_indexSpace, RefCountedPtr< AMRLevelOpFactory< LevelDataType > > &a_opFact) |
|
int | size () const |
| Returns the number of grid levels on which this integrator operates. More...
|
|
template<class LevelDataType, class FluxDataType, class FluxRegisterType>
class BaseLevelHeatSolver< LevelDataType, FluxDataType, FluxRegisterType >
This base class implements the 2nd-order implicit L0-stable time integration algorithm developed by Twizell, Gumel, and Arigu for solving elliptic equations. It relies upon linear algebraic operations defined in the underlying Helmholtz operators.
- Template Parameters
-
LevelDataType | The type used to store data at a grid level. This is usually LevelData<T>, where T is some cell-centered FArrayBox type. |
FluxDataType | The type used to store flux data at a grid level. This is usually an array box clas that stores fluxes. |
FluxRegisterType | The type used to store flux register data for interactions between different grid levels. This is usually a flux register class. |
template<class LevelDataType, class FluxDataType, class FluxRegisterType>
virtual void BaseLevelHeatSolver< LevelDataType, FluxDataType, FluxRegisterType >::updateSoln |
( |
LevelDataType & |
a_phiNew, |
|
|
LevelDataType & |
a_phiOld, |
|
|
LevelDataType & |
a_src, |
|
|
LevelData< FluxDataType > & |
a_flux, |
|
|
FluxRegisterType * |
a_fineFluxRegPtr, |
|
|
FluxRegisterType * |
a_crseFluxRegPtr, |
|
|
const LevelDataType * |
a_crsePhiOldPtr, |
|
|
const LevelDataType * |
a_crsePhiNewPtr, |
|
|
Real |
a_oldTime, |
|
|
Real |
a_crseOldTime, |
|
|
Real |
a_crseNewTime, |
|
|
Real |
a_dt, |
|
|
int |
a_level, |
|
|
bool |
a_zeroPhi = true , |
|
|
bool |
a_rhsAlreadyKappaWeighted = false , |
|
|
int |
a_fluxStartComponent = 0 |
|
) |
| |
|
pure virtual |
Integrates the helmholtz equation represented by this object, placing the new solution in a_phiNew.
- Parameters
-
a_phiNew | The new solution (the value of phi at time n + 1) will be stored here. |
a_phiOld | The old solution (the value of phi at time n). |
a_src | The source term on the right hand side of the Helmholtz equation. |
a_flux | This will store the flux computed at the current grid level during the solution of the Helmholtz equation. |
a_fineFluxRegPtr | A pointer to the flux register representing the finer grid level adjacent to this one, or NULL if there is no finer grid level. |
a_crseFluxRegPtr | A pointer to the flux register representing the coarser grid level adjacent to this one, or NULL if there is no coarser grid level. |
a_oldTime | The time at the beginning of the integration step at the current grid level. |
a_crseOldTime | The time at the beginning of the integration step at the coarser adjacent grid level. This parameter is ignored if there is no coarser grid level. |
a_crseNewTime | The time at the end of the integration step at the coarser adjacent grid level. This parameter is ignored if there is no coarser grid level. |
a_dt | The size of the integration step at the current grid level. |
a_level | The current grid level. |
a_zeroPhi | If set to true, a_phiNew will be set to zero before the integration takes place. Otherwise, a_phiNew is assumed to be an initial estimate for the solution in the iterative linear solve. |
a_fluxStartComponent | An index identifying the component at which flux data begins within a_fineFluxRegPtr and a_crseFluxRegPtr. |
Implemented in BaseLevelTGA< LevelDataType, FluxDataType, FluxRegisterType >, BaseLevelTGA< LevelData< EBCellFAB >, EBFluxFAB, EBFluxRegister >, BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >, BaseLevelBackwardEuler< LevelDataType, FluxDataType, FluxRegisterType >, BaseLevelCrankNicolson< LevelDataType, FluxDataType, FluxRegisterType >, BaseLevelBackwardEuler< LevelData< EBCellFAB >, EBFluxFAB, EBFluxRegister >, BaseLevelBackwardEuler< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >, and BaseLevelCrankNicolson< LevelData< EBCellFAB >, EBFluxFAB, EBFluxRegister >.
Referenced by BaseLevelHeatSolver< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::computeDiffusion(), and BaseLevelHeatSolver< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::~BaseLevelHeatSolver().
template<class LevelDataType, class FluxDataType, class FluxRegisterType>
virtual void BaseLevelHeatSolver< LevelDataType, FluxDataType, FluxRegisterType >::computeDiffusion |
( |
LevelDataType & |
a_diffusiveTerm, |
|
|
LevelDataType & |
a_phiOld, |
|
|
LevelDataType & |
a_src, |
|
|
LevelData< FluxDataType > & |
a_flux, |
|
|
FluxRegisterType * |
a_fineFluxRegPtr, |
|
|
FluxRegisterType * |
a_crseFluxRegPtr, |
|
|
const LevelDataType * |
a_crsePhiOldPtr, |
|
|
const LevelDataType * |
a_crsePhiNewPtr, |
|
|
Real |
a_oldTime, |
|
|
Real |
a_crseOldTime, |
|
|
Real |
a_crseNewTime, |
|
|
Real |
a_dt, |
|
|
int |
a_level, |
|
|
bool |
a_zeroPhi = true , |
|
|
bool |
a_rhsAlreadyKappaWeighted = false |
|
) |
| |
|
inlinevirtual |
Computes the time-centered diffusion term L(phi). This can be used to find contributions to the solution from diffusion. The diffusion term is computed by computing a finite difference approximation for d phi/dt using the updated and original values of phi and the time step. Most of the arguments given here are passed along to updateSoln and retain their significance therein.
- Parameters
-
a_diffusiveTerm | The diffusion term L(phi) will be stored here. |
a_phiOld | The old solution (the value of phi at time n). |
a_src | The source term on the right hand side of the Helmholtz equation (used to fine the new value of phi). |
a_fineFluxRegPtr | A pointer to the flux register representing the finer grid level adjacent to this one, or NULL if there is no finer grid level. |
a_crseFluxRegPtr | A pointer to the flux register representing the coarser grid level adjacent to this one, or NULL if there is no coarser grid level. |
a_crsePhiOldTime | A pointer to the value of phi at the beginning of the integration step at the coarser grid level adjacent to this one, or NULL if there is no coarser grid level. |
a_crsePhiNewTime | A pointer to the value of phi at the end of the integration step at the coarser grid level adjacent to this one, or NULL if there is no coarser grid level. |
a_oldTime | The time at the beginning of the integration step at the current grid level. |
a_crseOldTime | The time at the beginning of the integration step at the coarser adjacent grid level. This parameter is ignored if there is no coarser grid level. |
a_crseNewTime | The time at the end of the integration step at the coarser adjacent grid level. This parameter is ignored if there is no coarser grid level. |
a_dt | The size of the integration step at the current grid level. |
a_level | The current grid level. |
a_zeroPhi | If set to true, the new value of phi will be set to zero before the integration takes place. Otherwise, it will be set to the value in a_diffusiveTerm. |
Reimplemented in BaseLevelTGA< LevelDataType, FluxDataType, FluxRegisterType >, BaseLevelTGA< LevelData< EBCellFAB >, EBFluxFAB, EBFluxRegister >, and BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >.
template<class LevelDataType, class FluxDataType, class FluxRegisterType>
void BaseLevelHeatSolver< LevelDataType, FluxDataType, FluxRegisterType >::applyHelm |
( |
LevelDataType & |
a_ans, |
|
|
const LevelDataType & |
a_phi, |
|
|
const LevelDataType * |
a_phiC, |
|
|
int |
a_level, |
|
|
Real |
a_mu, |
|
|
Real |
a_dt, |
|
|
bool |
a_homogeneousBC |
|
) |
| |
|
inline |
Applies the Helmholtz operator to the solution a_phi at the given grid level. This will set a_ans to (I + a_mu * a_dt * L(a_phi).
- Parameters
-
a_ans | This will store the value of the Helmholtz operator applied to a_phi. |
a_phi | The value of the solution to which the Helmholtz operator is applied. |
a_phiC | A pointer to the value of the solution at the coarser adjacent grid level, or NULL if there is no coarser grid level. |
a_level | The grid level at which the Helmholtz operator is applied. |
a_mu | A number between 0 and 1 that defines the time within the integration step at which the Helmholtz operator is evaluated. |
a_dt | The size of the integration step at the current grid level. |
a_applyBC | If set to true, inhomogeneous boundary conditions will be applied to the Helmholtz operator (both at the domain boundary and at the coarse-fine grid interfaces). If set to false, homogeneous boundary conditions will applied to the Helmholtz operator at these locations. protected because not flexible (alpha == 1) |
Referenced by BaseLevelCrankNicolson< LevelData< EBCellFAB >, EBFluxFAB, EBFluxRegister >::updateSoln(), BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSolnWithTimeDependentOp(), and BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSolnWithTimeIndependentOp().
template<class LevelDataType, class FluxDataType, class FluxRegisterType>
void BaseLevelHeatSolver< LevelDataType, FluxDataType, FluxRegisterType >::incrementFlux |
( |
LevelData< FluxDataType > & |
a_diffusiveFlux, |
|
|
LevelDataType & |
a_phi, |
|
|
int |
a_level, |
|
|
Real |
a_mu, |
|
|
Real |
a_dt, |
|
|
Real |
a_sign, |
|
|
bool |
a_setToZero |
|
) |
| |
|
inline |
Adds flux contributions from the Helmholtz operator at the current grid level.
- Parameters
-
a_diffusiveFlux | The flux to which Helmholtz operator flux contributions will be accumulated. |
a_phi | The solution at the current grid level to which the Helmholtz operator is applied in order to evaluate the flux. |
a_level | The grid level at which the flux is computed. |
a_mu | A number between 0 and 1 that defines the time within the integration step at which the Helmholtz operator is evaluated. |
a_dt | The size of the integration step at the current grid level. |
a_sign | A factor applied to the derivative term in the Helmholtz equation. This allows the sign to be made positive or negative as is necessary for the flux calculation. |
a_setToZero | If true, a_diffusiveFlux will be set to zero before the fluxes are accumulated to it. Otherwise, fluxes will be added to the existing value. |
Referenced by BaseLevelBackwardEuler< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSoln(), BaseLevelCrankNicolson< LevelData< EBCellFAB >, EBFluxFAB, EBFluxRegister >::updateSoln(), BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSolnWithTimeDependentOp(), and BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSolnWithTimeIndependentOp().
template<class LevelDataType, class FluxDataType, class FluxRegisterType>
void BaseLevelHeatSolver< LevelDataType, FluxDataType, FluxRegisterType >::solveHelm |
( |
LevelDataType & |
a_phi, |
|
|
LevelDataType & |
a_phiC, |
|
|
LevelDataType & |
a_rhs, |
|
|
int |
a_level, |
|
|
Real |
a_mu, |
|
|
Real |
a_dt, |
|
|
bool |
a_zeroPhi = true |
|
) |
| |
|
inline |
Solves the Helmholtz equation (I - a_mu * a_dt * L(a_phi) = a_rhs for a_phi. Here it is assumed that a solution a_phiC exists on a coarser grid level.
- Parameters
-
a_phi | This will store the solution to the Helmholtz equation on the given grid level. |
a_phiC | The value of the solution given at the coarser adjacent grid level. If there is no coarser grid level, this parameter is ignored. |
a_rhs | The right hand side of the Helmholtz equation at the given grid level. |
a_level | The grid level at which the Helmholtz equation is solved. |
a_mu | A number between 0 and 1 that defines the time within the integration step at which the Helmholtz operator is evaluated. |
a_dt | The size of the integration step at the current grid level. |
a_zeroPhi | If set to true, a_phi will be set to zero before the (iterative) solve takes place. Otherwise, a_phi is assumed to be an initial estimate for the solution. |
Referenced by BaseLevelBackwardEuler< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSoln(), BaseLevelCrankNicolson< LevelData< EBCellFAB >, EBFluxFAB, EBFluxRegister >::updateSoln(), BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSolnWithTimeDependentOp(), and BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSolnWithTimeIndependentOp().
template<class LevelDataType, class FluxDataType, class FluxRegisterType>
void BaseLevelHeatSolver< LevelDataType, FluxDataType, FluxRegisterType >::timeInterp |
( |
LevelDataType & |
a_data, |
|
|
const LevelDataType & |
a_oldData, |
|
|
const LevelDataType & |
a_newData, |
|
|
Real |
a_time, |
|
|
Real |
a_oldTime, |
|
|
Real |
a_newTime, |
|
|
int |
a_level |
|
) |
| |
|
inlineprotected |
Interpolates a given quantity linearly in time using its beginning- and end-of-step values and placing the result in a_data.
- Parameters
-
a_oldData | The value of the quantity at the beginning of the time step at the grid level identified by a_level. |
a_newData | The value of the quantity at the end of the time step at the grid level identified by a_level. |
a_time | The time at which the quantity is to be interpolated within the integration step. |
a_oldTime | The beginning of the integration step at the grid level identified by a_level. |
a_newTime | The end of the integration step at the grid level identified by a_level. |
a_level | An index identifying the grid level at which the interpolation takes place. |
Referenced by BaseLevelBackwardEuler< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSoln(), BaseLevelCrankNicolson< LevelData< EBCellFAB >, EBFluxFAB, EBFluxRegister >::updateSoln(), BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSolnWithTimeDependentOp(), and BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >::updateSolnWithTimeIndependentOp().