Chombo + EB  3.2
MultilevelLinearOp.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 // BVS, June 18, 2003
12 
13 // We can assume that template class T has null construction.
14 
15 #ifndef _MULTILEVELLINEAROP_H_
16 #define _MULTILEVELLINEAROP_H_
17 
18 #include <cmath>
19 
20 #include "RefCountedPtr.H"
21 #include "Vector.H"
22 #include "RealVect.H"
23 #include "LevelData.H"
24 #include "LinearSolver.H"
25 #include "AMRMultiGrid.H"
26 #include "BiCGStabSolver.H"
27 
28 #include "NamespaceHeader.H"
29 
30 ///
31 /**
32  Operator class for Linear solver for solving L(phi) = rhs for the case
33  where <T> is a multilevel (AMR) data structure; as such, T is actually
34  a Vector<T>. The MultilevelLinearOp is also
35  defined with a Vector<AMRLevelOp> which applies the operator in an AMR
36  way. For now, this is hardwired for T = LevelData<FArrayBox>, but we can
37  hopefully get back redo this once I get the kinks worked out...
38  */
39 template <class T>
40 class MultilevelLinearOp : public LinearOp< Vector<LevelData<T>* > >
41 {
42 
43 public:
44  /// default constructor (sets defaults for preconditioner)
46 
47  /// define function -- opFactory should be able to define all required levels
48  /**
49  */
50  virtual void define(const Vector<DisjointBoxLayout>& a_vectGrids,
51  const Vector<int>& a_refRatios,
52  const Vector<ProblemDomain>& a_domains,
53  const Vector<RealVect>& a_vectDx,
55  int a_lBase);
56 
57  ///
58  virtual ~MultilevelLinearOp();
59 
60  /// compute residual using AMRLevelOps AMRResidual functionality
61  /**
62  a_lhs = L(a_phi) - a_rhs.
63  If a_homogeneous is true, evaluate the operator using
64  homogeneous boundary conditions.
65  */
66  virtual void residual( Vector<LevelData<T>* >& a_lhs,
67  const Vector<LevelData<T>* >& a_phi,
68  const Vector<LevelData<T>* >& a_rhs,
69  bool a_homogeneous = false);
70 
71  /// Apply preconditioner
72  /**
73  Given the current state of the residual the correction,
74  apply your preconditioner to a_cor.
75  */
76  virtual void preCond(Vector<LevelData<T>* >& a_cor,
77  const Vector<LevelData<T>* >& a_residual);
78 
79  ///
80  /**
81  In the context of solving L(phi) = rhs, set a_lhs = L(a_phi).
82  If a_homogeneous is true,
83  evaluate the operator using homogeneous boundary conditions.
84  */
85  virtual void applyOp(Vector<LevelData<T>* >& a_lhs,
86  const Vector<LevelData<T>* >& a_phi,
87  bool a_homogeneous = false);
88 
89  ///
90  /**
91  Create data holder a_lhs that mirrors a_rhs. You do not need
92  to copy the data of a_rhs, just make a holder the same size.
93  */
94  virtual void create(Vector<LevelData<T>* >& a_lhs,
95  const Vector<LevelData<T>* >& a_rhs);
96 
97  ///
98  /**
99  Clean up data holder before it goes out of scope. This is necessary
100  because create calls new.
101  */
102  virtual void clear(Vector<LevelData<T>* >& a_lhs);
103 
104  ///
105  /**
106  Set a_lhs equal to a_rhs.
107  */
108  virtual void assign(Vector<LevelData<T>* >& a_lhs,
109  const Vector<LevelData<T>* >& a_rhs);
110 
111  ///
112  /**
113  Compute and return the dot product of a_1 and a_2. In most
114  contexts, this means return the sum over all data points of a_1*a_2.
115  */
116  virtual Real dotProduct(const Vector<LevelData<T>* >& a_1,
117  const Vector<LevelData<T>* >& a_2);
118 
119  ///
120  /**
121  Increment by scaled amount (a_lhs += a_scale*a_x).
122  */
123  virtual void incr (Vector<LevelData<T>* >& a_lhs,
124  const Vector<LevelData<T>* >& a_x,
125  Real a_scale);
126 
127  ///
128  /**
129  Set input to a scaled sum (a_lhs = a_a*a_x + a_b*a_y).
130  */
131  virtual void axby(Vector<LevelData<T>* >& a_lhs,
132  const Vector<LevelData<T>* >& a_x,
133  const Vector<LevelData<T>* >& a_y,
134  Real a_a,
135  Real a_b);
136 
137  ///
138  /**
139  Multiply the input by a given scale (a_lhs *= a_scale).
140  */
141  virtual void scale(Vector<LevelData<T>* >& a_lhs,
142  const Real& a_scale);
143 
144  ///
145  /**
146  Return the norm of a_rhs.
147  a_ord == 0 max norm, a_ord == 1 sum(abs(a_rhs)), else, L(a_ord) norm.
148  */
149  virtual Real norm(const Vector<LevelData<T>* >& a_rhs,
150  int a_ord);
151 
152  ///
153  /**
154  Set a_lhs to zero.
155  */
156  virtual void setToZero(Vector<LevelData<T>* >& a_lhs);
157 
158  virtual void write(const Vector<LevelData<T>* >* a_data,
159  const char* a_filename);
160 
161  /// if true, preCond(...) function uses multigrid v-cycle(s)
163 
164  /// number of mg v-cycles for each preCond call (if using mg preconditioner)
166 
167  /// number of smoothings in the multigrid preconditioner
169 
170  /// max depth for multigrid preconditioner (if using mg preconditioner)
172 
173  ///
175 
176  ///
178 
179  ///
181 
182  // cell spacing
184 
185  /// vector of level operators
187 
188  ///
189  int m_lBase;
190 
191  /// Solver which is used by the preCond function
192  /**
193  if m_use_multigrid_preconditioner is true, this solver is
194  used to apply a AMR multigrid v-cycle as the preconditioner
195  */
197 
198  /// bottom solver for m_preCondSolver
200 
201  /// has the preconditioner solver been initialized?
203 };
204 
205 #include "NamespaceFooter.H"
206 
207 #include "MultilevelLinearOpI.H"
208 #endif /*_MULTILEVELLINEAROP_H_*/
LinearSolver< LevelData< T > > * m_precondBottomSolverPtr
bottom solver for m_preCondSolver
Definition: MultilevelLinearOp.H:199
virtual void incr(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_x, Real a_scale)
Definition: MultilevelLinearOpI.H:525
Vector< int > m_refRatios
Definition: MultilevelLinearOp.H:177
Definition: LinearSolver.H:28
virtual Real norm(const Vector< LevelData< T > * > &a_rhs, int a_ord)
Definition: MultilevelLinearOpI.H:588
A reference-counting handle class.
Definition: RefCountedPtr.H:173
AMRMultiGrid< LevelData< T > > m_preCondSolver
Solver which is used by the preCond function.
Definition: MultilevelLinearOp.H:196
virtual void write(const Vector< LevelData< T > * > *a_data, const char *a_filename)
Definition: MultilevelLinearOpI.H:670
bool m_isPrecondSolverInitialized
has the preconditioner solver been initialized?
Definition: MultilevelLinearOp.H:202
Vector< ProblemDomain > m_domains
Definition: MultilevelLinearOp.H:180
Definition: AMRMultiGrid.H:308
virtual void axby(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_x, const Vector< LevelData< T > * > &a_y, Real a_a, Real a_b)
Definition: MultilevelLinearOpI.H:544
virtual ~MultilevelLinearOp()
Definition: MultilevelLinearOpI.H:118
virtual void setToZero(Vector< LevelData< T > * > &a_lhs)
Definition: MultilevelLinearOpI.H:655
MultilevelLinearOp()
default constructor (sets defaults for preconditioner)
Definition: MultilevelLinearOpI.H:21
virtual void scale(Vector< LevelData< T > * > &a_lhs, const Real &a_scale)
Definition: MultilevelLinearOpI.H:568
int m_num_mg_smooth
number of smoothings in the multigrid preconditioner
Definition: MultilevelLinearOp.H:168
new code
Definition: BoxLayoutData.H:170
int m_lBase
Definition: MultilevelLinearOp.H:189
int m_preCondSolverDepth
max depth for multigrid preconditioner (if using mg preconditioner)
Definition: MultilevelLinearOp.H:171
double Real
Definition: REAL.H:33
Vector< RefCountedPtr< AMRLevelOp< LevelData< T > > > > m_vectOperators
vector of level operators
Definition: MultilevelLinearOp.H:186
int m_num_mg_iterations
number of mg v-cycles for each preCond call (if using mg preconditioner)
Definition: MultilevelLinearOp.H:165
virtual void clear(Vector< LevelData< T > * > &a_lhs)
Definition: MultilevelLinearOpI.H:402
virtual void assign(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_rhs)
Definition: MultilevelLinearOpI.H:422
virtual void applyOp(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_phi, bool a_homogeneous=false)
Definition: MultilevelLinearOpI.H:319
virtual void residual(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_phi, const Vector< LevelData< T > * > &a_rhs, bool a_homogeneous=false)
compute residual using AMRLevelOps AMRResidual functionality
Definition: MultilevelLinearOpI.H:137
Definition: LinearSolver.H:156
virtual void preCond(Vector< LevelData< T > * > &a_cor, const Vector< LevelData< T > * > &a_residual)
Apply preconditioner.
Definition: MultilevelLinearOpI.H:217
virtual void define(const Vector< DisjointBoxLayout > &a_vectGrids, const Vector< int > &a_refRatios, const Vector< ProblemDomain > &a_domains, const Vector< RealVect > &a_vectDx, RefCountedPtr< AMRLevelOpFactory< LevelData< T > > > &a_opFactory, int a_lBase)
define function – opFactory should be able to define all required levels
Definition: MultilevelLinearOpI.H:35
Definition: MultilevelLinearOp.H:40
virtual Real dotProduct(const Vector< LevelData< T > * > &a_1, const Vector< LevelData< T > * > &a_2)
Definition: MultilevelLinearOpI.H:445
Definition: AMRMultiGrid.H:233
virtual void create(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_rhs)
Definition: MultilevelLinearOpI.H:379
bool m_use_multigrid_preconditioner
if true, preCond(...) function uses multigrid v-cycle(s)
Definition: MultilevelLinearOp.H:162
Vector< DisjointBoxLayout > m_vectGrids
Definition: MultilevelLinearOp.H:174
Vector< RealVect > m_vectDx
Definition: MultilevelLinearOp.H:183