Chombo + EB  3.0
MultilevelLinearOpI.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 _MULTILEVELLINEAROPI_H_
12 #define _MULTILEVELLINEAROPI_H_
13 
14 #include "AMRIO.H"
15 
16 #include "NamespaceHeader.H"
17 
18 // set a couple of defaults here (do before define so that
19 // they can be overridden before calling define)
20 template<class T>
22 {
23  // set some default values
24  m_use_multigrid_preconditioner = true;
25  m_num_mg_iterations = 1;
26  m_num_mg_smooth = 4;
27  // default is to coarsen all the way down
28  m_preCondSolverDepth = -1;
29  m_precondBottomSolverPtr = NULL;
30 }
31 
32 /// define function
33 template<class T>
34 void
36  const Vector<int>& a_refRatios,
37  const Vector<ProblemDomain>& a_domains,
38  const Vector<RealVect>& a_vectDx,
40  int a_lBase)
41 {
42  CH_TIME("MultilevelLinearOp::define");
43 
44  int numLevels = a_vectGrids.size();
45 
46  // couple of sanity checks...
47  CH_assert (a_lBase >= 0);
48  CH_assert (a_lBase < numLevels);
49  CH_assert ( a_domains.size() == numLevels);
50  // since you technically only need numLevels -1 refinement ratios...
51  CH_assert (a_refRatios.size() >= numLevels -1);
52 
53 
54  m_lBase = a_lBase;
55 
56  m_vectGrids = a_vectGrids;
57  m_refRatios = a_refRatios;
58  m_domains = a_domains;
59  m_vectDx = a_vectDx;
60  m_vectOperators.resize(numLevels);
61 
62  // need to at least define the levels to lBase -1 (for coarser BC's)
63  int coarsestLevel = Max(0, m_lBase-1);
64  for (int level=coarsestLevel; level<numLevels; level++)
65  {
66  RefCountedPtr<AMRLevelOp<LevelData<T> > > op = RefCountedPtr<AMRLevelOp<LevelData<T> > >(a_opFactory->AMRnewOp(m_domains[level]));
67  m_vectOperators[level] = op;
68  }
69 
70  // finally, define AMRMultigrid preconditioner if required
71  if (m_use_multigrid_preconditioner)
72  {
73  // preconditioner requires a bottom smoother
75 
76  m_precondBottomSolverPtr = newPtr;
77  int bottomSolverVerbosity = 1;
78 
79  newPtr->m_verbosity = bottomSolverVerbosity;
80 
81  if (m_preCondSolverDepth >= 0)
82  {
83  m_preCondSolver.m_maxDepth = m_preCondSolverDepth;
84  }
85  m_preCondSolver.define(m_domains[0],
86  *a_opFactory,
87  m_precondBottomSolverPtr,
88  numLevels);
89 
90  // most of these have no real meaning since we're only doina a
91  // single V-cycle, rather than a full solve
92  Real solverEps = 1.0e-7;
93  int maxIterations = 1;
94  Real hang = 1.0e-7;
95  Real normThresh = 1.0e-30;
96 
97  int num_mg = 1;
98  m_preCondSolver.setSolverParameters(m_num_mg_smooth,
99  m_num_mg_smooth,
100  m_num_mg_smooth,
101  num_mg,
102  maxIterations,
103  solverEps,
104  hang,
105  normThresh);
106 
107  // set preconditioner solver to be _really_ quiet...
108  m_preCondSolver.m_verbosity = 1;
109 
110  // AMRMultiGrid::init has not yet been called
111  // (will call later, during actual solve)
112  m_isPrecondSolverInitialized = false;
113 
114  }
115 
116 
117 }
118 
119  ///
120 template<class T>
122 {
123  CH_TIME("MultilevelLinearOp::~MultilevelLinearOp");
124 
125  if (m_precondBottomSolverPtr != NULL)
126  {
127  delete m_precondBottomSolverPtr;
128  m_precondBottomSolverPtr = NULL;
129  }
130 }
131 
132 ///
133 /**
134  Say you are solving L(phi) = rhs. Make a_lhs = L(a_phi) - a_rhs.
135  If a_homogeneous is true, evaluate the operator using homogeneous
136  boundary conditions.
137 */
138 template<class T>
139 void
141  const Vector<LevelData<T>* >& a_phi,
142  const Vector<LevelData<T>* >& a_rhs,
143  bool a_homogeneous )
144 {
145  CH_TIME("MultilevelLinearOp::residual");
146 
147  // do this using the operators in AMRLevelOp
148  LevelData<T> dummyT;
149 
150  int numLevels = a_lhs.size();
151  CH_assert (a_phi.size() == numLevels);
152  CH_assert (a_rhs.size() == numLevels);
153  CH_assert (m_vectOperators.size() >= numLevels);
154 
155  int maxLevel = numLevels -1;
156  for (int level = m_lBase; level < numLevels; level++)
157  {
158  const LevelData<T>& thisPhi = *(a_phi[level]);
159  LevelData<T>& thisResid = *(a_lhs[level]);
160  const LevelData<T>& thisRHS = *(a_rhs[level]);
161  // theres gotta be a cleaner way to do this...
162  if (level >0)
163  {
164  const LevelData<T>& phiCrse = *a_phi[level-1];
165 
166  if (level < maxLevel)
167  {
168  const LevelData<T>& phiFine = *a_phi[level+1];
169  m_vectOperators[level]->AMRResidual(thisResid,
170  phiFine,
171  thisPhi,
172  phiCrse,
173  thisRHS,
174  a_homogeneous,
175  (m_vectOperators[level+1].operator->()));
176  }
177  else
178  {
179  // no finer level exists
180  m_vectOperators[level]->AMRResidualNF(thisResid,
181  thisPhi,
182  phiCrse,
183  thisRHS,
184  a_homogeneous);
185  } // end if no finer level
186 
187  } // end if a coarser level exists
188  else
189  {
190  // we're on level 0
191  if (level < maxLevel)
192  {
193  const LevelData<T>& phiFine = *a_phi[level+1];
194  m_vectOperators[level]->AMRResidualNC(thisResid,
195  phiFine,
196  thisPhi,
197  thisRHS,
198  a_homogeneous,
199  (m_vectOperators[level+1].operator->()));
200  }
201  else
202  {
203  // no finer level exists
204  m_vectOperators[level]->residual(thisResid,
205  thisPhi,
206  thisRHS,
207  a_homogeneous);
208  } // end if no finer level
209  } // end if on level 0
210  } // end loop over levels
211 }
212 
213 
214 
215 ///
216 /**
217  Given the current state of the residual the correction, apply
218  your preconditioner to a_cor.
219 */
220 template<class T>
221 void
223  const Vector<LevelData<T>* >& a_residual)
224 
225 {
226  CH_TIME("MultilevelLinearOp::preCond");
227 
228  // sanity checks
229  int numLevels = a_cor.size();
230  CH_assert (a_residual.size() == numLevels);
231  CH_assert (m_vectOperators.size() <= numLevels);
232 
233  // use AMR multigrid for a preconditioner
234  if (m_use_multigrid_preconditioner)
235  {
236  CH_TIME("MultilevelLinearOp::preCond::Multigrid");
237 
238  Vector<LevelData<T>* > localCorr;
239  Vector<LevelData<T>* > localResid;
240  create(localCorr, a_cor);
241  create(localResid, a_residual);
242 
243  // this is kinda silly, but I need to get rid of RefCountedPtr
244  // part of things
245  Vector<LevelData<T>* > tempCorr(numLevels, NULL);
246  Vector<LevelData<T>* > tempResid(numLevels, NULL);
247 
248  {
249  CH_TIME("MultilevelLinearOp::preCond::Multigrid::NewStorage");
250 
251  // in place. Otherwise, need to allocate new storage for them
252  for (int lev=0; lev<numLevels; lev++)
253  {
254  tempCorr[lev] = localCorr[lev];
255  tempResid[lev] = const_cast<LevelData<T>*>(localResid[lev]);
256  }
257  }
258 
259  int finestLevel = numLevels -1;
260  {
261  CH_TIME("MultilevelLinearOp::preCond::Multigrid::Initialize");
262 
263  if (!m_isPrecondSolverInitialized)
264  {
265  // because we're not calling solve, need to call init directly
266  m_preCondSolver.init(tempCorr, tempResid, finestLevel,
267  m_lBase);
268  // also need to set bottom solver
269  m_preCondSolver.setBottomSolver(finestLevel, m_lBase);
270  m_isPrecondSolverInitialized = true;
271  }
272  }
273 
274  for (int iter = 0; iter<m_num_mg_iterations; iter++)
275  {
276  CH_TIME("MultilevelLinearOp::preCond::Multigrid::Iterate");
277 
278  // use AMRVCycle instead of solve because we're already
279  // in residual-correction form, so we want to use
280  // homogeneous form of boundary conditions.
281 
282  // however, AMRVCycle requires that the initial correction be zero
283  // to do this properly, recompute residual before calling V-Cycle
284 
285  bool homogeneous = true;
286  residual(localResid, a_cor, a_residual, homogeneous);
287  setToZero(localCorr);
288 
289  m_preCondSolver.AMRVCycle(tempCorr, tempResid,
290  finestLevel, finestLevel,
291  m_lBase);
292 
293  // now increment a_cor with local correction
294  incr(a_cor, localCorr, 1.0);
295  }
296 
297  {
298  CH_TIME("MultilevelLinearOp::preCond::Multigrid::Clear");
299 
300  clear(localCorr);
301  clear(localResid);
302  }
303  }
304  else
305  {
306  CH_TIME("MultilevelLinearOp::preCond::LevelBased");
307 
308  // just use whatever the AMRLevelOp provides for levels finer than lBase
309  for (int level = m_lBase; level < numLevels; level++)
310  {
311  m_vectOperators[level]->preCond(*a_cor[level], *a_residual[level]);
312  }
313  }
314 }
315 
316 
317 ///
318 /**
319  In the context of solving L(phi) = rhs, set a_lhs = L(a_phi).
320  If a_homogeneous is true,
321  evaluate the operator using homogeneous boundary conditions.
322 */
323 template<class T>
324 void
326  const Vector<LevelData<T>* >& a_phi,
327  bool a_homogeneous)
328 {
329  CH_TIME("MultilevelLinearOp::applyOp");
330 
331  int numLevels = a_phi.size();
332  int maxLevel = numLevels-1;
333  LevelData<T> dummyLDF;
334  for (int level = m_lBase; level<numLevels; level++)
335  {
336  if (level == 0)
337  {
338  if (level < maxLevel)
339  {
340  m_vectOperators[level]->AMROperatorNC(*a_lhs[level],
341  *a_phi[level+1],
342  *a_phi[level],
343  a_homogeneous,
344  m_vectOperators[level+1].operator->());
345  }
346  else
347  {
348  // this is the single-level case...
349  m_vectOperators[level]->applyOp(*a_lhs[level],
350  *a_phi[level],
351  a_homogeneous);
352  }
353  } // end if level is 0
354  else
355  {
356  // coarser level exists
357  if (level < maxLevel)
358  {
359  m_vectOperators[level]->AMROperator(*a_lhs[level],
360  *a_phi[level+1],
361  *a_phi[level],
362  *a_phi[level-1],
363  a_homogeneous,
364  m_vectOperators[level+1].operator->());
365  }
366  else
367  {
368  m_vectOperators[level]->AMROperatorNF(*a_lhs[level],
369  *a_phi[level],
370  *a_phi[level-1],
371  a_homogeneous);
372  }
373  } // end if coarser level exists
374  } // end loop over levels
375 
376 }
377 
378 ///
379 /**
380  Create data holder a_lhs that mirrors a_rhs. You do not need
381  to copy the data of a_rhs, just make a holder the same size.
382 */
383 template<class T>
384 void
386  const Vector<LevelData<T>* >& a_rhs)
387 {
388  CH_TIME("MultilevelLinearOp::create");
389 
390  a_lhs.resize(a_rhs.size());
391  // need to create an lBase-1 just in case it's needed for C/F BCs
392  int coarsestLevel = Max(0, m_lBase-1);
393  for (int level =coarsestLevel; level < a_rhs.size(); level++)
394  {
395  a_lhs[level] = new LevelData<T>;
396  m_vectOperators[level]->create(*a_lhs[level],
397  *a_rhs[level]);
398  }
399 }
400 
401 
402 ///
403 /**
404  Clean up data holder before it goes out of scope. This is necessary
405  because create calls new.
406 */
407 template<class T>
408 void
410 {
411  CH_TIME("MultilevelLinearOp::clear");
412 
413  for (int level =0; level < a_lhs.size(); level++)
414  {
415  if (a_lhs[level] != NULL)
416  {
417  delete a_lhs[level];
418  a_lhs[level] = NULL;
419  }
420  }
421 }
422 
423 
424 ///
425 /**
426  Set a_lhs equal to a_rhs.
427 */
428 template<class T>
429 void
431  const Vector<LevelData<T>* >& a_rhs)
432 {
433  CH_TIME("MultilevelLinearOp::assign");
434 
435  CH_assert (a_lhs.size() == a_rhs.size());
436  // only do this for lBase and finer levels
437  for (int level = m_lBase; level < a_lhs.size(); level++)
438  {
439  if (!(a_lhs[level] ==NULL) && !(a_rhs[level] == NULL))
440  {
441  m_vectOperators[level]->assign(*a_lhs[level], *a_rhs[level]);
442  }
443  }
444 }
445 
446 ///
447 /**
448  Compute and return the dot product of a_1 and a_2. In most
449  contexts, this means return the sum over all data points of a_1*a_2.
450 */
451 template<class T>
452 Real
454  const Vector<LevelData<T>* >& a_2)
455 {
456  CH_TIME("MultilevelLinearOp::dotProduct");
457 
458  // want to do this in an AMR way, so need to generate a temporary
459  Vector<LevelData<T>* > temp1, temp2;
460  create(temp1, a_1);
461  create (temp2, a_2);
462 
463  // first set to zero, then call assign, since assign only sets
464  // valid regions (that way ghost cells are set to zero)
465 
466  setToZero(temp1);
467  setToZero(temp2);
468 
469  assign(temp1, a_1);
470  assign(temp2, a_2);
471 
472  // now set covered regions to zero
473  for (int level =m_lBase; level<temp1.size()-1; level++)
474  {
475  LevelData<T>& temp1Level = *temp1[level];
476  LevelData<T>& temp2Level = *temp2[level];
477 
478  CH_assert(temp1[level]->getBoxes() == temp2[level]->getBoxes());
479  CH_assert(temp1[level+1] != NULL);
480 
481  int nRefFine = m_refRatios[level];
482  const DisjointBoxLayout& finerGrids = temp1[level+1]->getBoxes();
483  const DisjointBoxLayout& levelGrids = temp1[level]->getBoxes();
484  DataIterator levelDit = levelGrids.dataIterator();
485  LayoutIterator finerLit = finerGrids.layoutIterator();
486 
487  for (levelDit.begin(); levelDit.ok(); ++levelDit)
488  {
489  const Box& thisBox = levelGrids[levelDit];
490  for (finerLit.begin(); finerLit.ok(); ++finerLit)
491  {
492  Box testBox = finerGrids[finerLit];
493  testBox.coarsen(nRefFine);
494  testBox &= thisBox;
495  if (!testBox.isEmpty())
496  {
497  temp1Level[levelDit].setVal(0.0, testBox,
498  0, temp1Level.nComp());
499 
500  temp2Level[levelDit].setVal(0.0, testBox,
501  0, temp2Level.nComp());
502  }
503  } // end loop over finer boxes
504  } // end loop over boxes on this level
505  } // end loop over levels for setting covered regions to zero;
506 
507  // now loop over levels and call AMRLevelOp dotProduct
508  Real prod = 0.0;
509  for (int level =m_lBase; level<temp1.size(); level++)
510  {
511  Real levelProd = m_vectOperators[level]->dotProduct(*temp1[level],
512  *temp2[level]);
513 
514  // incorporate scaling for AMR dot products
515  RealVect dxLev = m_vectDx[level];
516  Real scale = D_TERM(dxLev[0], *dxLev[1], *dxLev[2]);
517  levelProd *= scale;
518  prod += levelProd;
519  }
520 
521  clear(temp1);
522  clear (temp2);
523 
524  return prod;
525 }
526 
527 ///
528 /**
529  Increment by scaled amount (a_lhs += a_scale*a_x).
530 */
531 template<class T>
532 void
534  const Vector<LevelData<T>* >& a_x,
535  Real a_scale)
536 {
537  CH_TIME("MultilevelLinearOp::incr");
538 
539  // int this case, only operate on lBase and finer
540  for (int level=m_lBase; level<a_lhs.size(); level++)
541  {
542  m_vectOperators[level]->incr(*a_lhs[level], *a_x[level], a_scale);
543  }
544 }
545 
546 ///
547 /**
548  Set input to a scaled sum (a_lhs = a_a*a_x + a_b*a_y).
549 */
550 template<class T>
551 void
553  const Vector<LevelData<T>* >& a_x,
554  const Vector<LevelData<T>* >& a_y,
555  Real a_a,
556  Real a_b)
557 {
558  CH_TIME("MultilevelLinearOp::axby");
559 
560  // only do this for lBase and finer
561  for (int level=m_lBase; level<a_lhs.size(); ++level)
562  {
563  m_vectOperators[level]->axby(*a_lhs[level],
564  *a_x[level],
565  *a_y[level],
566  a_a, a_b);
567  }
568 }
569 
570 ///
571 /**
572  Multiply the input by a given scale (a_lhs *= a_scale).
573 */
574 template<class T>
575 void
577  const Real& a_scale)
578 {
579  CH_TIME("MultilevelLinearOp::scale");
580 
581  // only do this for lBase and finer
582  for (int level =m_lBase; level<a_lhs.size(); level++)
583  {
584  m_vectOperators[level]->scale(*a_lhs[level], a_scale);
585  }
586 
587 }
588 
589 ///
590 /**
591  Return the norm of a_rhs.
592  a_ord == 0 max norm, a_ord == 1 sum(abs(a_rhs)), else, L(a_ord) norm.
593 */
594 template<class T>
595 Real
597  int a_ord)
598 {
599  CH_TIME("MultilevelLinearOp::norm");
600 
601  // now loop over levels and call AMRLevelOp AMRNorm
602  int maxLevel = a_rhs.size()-1;
603  Real thisNorm = 0.0;
604  for (int level = m_lBase; level<a_rhs.size(); level++)
605  {
606  Real normLevel;
607 
608  if (level < maxLevel)
609  {
610  // finer level exists
611  int refRatio = m_refRatios[level];
612  normLevel = m_vectOperators[level]->AMRNorm(*a_rhs[level],
613  *a_rhs[level+1],
614  refRatio,
615  a_ord);
616  }
617  else
618  {
619  // no finer level exists
620  LevelData<T> tempLDF;
621  int refRatio = -1;
622  normLevel = m_vectOperators[level]->AMRNorm(*a_rhs[level],
623  tempLDF,
624  refRatio,
625  a_ord);
626  }
627 
628 
629  // now need to do a bit of rescaling to make things work out right
630  if (a_ord > 2)
631  {
632  MayDay::Error("MultilevelLinearOp::norm -- undefined for order > 2");
633  }
634 
635  if (a_ord == 2)
636  {
637  normLevel = normLevel*normLevel;
638  }
639 
640  if (a_ord != 0)
641  {
642  thisNorm += normLevel;
643  }
644  else if (normLevel > thisNorm)
645  {
646  thisNorm = normLevel;
647  }
648  } // end
649 
650  if (a_ord == 2)
651  {
652  thisNorm = sqrt(thisNorm);
653  }
654 
655  return thisNorm;
656 }
657 
658 ///
659 /**
660  Set a_lhs to zero.
661 */
662 template<class T>
663 void
665 {
666  CH_TIME("MultilevelLinearOp::setToZero");
667 
668  // do this for lBase -1 if necessary, to account for cases where correction
669  // lBase-1 needs to be zero
670  int coarsestLevel = Max(0, m_lBase-1);
671  for (int level = coarsestLevel; level< a_lhs.size(); level++)
672  {
673  m_vectOperators[level]->setToZero(*a_lhs[level]);
674  }
675 }
676 
677 template<class T>
678 void
680  const char* a_filename)
681 {
682 #ifdef CH_USE_HDF5
683  writeVectorLevelName(a_data, &m_refRatios, a_filename);
684 #else
685  MayDay::Warning("MultilevelLinearOp<T>::write unimplemented since CH_USE_HDF5 undefined");
686 #endif
687 }
688 
689 #include "NamespaceFooter.H"
690 
691 #endif
Box & coarsen(int refinement_ratio)
coarsening
virtual void incr(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_x, Real a_scale)
Definition: MultilevelLinearOpI.H:533
virtual Real norm(const Vector< LevelData< T > * > &a_rhs, int a_ord)
Definition: MultilevelLinearOpI.H:596
A reference-counting handle class.
Definition: RefCountedPtr.H:66
virtual void write(const Vector< LevelData< T > * > *a_data, const char *a_filename)
Definition: MultilevelLinearOpI.H:679
#define CH_assert(cond)
Definition: CHArray.H:37
int m_verbosity
Definition: BiCGStabSolver.H:78
void begin()
initialize this iterator to the first Box in its Layout
Definition: LayoutIterator.H:115
virtual bool ok() const
return true if this iterator is still in its Layout
Definition: LayoutIterator.H:110
IndexTM< T, N > scale(const IndexTM< T, N > &a_p, T a_s)
Definition: IndexTMI.H:384
Definition: DataIterator.H:140
An Iterator based on a BoxLayout object.
Definition: LayoutIterator.H:38
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:552
void resize(unsigned int isize)
Definition: Vector.H:323
virtual ~MultilevelLinearOp()
Definition: MultilevelLinearOpI.H:121
virtual void setToZero(Vector< LevelData< T > * > &a_lhs)
Definition: MultilevelLinearOpI.H:664
#define CH_TIME(name)
Definition: CH_Timer.H:59
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:576
Definition: BoxLayoutData.H:136
double Real
Definition: REAL.H:33
A BoxLayout that has a concept of disjointedness.
Definition: DisjointBoxLayout.H:31
virtual void clear(Vector< LevelData< T > * > &a_lhs)
Definition: MultilevelLinearOpI.H:409
size_t size() const
Definition: Vector.H:177
static void Error(const char *const a_msg=m_nullString, int m_exitCode=CH_DEFAULT_ERROR_CODE)
Print out message to cerr and exit with the specified exit code.
virtual void assign(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_rhs)
Definition: MultilevelLinearOpI.H:430
virtual void applyOp(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_phi, bool a_homogeneous=false)
Definition: MultilevelLinearOpI.H:325
static void Warning(const char *const a_msg=m_nullString)
Print out message to cerr and continue.
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:140
A Rectangular Domain on an Integer Lattice.
Definition: Box.H:465
A Real vector in SpaceDim-dimensional space.
Definition: RealVect.H:41
virtual void preCond(Vector< LevelData< T > * > &a_cor, const Vector< LevelData< T > * > &a_residual)
Apply preconditioner.
Definition: MultilevelLinearOpI.H:222
int nComp() const
Definition: BoxLayoutData.H:258
T Max(const T &a_a, const T &a_b)
Definition: Misc.H:39
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
virtual Real dotProduct(const Vector< LevelData< T > * > &a_1, const Vector< LevelData< T > * > &a_2)
Definition: MultilevelLinearOpI.H:453
void writeVectorLevelName(const Vector< LevelData< FArrayBox > *> *a_dataPtr, const Vector< int > *a_refRatios, const char *a_filename)
LayoutIterator layoutIterator() const
Iterator that processes through ALL the boxes in a BoxLayout.
Definition: AMRMultiGrid.H:231
DataIterator dataIterator() const
Parallel iterator.
bool isEmpty() const
{ Comparison Functions}
Definition: Box.H:1863
virtual void create(Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_rhs)
Definition: MultilevelLinearOpI.H:385
Definition: BiCGStabSolver.H:26