#include <LevelTGA.H>

The LevelTGA class implements the Runge-Kutta-based approach to solving the heat equation due to Twizell, Gumel, and Arigu.
Public Member Functions | |
| LevelTGA (const Vector< DisjointBoxLayout > &a_grids, const Vector< int > &a_refRat, const ProblemDomain &a_level0Domain, RefCountedPtr< AMRLevelOpFactory< LevelData< FArrayBox > > > &a_opFact, const RefCountedPtr< AMRMultiGrid< LevelData< FArrayBox > > > &a_solver) | |
| full constructor | |
| virtual | ~LevelTGA () |
| destructor | |
| void | computeDiffusion (LevelData< FArrayBox > &a_DiffusiveTerm, LevelData< FArrayBox > &a_phiOld, LevelData< FArrayBox > &a_src, LevelFluxRegister *a_FineFluxRegPtr, LevelFluxRegister *a_crseFluxRegPtr, const LevelData< FArrayBox > *a_crsePhiOldPtr, const LevelData< FArrayBox > *a_crsePhiNewPtr, Real a_oldTime, Real a_crseOldTime, Real a_crseNewTime, Real a_dt, int a_level) |
| computes time-centered diffusive term for explicit updating | |
| void | updateSoln (LevelData< FArrayBox > &a_phiNew, LevelData< FArrayBox > &a_phiOld, LevelData< FArrayBox > &a_src, LevelFluxRegister *a_fineFluxRegPtr, LevelFluxRegister *a_crseFluxRegPtr, const LevelData< FArrayBox > *a_crsePhiOldPtr, const LevelData< FArrayBox > *a_crsePhiNewPtr, Real a_oldTime, Real a_crseOldTime, Real a_crseNewTime, Real a_dt, int a_level, bool a_zeroPhi=true, int a_fluxStartComponent=0) |
| do solve | |
| void | setSourceGhostCells (LevelData< FArrayBox > &a_src, const DisjointBoxLayout &a_grids, int a_level) |
| LevelTGA::LevelTGA | ( | const Vector< DisjointBoxLayout > & | a_grids, | |
| const Vector< int > & | a_refRat, | |||
| const ProblemDomain & | a_level0Domain, | |||
| RefCountedPtr< AMRLevelOpFactory< LevelData< FArrayBox > > > & | a_opFact, | |||
| const RefCountedPtr< AMRMultiGrid< LevelData< FArrayBox > > > & | a_solver | |||
| ) | [inline] |
full constructor
| virtual LevelTGA::~LevelTGA | ( | ) | [inline, virtual] |
destructor
| void LevelTGA::computeDiffusion | ( | LevelData< FArrayBox > & | a_DiffusiveTerm, | |
| LevelData< FArrayBox > & | a_phiOld, | |||
| LevelData< FArrayBox > & | a_src, | |||
| LevelFluxRegister * | a_FineFluxRegPtr, | |||
| LevelFluxRegister * | a_crseFluxRegPtr, | |||
| const LevelData< FArrayBox > * | a_crsePhiOldPtr, | |||
| const LevelData< FArrayBox > * | a_crsePhiNewPtr, | |||
| Real | a_oldTime, | |||
| Real | a_crseOldTime, | |||
| Real | a_crseNewTime, | |||
| Real | a_dt, | |||
| int | a_level | |||
| ) |
computes time-centered diffusive term for explicit updating
compute time-centered L(phi) for use in subsequent update operations. In this case, we do solve for phiNew, then subtract source and old phi back out to get L(phi).
| void LevelTGA::updateSoln | ( | LevelData< FArrayBox > & | a_phiNew, | |
| LevelData< FArrayBox > & | a_phiOld, | |||
| LevelData< FArrayBox > & | a_src, | |||
| LevelFluxRegister * | a_fineFluxRegPtr, | |||
| LevelFluxRegister * | a_crseFluxRegPtr, | |||
| const LevelData< FArrayBox > * | a_crsePhiOldPtr, | |||
| const LevelData< FArrayBox > * | a_crsePhiNewPtr, | |||
| Real | a_oldTime, | |||
| Real | a_crseOldTime, | |||
| Real | a_crseNewTime, | |||
| Real | a_dt, | |||
| int | a_level, | |||
| bool | a_zeroPhi = true, |
|||
| int | a_fluxStartComponent = 0 | |||
| ) |
do solve
update phi to a_phiNew
| void LevelTGA::setSourceGhostCells | ( | LevelData< FArrayBox > & | a_src, | |
| const DisjointBoxLayout & | a_dbl, | |||
| int | a_level | |||
| ) | [virtual] |
Sets the value of the source term in the Helmholtz equation on ghost cells. This method should be overridden by subclasses of BaseLevelTGA.
| a_src | The source term in the Helmholtz equation whose ghost cell values are to be set. | |
| a_dbl | The disjoint box layout that indexes a_src. | |
| a_level | The grid level at which the ghost cell values are to be set. |
Implements BaseLevelTGA< LevelData< FArrayBox >, FluxBox, LevelFluxRegister >.
1.5.5