Chombo + EB + MF  3.2
LevelTGA.H
Go to the documentation of this file.
1 #ifdef CH_LANG_CC
2 /*
3  * _______ __
4  * / ___/ / ___ __ _ / / ___
5  * / /__/ _ \/ _ \/ V \/ _ \/ _ \
6  * \___/_//_/\___/_/_/_/_.__/\___/
7  * Please refer to Copyright.txt, in Chombo's root directory.
8  */
9 #endif
10 
11 #ifndef _LEVELTGA_H_
12 #define _LEVELTGA_H_
13 
14 #include <iostream>
15 #include <math.h>
16 #include "SPACE.H"
17 #include <stdlib.h>
18 #include <REAL.H>
19 #include <Box.H>
20 #include <DisjointBoxLayout.H>
21 #include <LevelData.H>
22 #include <FArrayBox.H>
23 #include "FluxBox.H"
24 #include "LevelFluxRegister.H"
25 #include <ProblemDomain.H>
26 #include "BaseLevelTGA.H"
27 #include "BaseLevelBackwardEuler.H"
28 #include "NamespaceHeader.H"
29 
30 /// Implements TGA algorithm to update solution to heat equation
31 /** The LevelTGA class implements the Runge-Kutta-based approach
32  to solving the heat equation due to Twizell, Gumel, and Arigu.
33 */
34 class LevelTGA : public BaseLevelTGA<LevelData<FArrayBox>, FluxBox, LevelFluxRegister>
35 {
36 
37 public:
38 
39  /// full constructor
40  /**
41  */
43  const Vector<int>& a_refRat,
44  const ProblemDomain& a_level0Domain,
47  :BaseLevelTGA<LevelData<FArrayBox>, FluxBox, LevelFluxRegister>(a_grids, a_refRat, a_level0Domain, a_opFact, a_solver)
48  {
49  ;
50  }
51 
52  /// destructor
53  virtual ~LevelTGA()
54  {
55  }
56 
57  /// computes time-centered diffusive term for explicit updating
58  /** compute time-centered L(phi) for use in subsequent update
59  operations. In this case, we do solve for phiNew, then
60  subtract source and old phi back out to get L(phi).
61  */
62  void computeDiffusion(LevelData<FArrayBox>& a_DiffusiveTerm,
63  LevelData<FArrayBox>& a_phiOld,
64  LevelData<FArrayBox>& a_src,
65  LevelFluxRegister* a_FineFluxRegPtr,
66  LevelFluxRegister* a_crseFluxRegPtr,
67  const LevelData<FArrayBox>* a_crsePhiOldPtr,
68  const LevelData<FArrayBox>* a_crsePhiNewPtr,
69  Real a_oldTime,
70  Real a_crseOldTime,
71  Real a_crseNewTime,
72  Real a_dt,
73  int a_level);
74 
75  /// do solve
76  /** update phi to a_phiNew
77  */
78  void updateSoln(LevelData<FArrayBox>& a_phiNew,
79  LevelData<FArrayBox>& a_phiOld,
80  LevelData<FArrayBox>& a_src,
81  LevelFluxRegister* a_fineFluxRegPtr,
82  LevelFluxRegister* a_crseFluxRegPtr,
83  const LevelData<FArrayBox>* a_crsePhiOldPtr,
84  const LevelData<FArrayBox>* a_crsePhiNewPtr,
85  Real a_oldTime,
86  Real a_crseOldTime,
87  Real a_crseNewTime,
88  Real a_dt,
89  int a_level,
90  bool a_zeroPhi = true,
91  int a_fluxStartComponent = 0);
92 
94  const DisjointBoxLayout& a_grids,
95  int a_level);
96 protected:
97 
98 };
99 
100 /// Implements BackwardEuler algorithm to update solution to heat equation
101 /** The LevelBackwardEuler class implements the Runge-Kutta-based approach
102  to solving the heat equation due to Twizell, Gumel, and Arigu.
103 */
104 class LevelBackwardEuler : public BaseLevelBackwardEuler<LevelData<FArrayBox>, FluxBox, LevelFluxRegister>
105 {
106 
107 public:
108 
109  /// full constructor
110  /**
111  */
113  const Vector<int>& a_refRat,
114  const ProblemDomain& a_level0Domain,
116  const RefCountedPtr<AMRMultiGrid<LevelData<FArrayBox> > >& a_solver)
117  :BaseLevelBackwardEuler<LevelData<FArrayBox>, FluxBox, LevelFluxRegister>(a_grids, a_refRat, a_level0Domain, a_opFact, a_solver)
118  {
119  ;
120  }
121 
122  /// destructor
124  {
125  }
126 
127  /// computes time-centered diffusive term for explicit updating
128  /** compute time-centered L(phi) for use in subsequent update
129  operations. In this case, we do solve for phiNew, then
130  subtract source and old phi back out to get L(phi).
131  */
132  void computeDiffusion(LevelData<FArrayBox>& a_DiffusiveTerm,
133  LevelData<FArrayBox>& a_phiOld,
134  LevelData<FArrayBox>& a_src,
135  LevelFluxRegister* a_FineFluxRegPtr,
136  LevelFluxRegister* a_crseFluxRegPtr,
137  const LevelData<FArrayBox>* a_crsePhiOldPtr,
138  const LevelData<FArrayBox>* a_crsePhiNewPtr,
139  Real a_oldTime,
140  Real a_crseOldTime,
141  Real a_crseNewTime,
142  Real a_dt,
143  int a_level);
144 
145  /// do solve
146  /** update phi to a_phiNew
147  */
148  void updateSoln(LevelData<FArrayBox>& a_phiNew,
149  LevelData<FArrayBox>& a_phiOld,
150  LevelData<FArrayBox>& a_src,
151  LevelFluxRegister* a_fineFluxRegPtr,
152  LevelFluxRegister* a_crseFluxRegPtr,
153  const LevelData<FArrayBox>* a_crsePhiOldPtr,
154  const LevelData<FArrayBox>* a_crsePhiNewPtr,
155  Real a_oldTime,
156  Real a_crseOldTime,
157  Real a_crseNewTime,
158  Real a_dt,
159  int a_level,
160  bool a_zeroPhi = true,
161  int a_fluxStartComponent = 0);
162 
164  const DisjointBoxLayout& a_grids,
165  int a_level);
166 protected:
167 
168 };
169 
170 #include "NamespaceFooter.H"
171 
172 #endif
A reference-counting handle class.
Definition: RefCountedPtr.H:173
virtual ~LevelTGA()
destructor
Definition: LevelTGA.H:53
A class to facilitate interaction with physical boundary conditions.
Definition: ProblemDomain.H:141
LevelBackwardEuler(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
Definition: LevelTGA.H:112
Definition: BaseLevelBackwardEuler.H:44
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
Definition: LevelTGA.H:42
void setSourceGhostCells(LevelData< FArrayBox > &a_src, const DisjointBoxLayout &a_grids, int a_level)
Definition: AMRMultiGrid.H:308
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
A FArrayBox-like container for face-centered fluxes.
Definition: FluxBox.H:22
Implements BackwardEuler algorithm to update solution to heat equation.
Definition: LevelTGA.H:104
double Real
Definition: REAL.H:33
Definition: BaseLevelTGA.H:44
A BoxLayout that has a concept of disjointedness.
Definition: DisjointBoxLayout.H:30
Implements TGA algorithm to update solution to heat equation.
Definition: LevelTGA.H:34
LevelFluxRegister-A class to encapsulate a levels worth of flux registers.
Definition: LevelFluxRegister.H:29
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
Definition: FArrayBox.H:45
Definition: AMRMultiGrid.H:233
virtual ~LevelBackwardEuler()
destructor
Definition: LevelTGA.H:123