Chombo + EB  3.2
EBAMRPoissonOp.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 _EBAMRPOISSONOP_H_
12 #define _EBAMRPOISSONOP_H_
13 
14 #include "REAL.H"
15 #include "Box.H"
16 #include "FArrayBox.H"
17 #include "Vector.H"
18 #include <map>
19 #include "RefCountedPtr.H"
20 
21 #include "AMRMultiGrid.H"
22 
23 #include "EBIndexSpace.H"
24 #include "EBCellFAB.H"
25 #include "EBCellFactory.H"
26 ///only use the non-aggregated for testing
27 #define USE_NONAGG 0
28 
29 #if USE_NONAGG==1
30 #define EBSTENCIL_T NonAggregatedEBStencil
31 #else
32 #define EBSTENCIL_T EBStencil
33 #endif
34 
35 #include "EBStencil.H"
36 #include "NonAggregatedEBStencil.H"
37 
38 #include "EBLevelDataOps.H"
39 #include "BaseEBBC.H"
40 #include "BaseDomainBC.H"
41 #include "CFIVS.H"
42 #include "EBFastFR.H"
43 #include "EBMGAverage.H"
44 #include "EBMGInterp.H"
45 #include "PolyGeom.H"
46 #include "EBQuadCFInterp.H"
47 #include "EBLevelGrid.H"
48 #include "AMRTGA.H"
49 #include "EBAMRIO.H"
50 #include "AMRPoissonOp.H"
51 #include "CFRegion.H"
52 #include "NamespaceHeader.H"
53 
54 #ifdef CH_USE_PETSC
55 #include "petsc.h"
56 #include "petscmat.h"
57 #include "petscksp.h"
58 #include "petscviewer.h"
59 #endif
60 
61 #if CH_SPACEDIM==2
62 #define EBAMRPO_NUMSTEN 4
63 #elif CH_SPACEDIM==3
64 #define EBAMRPO_NUMSTEN 8
65 #else
67 {
68  THIS_WILL_ONLY_COMPILE_WHEN_CH_SPACEDIM_IS_2_OR_3;
69 }
70 #endif
71 
72 ///
73 /**
74  Operator to solve (alpha + beta lapl)phi = rhs. This follows the AMRLevelOp interface.
75 */
76 class EBAMRPoissonOp: public LevelTGAHelmOp<LevelData<EBCellFAB>, EBFluxFAB >
77 {
78 public:
79  void
80  getVoFStencil(VoFStencil& a_stenc,
81  const VolIndex& a_vof,
82  const DataIndex& a_dit) const;
83 
84 #ifdef CH_USE_HDF5
85  ///
86  virtual void dumpAMR(Vector<LevelData<EBCellFAB>*>& a_data, string name)
87  {
88  writeEBAMRname(&a_data, name.c_str());
89  }
90 
91  virtual void dumpLevel(LevelData<EBCellFAB>& a_data, string name)
92  {
93  writeEBLevelname(&a_data, name.c_str());
94  }
95 #endif
96 
97 #ifdef CH_USE_PETSC
98  ///Fill a_petsc_mat with a matrix that describes the operator stencil
99  int getPetscMatrix(Mat& a_petsc_mat);
100 
101  ///put the real components into comp 0, the imaginary ones into comp 1
102  int getLevelDataFromPetscVector(LevelData<EBCellFAB>& a_data, const Vec& a_petsc_vec_real, const Vec& a_petsc_vec_imag );
103 
104  ///
105  /**
106  get the matrix indicies associated with every point on the level
107  Not really going to work if level does not cover the domain
108  */
109  static int getMatrixIndexingLD(LevelData<BaseEBCellFAB<int> >& a_gids, int & a_data,
110  const EBLevelGrid& a_eblg,
111  const IntVect & a_ghostCellsPhi,
112  const bool & a_hasCoar);
113 #endif
114  static void
116  const ProblemDomain & a_domainCoar,
117  const EBIndexSpace * const a_ebisPtr,
118  const int & a_maxBoxSize);
119 
120  ///
121  /**
122  version that does not fill ebislcoar
123  */
124  static bool getCoarserLayouts(DisjointBoxLayout& a_dblCoar,
125  ProblemDomain& a_domainCoar,
126  const DisjointBoxLayout& a_dblFine,
127  const EBISLayout& a_ebislFine,
128  const ProblemDomain& a_domainFine,
129  int a_refToCoar,
130  const EBIndexSpace* a_ebisPtr,
131  int a_maxBoxSize,
132  bool& a_layoutChanged,
133  int a_testRef = 2);
134 
135  static Real staticMaxNorm(const LevelData<EBCellFAB>& a_rhs, const EBLevelGrid& a_eblg);
136 
137  //for tga to reset stuff
138  virtual void setAlphaAndBeta(const Real& a_alpha,
139  const Real& a_beta);
140 
141  //another tgaism
142  virtual void diagonalScale(LevelData<EBCellFAB>& a_rhs,
143  bool a_kappaWeighted = true);
145  {
146  //no acoef here.
147  }
148 
149  virtual void kappaScale(LevelData<EBCellFAB> & a_rhs)
150  {
151  //since this is a constant coefficient operator with a=1,
152  //this means the same thing as diagonal scale
153  diagonalScale(a_rhs);
154 
155  }
156 
157  ///returns m_dx, such function is required by some LinearSolvers
158  Real dx() const
159  {
160  return m_dx[0];
161  }
162 
163  ///a leveltgaism
164  virtual void fillGrad(const LevelData<EBCellFAB>& a_phi)
165  {;}
166 
167  ///
168  /**
169  dump stencil as matrix to stdout
170  */
172  {
173  public:
175  int& a_phase)
176  {
177  m_vof = a_vof;
178  m_phase = a_phase;
179  }
180 
182  {
183  MayDay::Error("No weak construction of StencilIndex class.");
184  }
185 
187  {
188  m_vof =a_sin.m_vof;
189  m_phase = a_sin.m_phase;
190  return *this;
191  }
192 
193  VolIndex vof() const
194  {
195  return m_vof;
196  }
197 
198  int phase() const
199  {
200  return m_phase;
201  }
202 
203  bool operator!=(const StencilIndex& a_sin) const
204  {
205  return m_vof!= a_sin.m_vof || m_phase!=a_sin.m_phase;
206  }
207 
208  protected:
210  int m_phase;
211  };
212 
214  {
215  public:
216  bool operator() (const StencilIndex& a_s1, const StencilIndex& a_s2) const
217  {
218  /*
219  int p1 = a_s1.phase();
220  int p2 = a_s2.phase();
221  if (p1 == p2)
222  {
223  return (a_s1.vof()<a_s2.vof());
224  }
225  else
226  {
227  const IntVect& iv1 = a_s1.vof().gridIndex();
228  const IntVect& iv2 = a_s2.vof().gridIndex();
229  if (iv1 == iv2)
230  {
231  return (p1<p2);
232  }
233  else
234  {
235  return (iv1<iv2);
236  }
237  }
238  */
239  const IntVect& iv1 = a_s1.vof().gridIndex();
240  const IntVect& iv2 = a_s2.vof().gridIndex();
241  for (int idir=0; idir<SpaceDim; ++idir)
242  {
243  if (iv1[idir] != iv2[idir])
244  {
245  return (iv1[idir]<iv2[idir]);
246  }
247  }
248  int p1 = a_s1.phase();
249  int p2 = a_s2.phase();
250  if (p1 == p2)
251  {
252  int ci1 = a_s1.vof().cellIndex();
253  int ci2 = a_s2.vof().cellIndex();
254  return (ci1<ci2);
255  }
256  else
257  {
258  return (p1<p2);
259  }
260  }
261  };
262 
263  /*
264  Calculate stencil weights for each VoF. A mapping of each VolIndex
265  an integer is also performed. The first output is this mapping, as
266  a matrix that can be read into Matlab/Octave. The second output is
267  the stencil weight for each VoF pair, in matrix form.
268  The volume fraction is an additional output of the mapping matrix.
269  The matrices are output in the format:
270 
271  M_(RES)(imap, :) = [ gridIndex[0] gridIndex[1] (gridIndex[2]) alphaWeight kappa];
272  where 'RES' is the resolution at the EBPoissonOp's level, 'imap' is
273  the integer the VoF has been mapper to, 'alphaWeight' is the
274  weighting used to get good conditioning of the matrix, and kappa is
275  the VoF's volume fraction.
276 
277  L_(RES)(imap, jmap) = value;
278  where 'RES' is as above, 'imap' is the integer that the VoF whose
279  stencil is being computed has been mapped to, 'jmap' is the integer
280  for a VoF in that stencil. 'value' is the stencil weight.
281  */
282  void dumpStencilMatrix();
283 
284  /*
285  Calculates the stencil for VoFs at the domain boundary.
286  */
287  void getDomainFluxStencil( VoFStencil& a_stencil,
288  const VolIndex& a_vof,
289  const int a_comp,
290  const DataIndex& a_dit);
291  /*
292  A method for testing EBAMRPoissonOp::dumpStencilMatrix. Since it
293  uses calls to 'applyOp' to get the stencil weights, it is slow
294  but correct.
295  */
297 
298  virtual void getFlux(EBFluxFAB& a_flux,
299  const LevelData<EBCellFAB>& a_data,
300  const Box& a_grid,
301  const DataIndex& a_dit,
302  Real a_scale)
303  {
304  for (int idir = 0; idir < SpaceDim; idir++)
305  {
306  Box ghostedBox = a_grid;
307  ghostedBox.grow(1);
308  ghostedBox.grow(idir,-1);
309  ghostedBox &= m_eblg.getDomain();
310 
311  getFlux(a_flux[idir], a_data[a_dit], ghostedBox, a_grid,
312  m_eblg.getDomain(),
313  m_eblg.getEBISL()[a_dit], m_dx, idir);
314  }
315  }
317  {
318  return m_eblg;
319  }
321  {
322  return m_eblgCoarMG;
323  }
324  static void setOperatorTime(Real a_time)
325  {
326  s_time = a_time;
327  }
328 
329  //conforms to AMRTGA
330  virtual void setTime(Real a_time)
331  {
332  setOperatorTime(a_time);
333  }
334  ///
335  virtual ~EBAMRPoissonOp();
336 
337  ///
338  EBAMRPoissonOp();
339 
340  ///
341  /** a_residual = a_rhs - L(a_phiFine, a_phi) no coaser AMR level*/
342  void AMRResidualNC(LevelData<EBCellFAB>& a_residual,
343  const LevelData<EBCellFAB>& a_phiFine,
344  const LevelData<EBCellFAB>& a_phi,
345  const LevelData<EBCellFAB>& a_rhs,
346  bool a_homogeneousBC,
347  AMRLevelOp<LevelData<EBCellFAB> >* a_finerOp);
348 
349 
350  ///
351  /** apply AMR operator no coaser AMR level*/
352  void AMROperatorNC(LevelData<EBCellFAB>& a_LofPhi,
353  const LevelData<EBCellFAB>& a_phiFine,
354  const LevelData<EBCellFAB>& a_phi,
355  bool a_homogeneousBC,
356  AMRLevelOp<LevelData<EBCellFAB> >* a_finerOp);
357 
358  ///
359  /**
360  If you are approaching this operator from this interface, consider backing away and using
361  EBAMRPoissonOpFactory to generate these objects. Really.
362  a_eblgFine, : grid at finer level \\
363  a_eblg, : grid at this level \\
364  a_eblgCoar, : grid at coarser level \\
365  a_eblgCoarMG, : grid at intermediate multigrid level \\
366  a_domainBC, : domain boundary conditions at this level \\
367  a_ebBC: eb boundary conditions at this level \\
368  a_dx: grid spacing at this level \\
369  a_origin: offset to lowest corner of the domain \\
370  a_refToFine: refinement ratio to finer level \\
371  a_refToCoar: refinement ratio to coarser level \\
372  a_hasFiner: true if there is a finer AMR level, false otherwise. \\
373  a_hasCoarser: true if there is a coarser AMR level. \\
374  a_hasCoarserMG: true if there is a coarser MultiGrid level. \\
375  a_preCondIters: number of iterations to do for pre-conditioning \\
376  a_relaxType: 0 means point Jacobi, 1 is Gauss-Seidel. \\
377  a_alpha: coefficent of identity \\
378  a_beta: coefficient of laplacian.\\
379  a_ghostCellsPhi: Number of ghost cells in phi, correction\\
380  a_ghostCellsRhs: Number of ghost cells in RHS, residual, lphi\\
381  Ghost cell arguments are there for caching reasons. Once you set them, an error is thrown if
382  you send in data that does not match.
383  */
384 EBAMRPoissonOp(const EBLevelGrid & a_eblgFine,
385  const EBLevelGrid & a_eblg,
386  const EBLevelGrid & a_eblgCoar,
387  const EBLevelGrid & a_eblgCoarMG,
388  const RefCountedPtr<EBQuadCFInterp>& a_quadCFI,
389  const RefCountedPtr<BaseDomainBC>& a_domainBC,
390  const RefCountedPtr<BaseEBBC>& a_ebBC,
391  const RealVect& a_dx,
392  const RealVect& a_dxCoar,
393  const RealVect& a_origin,
394  const int& a_refToFine,
395  const int& a_refToCoar,
396  const bool& a_hasFine,
397  const bool& a_hasCoar,
398  const bool& a_hasMGObjects,
399  const bool& a_layoutChanged,
400  const int& a_numPreCondIters,
401  const int& a_relaxType,
402  const Real& a_alpha,
403  const Real& a_beta,
404  const IntVect& a_ghostCellsPhi,
405  const IntVect& a_ghostCellsRHS,
406  int a_testRef = 2);
407 
408  //MGOp operations. no finer or coarser
409 
410  ///
411  /**
412  */
413  virtual void residual(LevelData<EBCellFAB>& a_residual,
414  const LevelData<EBCellFAB>& a_phi,
415  const LevelData<EBCellFAB>& a_rhs,
416  bool a_homogeneousPhysBC=false);
417 
418  ///
419  /**
420  For debugging purpose, get the matrix of the op (could be change from time to time due to EBBC update) by calculating 0 - L(phi(i)) for i=0,1,...,n, each time phi(i)=1 for i-th cell and 0 for all other cells
421  */
422 
423  virtual void getOpMatrix(const LevelData<EBCellFAB>& a_phi,
424  const LevelData<EBCellFAB>& a_rhs);
425 
426  ///
427  /**
428  */
429  virtual void preCond(LevelData<EBCellFAB>& a_opPhi,
430  const LevelData<EBCellFAB>& a_phi);
431 
432  ///
433  /**
434  This function assumes that coarse-fine boundary condtions have
435  been dealt with.
436  */
437  virtual void applyOp(LevelData<EBCellFAB>& a_opPhi,
438  const LevelData<EBCellFAB>& a_phi,
439  const LevelData<EBCellFAB>* const a_phiCoarse,
440  const bool& a_homogeneousPhysBC,
441  const bool& a_homogeneousCFBC);
442 
443  /// virtual function called by LevelTGA
444  virtual void applyOpNoBoundary(LevelData<EBCellFAB>& a_opPhi,
445  const LevelData<EBCellFAB>& a_phi);
446 
447  ///
448  /**
449  get EB flux from leveldata if we have one. otherwise use m_ebBC
450  ignore input data in the case of homogeneous Phys BC
451  */
452  virtual void
453  applyOp(LevelData<EBCellFAB>& a_opPhi,
454  const LevelData<EBCellFAB>& a_phi,
455  const LevelData<EBCellFAB>* const a_phiCoar,
456  const bool& a_homogeneousPhysBC,
457  const bool& a_homogeneousCFBC,
458  const LevelData<BaseIVFAB<Real> >* const a_ebFluxBCLD);
459 
460  virtual void
461  applyOp(LevelData<EBCellFAB>& a_opPhi,
462  const LevelData<EBCellFAB>& a_phi,
463  const LevelData<EBCellFAB>* const a_phiCoar,
464  DataIterator& a_dit,
465  const bool& a_homogeneousPhysBC,
466  const bool& a_homogeneousCFBC,
467  const LevelData<BaseIVFAB<Real> >* const a_ebFluxBCLD);
468 
469  virtual void
471  const BaseFab<Real>& a_rhs,
472  const int& a_icolor,
473  const Real& a_weight,
474  const bool& a_homogeneousPhysBC,
475  const DataIndex& a_dit);
476 
477  virtual void
479  const EBCellFAB& a_rhs,
480  const int& a_icolor,
481  const bool& a_homogeneousPhysBC,
482  const DataIndex& a_dit);
483 
484  virtual void
486  const LevelData<EBCellFAB>& a_phiOld,
487  const LevelData<EBCellFAB>& a_rhs,
488  const int& a_icolor,
489  const Real& a_weight,
490  const bool& a_homogeneousPhysBC);
491 
492  virtual void
494  const LevelData<EBCellFAB>& a_phiOld,
495  const LevelData<EBCellFAB>& a_rhs,
496  const int& a_icolor,
497  const bool& a_homogeneousPhysBC);
498 
499  ///apply domainflux in multivariable mode
500  void
501  mvApplyDomainFlux(BaseFab<Real> & a_phiFAB,
502  const Box& a_grid,
503  const DataIndex& a_dit);
504  /***/
505  ///for taking the multi-variable laplacian (think source term in viscous flow)
506  /// evaluates beta*laplacian (with inhomogeneous bcs).
507  /***/
508  void
510  const EBCellFAB & a_phi,
511  const DataIndex & a_dit);
512  ///
513  /**
514  this is the linearop function. CFBC is set to homogeneous. phic is null
515  */
516  virtual void applyOp(LevelData<EBCellFAB>& a_opPhi,
517  const LevelData<EBCellFAB>& a_phi,
518  bool a_homogeneousPhysBC);
519 
520  /// no exchange of cf interp
521  void
523  const LevelData<EBCellFAB>& a_phi,
524  const LevelData<EBCellFAB>* const a_phiCoar,
525  DataIterator& a_dit,
526  const bool& a_homogeneousPhysBC,
527  const bool& a_homogeneousCFBC,
528  const LevelData<BaseIVFAB<Real> >* const a_ebFluxBCLD //only non null in multifluid
529  );
530 
531  ///
532  /**
533  */
534  virtual void create(LevelData<EBCellFAB>& a_lhs,
535  const LevelData<EBCellFAB>& a_rhs);
536 
537  ///
538  virtual void createCoarsened(LevelData<EBCellFAB>& a_lhs,
539  const LevelData<EBCellFAB>& a_rhs,
540  const int& a_refRat);
541 
542  virtual void buildCopier(Copier & a_copier,
543  const LevelData<EBCellFAB>& a_lhs,
544  const LevelData<EBCellFAB>& a_rhs);
545 
546 
547  virtual void assignCopier(LevelData<EBCellFAB> & a_lhs,
548  const LevelData<EBCellFAB>& a_rhs,
549  const Copier & a_copier);
550 
551  Real
552  AMRNorm(const LevelData<EBCellFAB>& a_coarResid,
553  const LevelData<EBCellFAB>& a_fineResid,
554  const int& a_refRat,
555  const int& a_ord);
556 
557  ///
558  /**
559  */
560  virtual void assign(LevelData<EBCellFAB>& a_lhs,
561  const LevelData<EBCellFAB>& a_rhs);
562 
563  ///copier definition was killing us.
564  virtual void
566  const LevelData<EBCellFAB>& a_rhs);
567 
568 
569  ///
570  /**
571  */
572  virtual Real dotProduct(const LevelData<EBCellFAB>& a_1,
573  const LevelData<EBCellFAB>& a_2);
574 
575  ///
576  /**
577  */
578  virtual void incr(LevelData<EBCellFAB>& a_lhs,
579  const LevelData<EBCellFAB>& a_x,
580  Real a_scale);
581 
582  ///
583  /**
584  */
585  virtual void axby(LevelData<EBCellFAB>& a_lhs,
586  const LevelData<EBCellFAB>& a_x,
587  const LevelData<EBCellFAB>& a_y,
588  Real a_a,
589  Real a_b);
590 
591  ///
592  /**
593  */
594  virtual void scale(LevelData<EBCellFAB>& a_lhs,
595  const Real& a_scale);
596 
597  ///
598  /**
599  */
600  virtual Real norm(const LevelData<EBCellFAB>& a_rhs,
601  int a_ord);
602 
603  ///
604  /**
605  */
606  virtual Real localMaxNorm(const LevelData<EBCellFAB>& a_rhs);
607 
608  ///
609  /**
610  */
611  virtual void setToZero(LevelData<EBCellFAB>& a_lhs);
612 
613  ///
614  /**
615  */
616  virtual void setVal(LevelData<EBCellFAB>& a_lhs, const Real& a_value);
617 
618  ///
619  /**
620  */
621  virtual void createCoarser(LevelData<EBCellFAB>& a_coarse,
622  const LevelData<EBCellFAB>& a_fine,
623  bool a_ghosted);
624 
625  ///
626  /**
627  */
628  virtual void relax(LevelData<EBCellFAB>& a_e,
629  const LevelData<EBCellFAB>& a_residual,
630  int a_iterations);
631 
632  ///
633  /**
634  Calculate restricted residual:
635  a_resCoarse[2h] = I[h->2h] (a_rhsFine[h] - L[h](a_phiFine[h]))
636  */
637  virtual void restrictResidual(LevelData<EBCellFAB>& a_resCoarse,
638  LevelData<EBCellFAB>& a_phiFine,
639  const LevelData<EBCellFAB>& a_rhsFine);
640 
641  ///
642  /**
643  Correct the fine solution based on coarse correction:
644  a_phiThisLevel += I[2h->h] (a_correctCoarse)
645  */
646  virtual void prolongIncrement(LevelData<EBCellFAB>& a_phiThisLevel,
647  const LevelData<EBCellFAB>& a_correctCoarse);
648 
649  ///
650  /** Refinement ratio between this level and coarser level.
651  Returns 1 when there are no coarser AMRLevelOp objects */
652  virtual int refToCoarser();
653 
654  ///
655  /** Refinement ratio between this level and coarser level.
656  Returns 1 when there are no coarser AMRLevelOp objects */
657  virtual int refToFiner();
658 
659  ///
660  /** a_residual = a_rhs - L(a_phi, a_phiFine, a_phiCoarse) */
661  virtual void AMRResidual(LevelData<EBCellFAB>& a_residual,
662  const LevelData<EBCellFAB>& a_phiFine,
663  const LevelData<EBCellFAB>& a_phi,
664  const LevelData<EBCellFAB>& a_phiCoarse,
665  const LevelData<EBCellFAB>& a_rhs,
666  bool a_homogeneousBC,
667  AMRLevelOp<LevelData<EBCellFAB> >* a_finerOp);
668 
669  ///
670  /** a_residual = a_rhs - L(a_phi, a_phiCoarse) */
671  virtual void AMRResidualNF(LevelData<EBCellFAB>& a_residual,
672  const LevelData<EBCellFAB>& a_phi,
673  const LevelData<EBCellFAB>& a_phiCoarse,
674  const LevelData<EBCellFAB>& a_rhs,
675  bool a_homogeneousBC);
676 
677 
678  ///
679  /** a_residual = a_rhs - L(a_phi, a_phiFine, a_phiCoarse) */
680  virtual void AMROperator(LevelData<EBCellFAB>& a_LofPhi,
681  const LevelData<EBCellFAB>& a_phiFine,
682  const LevelData<EBCellFAB>& a_phi,
683  const LevelData<EBCellFAB>& a_phiCoarse,
684  bool a_homogeneousBC,
685  AMRLevelOp<LevelData<EBCellFAB> >* a_finerOp);
686 
687  ///
688  /** a_residual = a_rhs - L(a_phi, a_phiCoarse) */
689  virtual void AMROperatorNF(LevelData<EBCellFAB>& a_LofPhi,
690  const LevelData<EBCellFAB>& a_phi,
691  const LevelData<EBCellFAB>& a_phiCoarse,
692  bool a_homogeneousBC);
693 
694  ///
695  /** a_resCoarse = I[h-2h] (a_residual - L(a_correction, a_coarseCorrection)) */
696  virtual void AMRRestrict(LevelData<EBCellFAB>& a_resCoarse,
697  const LevelData<EBCellFAB>& a_residual,
698  const LevelData<EBCellFAB>& a_correction,
699  const LevelData<EBCellFAB>& a_coarseCorrection,
700  bool a_skip_res = false );
701 
702  ///
703  /** a_correction += I[2h->h](a_coarseCorrection) */
704  virtual void AMRProlong(LevelData<EBCellFAB>& a_correction,
705  const LevelData<EBCellFAB>& a_coarseCorrection);
706 
707  ///
708  /** a_residual = a_residual - L(a_correction, a_coarseCorrection) */
709  virtual void AMRUpdateResidual(LevelData<EBCellFAB>& a_residual,
710  const LevelData<EBCellFAB>& a_correction,
711  const LevelData<EBCellFAB>& a_coarseCorrection);
712 
713  ///
714  /** a_residual = a_residual - L(a_correction, a_coarseCorrection)
715  used in multifluid */
716  void AMRUpdateResidual(LevelData<EBCellFAB>& a_residual,
717  const LevelData<EBCellFAB>& a_correction,
718  const LevelData<EBCellFAB>& a_coarseCorrection,
719  const LevelData<BaseIVFAB<Real> >* const a_ebFluxBCLD);
720 
721 
722  void reflux(LevelData<EBCellFAB>& a_residual,
723  const LevelData<EBCellFAB>& a_phiFine,
724  const LevelData<EBCellFAB>& a_phi,
725  AMRLevelOp<LevelData<EBCellFAB> >* a_finerOp);
726 
727 
728  void fast_reflux(LevelData<EBCellFAB>& a_residual,
729  const LevelData<EBCellFAB>& a_phiFine,
730  const LevelData<EBCellFAB>& a_phi,
731  AMRLevelOp<LevelData<EBCellFAB> >* a_finerOp);
732 
733  void setEBBC(const RefCountedPtr<BaseEBBC>& a_ebBC);
734 
735  //for debugging other operators
737  const LevelData<EBCellFAB>& a_res,
738  int a_iredblack);
739 
740  //for debugging other operators
742  const LevelData<EBCellFAB>& a_rhs);
743 
744  void levelJacobi(LevelData<EBCellFAB>& a_phi,
745  const LevelData<EBCellFAB>& a_rhs,
746  int a_iterations);
747 
749  const LevelData<EBCellFAB>& a_rhs);
750 
752  const LevelData<EBCellFAB>& a_resid,
753  const IntVect& color);
754 
755  void colorGS(LevelData<EBCellFAB>& a_phi,
756  const LevelData<EBCellFAB>& a_rhs,
757  const int& a_icolor);
758 
759  void levelGSRB(LevelData<EBCellFAB>& a_phi,
760  const LevelData<EBCellFAB>& a_rhs);
761 
762  void levelGSRB(LevelData<EBCellFAB>& a_phi,
763  const LevelData<EBCellFAB>& a_rhs,
764  const int a_color);
765 
767  const LevelData<EBCellFAB>& a_rhs);
768 
770  const LevelData<EBCellFAB>& a_phiOld,
771  const LevelData<EBCellFAB>& a_rhs,
772  const int& a_icolor);
773 
774  static void doLazyRelax(bool a_doLazyRelax);
775  static void doEBEllipticLoadBalance(bool a_doEBEllipticLoadBalance);
776  static void areaFracWeighted(bool a_areaFracWeighted);
777 
778  static void getDivFStencil(VoFStencil& a_vofStencil,
779  const VolIndex& a_vof,
780  const EBISBox& a_ebisBox,
781  const IntVectSet& a_cfivs,
782  const RealVect& a_dx,
783  bool doFaceInterp = true);
784 
785  void getDivFStencil(VoFStencil& a_vofStencil,
786  const VolIndex& a_vof,
787  const DataIndex& a_dit,
788  bool doFaceInterp) const;
789 
790  //CP recently added
791  void getVoFStencil(LayoutData<BaseIVFAB<VoFStencil> > const*& a_vofStencil)
792  {
793  //I'm defining const* to make sure the content in m_opStencil is not modified
794  a_vofStencil = &m_opStencil;
795  }
796 
797  void getAlphaDiagWeight(LayoutData<BaseIVFAB<Real> > const*& a_alphaDiagWeight)
798  {
799  a_alphaDiagWeight = &m_alphaDiagWeight;
800  }
801 
802  void getAlphaBeta(Real& a_alpha, Real& a_beta)
803  {
804  a_alpha = m_alpha;
805  a_beta = m_beta;
806  }
807  //CP recently added
809  {
810  //I'm defining const* to make sure the content in m_opStencil is not modified
811  return m_domainBC;
812  }
813 
814  void setListValue(const LevelData<EBCellFAB>& a_data, Real a_value);
815 
816  void setListValue(EBCellFAB& a_data, const Vector<VolIndex>& a_setList, Real a_value);
817 
818  void setRhsSetList(const LayoutData<Vector<VolIndex> >& a_list);
819 
821  {
822  return m_rhsSetList;
823  }
824 
826  {
827  return m_rhsSetList;
828  }
829 
830  static void getFluxStencil(VoFStencil& a_fluxStencil,
831  const FaceIndex& a_face,
832  const EBISBox& a_ebisBox,
833  const IntVectSet& a_cfivs,
834  const RealVect& a_dx,
835  bool a_doFaceInterp);
836 
837  static void getFaceCenteredFluxStencil(VoFStencil& a_fluxStencil,
838  const FaceIndex& a_face,
839  const RealVect& a_dx);
840 
841  void applyCFBCs(LevelData<EBCellFAB>& a_phi,
842  const LevelData<EBCellFAB>* const a_phiCoarse,
843  bool a_homogeneousCFBC,
844  bool a_doOnlyRegularInterp = false);
845 
846  ///
847  /**
848  This function computes: a_lhs = (1/diagonal)*a_rhs
849  It is used to initialize the preconditioner, and by
850  MFPoissonOp::levelJacobi.
851  Consider using one of the level Gauss-Seidel methods
852  instead of monkeying with this.
853  */
855  const LevelData<EBCellFAB>& a_rhs);
856 
857  static int s_numComps;
858  static int s_whichComp;
859 
860 protected:
861 
863  void defineMGObjects(const EBLevelGrid& a_eblgCoarMG);
864  void defineWithCoarser(const EBLevelGrid& a_eblgCoar, const int& a_refToCoar);
865  void defineWithFiner(const EBLevelGrid& a_eblgFine,
866  const int& a_refToFine);
867 
868  static bool s_turnOffBCs;
870  static bool s_areaFracWeighted;
871  void defineStencils();
872  void defineEBCFStencils();
873 
875 
878 
880 
886 
890 
894 
902  static Real s_time;
903  static bool s_doLazyRelax;
905  static bool s_doTrimEdges;
906  static bool s_doSetListValueResid; // if this variable is set to true, it will apply setListValue in residual()
910  bool m_hasFine;
912  bool m_hasCoar;
915 
916  bool m_hasEBCF;
917 
919  //restriction object
921  //prolongation object
923 
924  //stencils for operator evaluation
933  //weights that get multiplied by alpha
935  //weights that get multiplied by beta
937  //constant for using in place of weight in EBStencil->apply for inhom dom bcs
939 
940  //cache the irreg vofiterator
943 
946 
951  // LayoutData<BaseIVFAB<Real> > m_cacheIrregFluxDomLo[EBAMRPO_NUMSTEN][SpaceDim];
952  // LayoutData<BaseIVFAB<Real> > m_cacheIrregFluxDomHi[EBAMRPO_NUMSTEN][SpaceDim];
953 
954  // Coarse-fine stencils for homogeneous CFInterp
957 
958  //flux register with finer level
960 
961  //stuff to make EBCF go faster
964 
965 
966  //special mg objects for when we do not have
967  //a coarser level or when the refinement to coarse
968  //is greater than two
969  //flag for when we need special MG objects
972 
974 
975 
976  //stuff below is only defined if m_hasMGObjects==true
982 
983 private:
984  //internally useful but not for general consumption
985  //lots of hidden assumptions and all that
986 
987  void fast_incrementFRCoar(const LevelData<EBCellFAB>& a_phiFine,
988  const LevelData<EBCellFAB>& a_phi);
989 
990  void fast_incrementFRFine(const LevelData<EBCellFAB>& a_phiFine,
991  const LevelData<EBCellFAB>& a_phi,
992  AMRLevelOp<LevelData<EBCellFAB> >* a_finerOp);
993 
994  ///this one gets called by base level tga (called by the public getFlux)
995  void getFlux(EBFaceFAB& a_flux,
996  const EBCellFAB& a_phi,
997  const Box& a_ghostedBox,
998  const Box& a_fabBox,
999  const ProblemDomain& a_domainBox,
1000  const EBISBox& a_ebisBox,
1001  const RealVect& a_dx,
1002  const int& a_idir);
1003 
1004  ///this one is internal (called by refluxing)
1005  void getFluxEBCF(EBFaceFAB& a_flux,
1006  const EBCellFAB& a_phi,
1007  const Box& a_ghostedBox,
1008  Vector<FaceIndex>& a_faceitEBCF,
1009  Vector<VoFStencil>& a_ebcfsten,
1010  const RealVect& a_dx);
1011 
1012  ///this one is internal (called by refluxing)
1013  void getFluxRegO(EBFaceFAB& a_flux,
1014  const EBCellFAB& a_phi,
1015  const Box& a_ghostedBox,
1016  const RealVect& a_dx);
1017 
1018  void getOpVoFStencil(VoFStencil& a_stencil,
1019  const EBISBox& a_ebisbox,
1020  const VolIndex& a_vof);
1021 
1022  void getOpVoFStencil(VoFStencil& a_stencil,
1023  const int& a_dir,
1024  const Vector<VolIndex>& a_allMonotoneVoFs,
1025  const EBISBox& a_ebisbox,
1026  const VolIndex& a_vof,
1027  const bool& a_lowOrder);
1028 
1029 
1030  void getOpFaceStencil(VoFStencil& a_stencil,
1031  const Vector<VolIndex>& a_allMonotoneVofs,
1032  const EBISBox& a_ebisbox,
1033  const VolIndex& a_vof,
1034  int a_dir,
1035  const Side::LoHiSide& a_side,
1036  const FaceIndex& a_face,
1037  const bool& a_lowOrder);
1038 
1039  void levelJacobi(LevelData<EBCellFAB>& a_phi,
1040  const LevelData<EBCellFAB>& a_rhs);
1041 
1043 
1044  void applyHomogeneousCFBCs(EBCellFAB& a_phi,
1045  const DataIndex& a_datInd,
1046  int a_idir,
1047  Side::LoHiSide a_hiorlo);
1048 
1049  Real getRadius(const FaceIndex& a_face, const RealVect& a_centroid);
1050 
1051  ///
1052  /** applyOp() on the regular cells for all directions
1053  opPhi comes in holding alpha*phi. this adds on beta*lapl phi*/
1054  void applyOpRegularAllDirs(Box * a_loBox,
1055  Box * a_hiBox,
1056  int * a_hasLo,
1057  int * a_hasHi,
1058  Box & a_curOpPhiBox,
1059  Box & a_curPhiBox,
1060  int a_nComps,
1061  BaseFab<Real> & a_curOpPhiFAB,
1062  const BaseFab<Real> & a_curPhiFAB,
1063  bool a_homogeneousPhysBC,
1064  const DataIndex& a_dit,
1065  const Real& a_beta);
1066 
1067  void applyDomainFlux(Box * a_loBox,
1068  Box * a_hiBox,
1069  int * a_hasLo,
1070  int * a_hasHi,
1071  Box & a_curPhiBox,
1072  int a_nComps,
1073  BaseFab<Real> & a_phiFAB,
1074  bool a_homogeneousPhysBC,
1075  const DataIndex& a_dit,
1076  const Real& a_beta);
1077 
1078 private:
1079  //copy constructor and operator= disallowed for all the usual reasons
1081  {
1082  MayDay::Error("invalid operator");
1083  }
1084 
1085  void operator=(const EBAMRPoissonOp& a_opin)
1086  {
1087  MayDay::Error("invalid operator");
1088  }
1089 };
1090 
1091 
1092 #include "NamespaceFooter.H"
1093 #endif
bool m_hasCoar
Definition: EBAMRPoissonOp.H:912
virtual void assignCopier(LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs, const Copier &a_copier)
EBMGAverage m_ebAverage
Definition: EBAMRPoissonOp.H:920
Vector< IntVect > m_colors
Definition: EBAMRPoissonOp.H:973
void levelJacobi(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs, int a_iterations)
virtual void AMRResidual(LevelData< EBCellFAB > &a_residual, const LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_phiCoarse, const LevelData< EBCellFAB > &a_rhs, bool a_homogeneousBC, AMRLevelOp< LevelData< EBCellFAB > > *a_finerOp)
LayoutData< VoFIterator > m_vofItIrregDomHi[SpaceDim]
Definition: EBAMRPoissonOp.H:945
virtual void relax(LevelData< EBCellFAB > &a_e, const LevelData< EBCellFAB > &a_residual, int a_iterations)
int m_testRef
Definition: EBAMRPoissonOp.H:874
virtual void create(LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs)
RealVect m_dxFine
Definition: EBAMRPoissonOp.H:891
virtual void AMRResidualNF(LevelData< EBCellFAB > &a_residual, const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_phiCoarse, const LevelData< EBCellFAB > &a_rhs, bool a_homogeneousBC)
static void setOperatorTime(Real a_time)
Definition: EBAMRPoissonOp.H:324
EBFastFR-A class to encapsulate a levels worth of flux registers.
Definition: EBFastFR.H:39
static int getMatrixIndexingLD(LevelData< BaseEBCellFAB< int > > &a_gids, int &a_data, const EBLevelGrid &a_eblg, const IntVect &a_ghostCellsPhi, const bool &a_hasCoar)
RealVect m_origin
Definition: EBAMRPoissonOp.H:907
EBISLayout getEBISL() const
Definition: EBLevelGrid.H:93
An irregular domain on an integer lattice.
Definition: IntVectSet.H:44
virtual void prolongIncrement(LevelData< EBCellFAB > &a_phiThisLevel, const LevelData< EBCellFAB > &a_correctCoarse)
bool m_hasFine
Definition: EBAMRPoissonOp.H:910
LayoutData< RefCountedPtr< EBSTENCIL_T > > m_opEBStencilInhomDomHi[SpaceDim]
Definition: EBAMRPoissonOp.H:928
static Real staticMaxNorm(const LevelData< EBCellFAB > &a_rhs, const EBLevelGrid &a_eblg)
A class to facilitate interaction with physical boundary conditions.
Definition: ProblemDomain.H:141
virtual void setToZero(LevelData< EBCellFAB > &a_lhs)
virtual void createCoarsened(LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs, const int &a_refRat)
void defineWithFiner(const EBLevelGrid &a_eblgFine, const int &a_refToFine)
LayoutData< BaseIVFAB< Real > > m_betaDiagWeight
Definition: EBAMRPoissonOp.H:936
Definition: EBIndexSpace.H:50
virtual void kappaScale(LevelData< EBCellFAB > &a_rhs)
for eb only. kappa weight the rhs but do not multiply by the identity coefficient ...
Definition: EBAMRPoissonOp.H:149
void getDomainFluxStencil(VoFStencil &a_stencil, const VolIndex &a_vof, const int a_comp, const DataIndex &a_dit)
EBFastFR m_fastFR
Definition: EBAMRPoissonOp.H:959
void defineWithCoarser(const EBLevelGrid &a_eblgCoar, const int &a_refToCoar)
Real getRadius(const FaceIndex &a_face, const RealVect &a_centroid)
virtual void GSColorAllIrregular(EBCellFAB &a_phi, const EBCellFAB &a_rhs, const int &a_icolor, const bool &a_homogeneousPhysBC, const DataIndex &a_dit)
void getFluxEBCF(EBFaceFAB &a_flux, const EBCellFAB &a_phi, const Box &a_ghostedBox, Vector< FaceIndex > &a_faceitEBCF, Vector< VoFStencil > &a_ebcfsten, const RealVect &a_dx)
this one is internal (called by refluxing)
void mvApplyDomainFlux(BaseFab< Real > &a_phiFAB, const Box &a_grid, const DataIndex &a_dit)
apply domainflux in multivariable mode
one dimensional dynamic array
Definition: Vector.H:53
virtual Real localMaxNorm(const LevelData< EBCellFAB > &a_rhs)
VolIndex m_vof
Definition: EBAMRPoissonOp.H:209
Definition: FaceIndex.H:28
Data that maintains a one-to-one mapping of T to the boxes in a BoxLayout.
Definition: BoxLayout.H:26
const ProblemDomain & getDomain() const
Definition: EBLevelGrid.H:132
virtual void restrictResidual(LevelData< EBCellFAB > &a_resCoarse, LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_rhsFine)
LayoutData< RefCountedPtr< EBSTENCIL_T > > m_colorEBStencil[EBAMRPO_NUMSTEN]
Definition: EBAMRPoissonOp.H:929
Definition: EBAMRPoissonOp.H:213
int m_refToCoar
Definition: EBAMRPoissonOp.H:909
void applyOpNoCFBCs(LevelData< EBCellFAB > &a_opPhi, const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > *const a_phiCoar, DataIterator &a_dit, const bool &a_homogeneousPhysBC, const bool &a_homogeneousCFBC, const LevelData< BaseIVFAB< Real > > *const a_ebFluxBCLD)
no exchange of cf interp
LayoutData< Vector< VolIndex > > & getRhsSetList()
Definition: EBAMRPoissonOp.H:820
A strange but true thing to make copying from one boxlayoutdata to another fast.
Definition: Copier.H:145
int getLevelDataFromPetscVector(LevelData< EBCellFAB > &a_data, const Vec &a_petsc_vec_real, const Vec &a_petsc_vec_imag)
put the real components into comp 0, the imaginary ones into comp 1
const IntVect m_ghostCellsRHS
Definition: EBAMRPoissonOp.H:877
Definition: EBISBox.H:46
static void doEBEllipticLoadBalance(bool a_doEBEllipticLoadBalance)
void writeEBAMRname(const Vector< LevelData< EBCellFAB > * > *a_dataPtr, const char *a_filename)
static int s_whichComp
Definition: EBAMRPoissonOp.H:858
int getPetscMatrix(Mat &a_petsc_mat)
Fill a_petsc_mat with a matrix that describes the operator stencil.
const IntVect & gridIndex() const
Definition: VolIndex.H:125
bool operator!=(const StencilIndex &a_sin) const
Definition: EBAMRPoissonOp.H:203
LayoutData< CFIVS > m_loCFIVS[SpaceDim]
Definition: EBAMRPoissonOp.H:955
LayoutData< BaseIVFAB< VoFStencil > > m_opStencil
Definition: EBAMRPoissonOp.H:925
Definition: EBLevelGrid.H:30
static bool s_doSetListValueResid
Definition: EBAMRPoissonOp.H:906
virtual Real dotProduct(const LevelData< EBCellFAB > &a_1, const LevelData< EBCellFAB > &a_2)
LayoutData< VoFIterator > m_vofItIrregColor[EBAMRPO_NUMSTEN]
Definition: EBAMRPoissonOp.H:942
EBISLayout m_ebislCoarMG
Definition: EBAMRPoissonOp.H:980
EBLevelGrid getEBLG()
Definition: EBAMRPoissonOp.H:316
static bool s_doInconsistentRelax
Definition: EBAMRPoissonOp.H:904
Definition: DataIterator.H:190
virtual void dumpAMR(Vector< LevelData< EBCellFAB > *> &a_data, string name)
Definition: EBAMRPoissonOp.H:86
const LayoutData< Vector< VolIndex > > & getListValue()
Definition: EBAMRPoissonOp.H:825
void writeEBLevelname(const LevelData< EBCellFAB > *a_dataPtr, const char *a_filename)
LevelData< EBCellFAB > m_resThisLevel
Definition: EBAMRPoissonOp.H:862
static bool s_turnOffBCs
Definition: EBAMRPoissonOp.H:868
int m_relaxType
Definition: EBAMRPoissonOp.H:914
void slowGSRBColor(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_res, int a_iredblack)
void mvBetaLaplacianGrid(EBCellFAB &a_lph, const EBCellFAB &a_phi, const DataIndex &a_dit)
virtual void applyOp(LevelData< EBCellFAB > &a_opPhi, const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > *const a_phiCoarse, const bool &a_homogeneousPhysBC, const bool &a_homogeneousCFBC)
Definition: EBFaceFAB.H:28
EBMGInterp m_ebInterp
Definition: EBAMRPoissonOp.H:922
void reflux(LevelData< EBCellFAB > &a_residual, const LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_phi, AMRLevelOp< LevelData< EBCellFAB > > *a_finerOp)
virtual void AMRRestrict(LevelData< EBCellFAB > &a_resCoarse, const LevelData< EBCellFAB > &a_residual, const LevelData< EBCellFAB > &a_correction, const LevelData< EBCellFAB > &a_coarseCorrection, bool a_skip_res=false)
void levelMultiColorGS(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs)
Copier m_exchangeCopier
Definition: EBAMRPoissonOp.H:918
virtual void GSColorAllRegularClone(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_phiOld, const LevelData< EBCellFAB > &a_rhs, const int &a_icolor, const Real &a_weight, const bool &a_homogeneousPhysBC)
void fast_incrementFRFine(const LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_phi, AMRLevelOp< LevelData< EBCellFAB > > *a_finerOp)
void setListValue(const LevelData< EBCellFAB > &a_data, Real a_value)
Real AMRNorm(const LevelData< EBCellFAB > &a_coarResid, const LevelData< EBCellFAB > &a_fineResid, const int &a_refRat, const int &a_ord)
int m_phase
Definition: EBAMRPoissonOp.H:210
Piecewise constant interpolation.
Definition: EBMGInterp.H:33
void operator=(const EBAMRPoissonOp &a_opin)
Definition: EBAMRPoissonOp.H:1085
int m_refToFine
Definition: EBAMRPoissonOp.H:908
const int SpaceDim
Definition: SPACE.H:38
static bool s_areaFracWeighted
Definition: EBAMRPoissonOp.H:870
virtual void scale(LevelData< EBCellFAB > &a_lhs, const Real &a_scale)
Definition: AMRMultiGrid.H:39
static Real s_time
Definition: EBAMRPoissonOp.H:902
virtual int refToCoarser()
virtual void assign(LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs)
VoF-centered stencil.
Definition: Stencils.H:60
RealVect m_invDx
Definition: EBAMRPoissonOp.H:895
RealVect m_invDx2
Definition: EBAMRPoissonOp.H:896
RealVect m_dx
Definition: EBAMRPoissonOp.H:892
void getAlphaDiagWeight(LayoutData< BaseIVFAB< Real > > const *&a_alphaDiagWeight)
Definition: EBAMRPoissonOp.H:797
Definition: EBAMRPoissonOp.H:76
A EBFaceFAB-like container for edge-centered fluxes.
Definition: EBFluxFAB.H:25
virtual void assignLocal(LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs)
copier definition was killing us.
void getOpFaceStencil(VoFStencil &a_stencil, const Vector< VolIndex > &a_allMonotoneVofs, const EBISBox &a_ebisbox, const VolIndex &a_vof, int a_dir, const Side::LoHiSide &a_side, const FaceIndex &a_face, const bool &a_lowOrder)
EBLevelGrid m_eblgCoarMG
Definition: EBAMRPoissonOp.H:884
virtual ~EBAMRPoissonOp()
Real m_dxScale
Definition: EBAMRPoissonOp.H:897
void getInvDiagRHS(LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs)
static void getFluxStencil(VoFStencil &a_fluxStencil, const FaceIndex &a_face, const EBISBox &a_ebisBox, const IntVectSet &a_cfivs, const RealVect &a_dx, bool a_doFaceInterp)
static bool getCoarserLayouts(DisjointBoxLayout &a_dblCoar, ProblemDomain &a_domainCoar, const DisjointBoxLayout &a_dblFine, const EBISLayout &a_ebislFine, const ProblemDomain &a_domainFine, int a_refToCoar, const EBIndexSpace *a_ebisPtr, int a_maxBoxSize, bool &a_layoutChanged, int a_testRef=2)
virtual int refToFiner()
Vector< Vector< RefCountedPtr< LayoutData< EBCellFAB > > > > m_cacheInhomDomBCLo
Definition: EBAMRPoissonOp.H:949
virtual void AMRProlong(LevelData< EBCellFAB > &a_correction, const LevelData< EBCellFAB > &a_coarseCorrection)
void AMRResidualNC(LevelData< EBCellFAB > &a_residual, const LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs, bool a_homogeneousBC, AMRLevelOp< LevelData< EBCellFAB > > *a_finerOp)
virtual void buildCopier(Copier &a_copier, const LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_rhs)
StencilIndex()
Definition: EBAMRPoissonOp.H:181
const char * name(const FArrayBox &a_dummySpecializationArg)
Definition: CH_HDF5.H:907
void getAlphaBeta(Real &a_alpha, Real &a_beta)
Definition: EBAMRPoissonOp.H:802
Definition: EBCellFAB.H:29
Real m_beta
Definition: EBAMRPoissonOp.H:900
LayoutData< VoFIterator > m_vofItIrreg
Definition: EBAMRPoissonOp.H:941
static int s_numComps
Definition: EBAMRPoissonOp.H:857
void dumpStencilMatrix()
void fast_incrementFRCoar(const LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_phi)
void defineEBCFStencils()
virtual void setTime(Real a_time)
Definition: EBAMRPoissonOp.H:330
EBMGInterp m_ebInterpMG
Definition: EBAMRPoissonOp.H:978
double Real
Definition: REAL.H:33
EBAMRPoissonOp(const EBAMRPoissonOp &a_opin)
Definition: EBAMRPoissonOp.H:1080
LayoutData< RefCountedPtr< EBSTENCIL_T > > m_opEBStencilInhomDomLo[SpaceDim]
Definition: EBAMRPoissonOp.H:927
static void doLazyRelax(bool a_doLazyRelax)
EBLevelGrid m_eblgCoar
Definition: EBAMRPoissonOp.H:883
bool m_hasMGObjects
Definition: EBAMRPoissonOp.H:970
void defineStencils()
Real m_bCoef
Definition: EBAMRPoissonOp.H:901
A BoxLayout that has a concept of disjointedness.
Definition: DisjointBoxLayout.H:30
LayoutData< VoFIterator > m_vofItIrregColorDomHi[EBAMRPO_NUMSTEN][SpaceDim]
Definition: EBAMRPoissonOp.H:948
Array defined at the VolIndexs of an Box in an EBIS.
Definition: BaseEBCellFAB.H:40
RefCountedPtr< BaseDomainBC > m_domainBC
Definition: EBAMRPoissonOp.H:887
LoHiSide
Definition: LoHiSide.H:27
static bool s_doEBEllipticLoadBalance
Definition: EBAMRPoissonOp.H:869
static void areaFracWeighted(bool a_areaFracWeighted)
void defineMGObjects(const EBLevelGrid &a_eblgCoarMG)
virtual void AMROperator(LevelData< EBCellFAB > &a_LofPhi, const LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_phiCoarse, bool a_homogeneousBC, AMRLevelOp< LevelData< EBCellFAB > > *a_finerOp)
virtual void AMROperatorNF(LevelData< EBCellFAB > &a_LofPhi, const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_phiCoarse, bool a_homogeneousBC)
EBLevelGrid getEBLGCoarMG()
Definition: EBAMRPoissonOp.H:320
LayoutData< RefCountedPtr< EBSTENCIL_T > > m_opEBStencil
Definition: EBAMRPoissonOp.H:926
virtual void fillGrad(const LevelData< EBCellFAB > &a_phi)
a leveltgaism
Definition: EBAMRPoissonOp.H:164
RefCountedPtr< EBQuadCFInterp > m_quadCFIWithCoar
Definition: EBAMRPoissonOp.H:879
LayoutData< BaseIVFAB< Real > > m_alphaDiagWeight
Definition: EBAMRPoissonOp.H:934
bool m_hasInterpAve
Definition: EBAMRPoissonOp.H:911
static bool s_doTrimEdges
Definition: EBAMRPoissonOp.H:905
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 Real norm(const LevelData< EBCellFAB > &a_rhs, int a_ord)
Vector< Vector< RefCountedPtr< LayoutData< EBCellFAB > > > > m_cacheInhomDomBCHi
Definition: EBAMRPoissonOp.H:950
virtual void applyOpNoBoundary(LevelData< EBCellFAB > &a_opPhi, const LevelData< EBCellFAB > &a_phi)
virtual function called by LevelTGA
virtual void preCond(LevelData< EBCellFAB > &a_opPhi, const LevelData< EBCellFAB > &a_phi)
LayoutData< Vector< FaceIndex > > m_faceitCoar[2 *SpaceDim]
Definition: EBAMRPoissonOp.H:962
RefCountedPtr< BaseEBBC > m_ebBC
Definition: EBAMRPoissonOp.H:888
void levelGSRB(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs)
virtual void setAlphaAndBeta(const Real &a_alpha, const Real &a_beta)
void setEBBC(const RefCountedPtr< BaseEBBC > &a_ebBC)
virtual void setVal(LevelData< EBCellFAB > &a_lhs, const Real &a_value)
LayoutData< RefCountedPtr< EBSTENCIL_T > > m_invDiagEBStencil
Definition: EBAMRPoissonOp.H:932
bool m_layoutChanged
Definition: EBAMRPoissonOp.H:971
Real m_alpha
Definition: EBAMRPoissonOp.H:898
A Rectangular Domain on an Integer Lattice.
Definition: Box.H:465
A Real vector in SpaceDim-dimensional space.
Definition: RealVect.H:41
const RefCountedPtr< BaseDomainBC > getDomainBC()
Definition: EBAMRPoissonOp.H:808
void getOpVoFStencil(VoFStencil &a_stencil, const EBISBox &a_ebisbox, const VolIndex &a_vof)
void AMROperatorNC(LevelData< EBCellFAB > &a_LofPhi, const LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_phi, bool a_homogeneousBC, AMRLevelOp< LevelData< EBCellFAB > > *a_finerOp)
LayoutData< Vector< VoFStencil > > m_stencilCoar[2 *SpaceDim]
Definition: EBAMRPoissonOp.H:963
void fast_reflux(LevelData< EBCellFAB > &a_residual, const LevelData< EBCellFAB > &a_phiFine, const LevelData< EBCellFAB > &a_phi, AMRLevelOp< LevelData< EBCellFAB > > *a_finerOp)
virtual void AMRUpdateResidual(LevelData< EBCellFAB > &a_residual, const LevelData< EBCellFAB > &a_correction, const LevelData< EBCellFAB > &a_coarseCorrection)
Definition: DataIndex.H:112
bool m_hasEBCF
Definition: EBAMRPoissonOp.H:916
Real dx() const
returns m_dx, such function is required by some LinearSolvers
Definition: EBAMRPoissonOp.H:158
int cellIndex() const
Definition: VolIndex.H:133
LayoutData< VoFIterator > m_vofItIrregDomLo[SpaceDim]
Definition: EBAMRPoissonOp.H:944
Piecewise constant interpolation.
Definition: EBMGAverage.H:31
Real m_aCoef
Definition: EBAMRPoissonOp.H:899
StencilIndex & operator=(const StencilIndex &a_sin)
Definition: EBAMRPoissonOp.H:186
void dumpReferenceStencilMatrix()
EBLevelGrid m_eblg
Definition: EBAMRPoissonOp.H:881
void colorGS(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs, const int &a_icolor)
LayoutData< Vector< VolIndex > > m_rhsSetList
Definition: EBAMRPoissonOp.H:889
LayoutData< VoFIterator > m_vofItIrregColorDomLo[EBAMRPO_NUMSTEN][SpaceDim]
Definition: EBAMRPoissonOp.H:947
An integer Vector in SpaceDim-dimensional space.
Definition: CHArray.H:42
ProblemDomain m_domainCoarMG
Definition: EBAMRPoissonOp.H:981
Definition: AMRTGA.H:159
void colorGSClone(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_phiOld, const LevelData< EBCellFAB > &a_rhs, const int &a_icolor)
LayoutData< CFIVS > m_hiCFIVS[SpaceDim]
Definition: EBAMRPoissonOp.H:956
void getVoFStencil(VoFStencil &a_stenc, const VolIndex &a_vof, const DataIndex &a_dit) const
LayoutData< RefCountedPtr< EBSTENCIL_T > > m_colorEBStencilDomHi[EBAMRPO_NUMSTEN][SpaceDim]
Definition: EBAMRPoissonOp.H:931
void THIS_IS_AN_ERROR_MESSAGE(void)
Definition: EBAMRPoissonOp.H:66
Volume of Fluid Index.
Definition: VolIndex.H:31
void applyCFBCs(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > *const a_phiCoarse, bool a_homogeneousCFBC, bool a_doOnlyRegularInterp=false)
void levelMultiColorGSClone(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs)
Box & grow(int i)
grow functions
Definition: Box.H:2247
EBMGAverage m_ebAverageMG
Definition: EBAMRPoissonOp.H:977
EBLevelGrid m_eblgCoarsenedFine
Definition: EBAMRPoissonOp.H:885
Definition: EBISLayout.H:39
virtual void getFlux(EBFluxFAB &a_flux, const LevelData< EBCellFAB > &a_data, const Box &a_grid, const DataIndex &a_dit, Real a_scale)
Definition: EBAMRPoissonOp.H:298
virtual void diagonalScale(LevelData< EBCellFAB > &a_rhs, bool a_kappaWeighted=true)
static bool s_doLazyRelax
Definition: EBAMRPoissonOp.H:903
void applyOpRegularAllDirs(Box *a_loBox, Box *a_hiBox, int *a_hasLo, int *a_hasHi, Box &a_curOpPhiBox, Box &a_curPhiBox, int a_nComps, BaseFab< Real > &a_curOpPhiFAB, const BaseFab< Real > &a_curPhiFAB, bool a_homogeneousPhysBC, const DataIndex &a_dit, const Real &a_beta)
void levelSlowRelax(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs)
virtual void residual(LevelData< EBCellFAB > &a_residual, const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs, bool a_homogeneousPhysBC=false)
StencilIndex(VolIndex &a_vof, int &a_phase)
Definition: EBAMRPoissonOp.H:174
const IntVect m_ghostCellsPhi
Definition: EBAMRPoissonOp.H:876
virtual void axby(LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_x, const LevelData< EBCellFAB > &a_y, Real a_a, Real a_b)
virtual void GSColorAllIrregularClone(LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_phiOld, const LevelData< EBCellFAB > &a_rhs, const int &a_icolor, const bool &a_homogeneousPhysBC)
virtual void GSColorAllRegular(BaseFab< Real > &a_phi, const BaseFab< Real > &a_rhs, const int &a_icolor, const Real &a_weight, const bool &a_homogeneousPhysBC, const DataIndex &a_dit)
void getFluxRegO(EBFaceFAB &a_flux, const EBCellFAB &a_phi, const Box &a_ghostedBox, const RealVect &a_dx)
this one is internal (called by refluxing)
void applyDomainFlux(Box *a_loBox, Box *a_hiBox, int *a_hasLo, int *a_hasHi, Box &a_curPhiBox, int a_nComps, BaseFab< Real > &a_phiFAB, bool a_homogeneousPhysBC, const DataIndex &a_dit, const Real &a_beta)
virtual void dumpLevel(LevelData< EBCellFAB > &a_data, string name)
Definition: EBAMRPoissonOp.H:91
VolIndex vof() const
Definition: EBAMRPoissonOp.H:193
static void getFaceCenteredFluxStencil(VoFStencil &a_fluxStencil, const FaceIndex &a_face, const RealVect &a_dx)
virtual void incr(LevelData< EBCellFAB > &a_lhs, const LevelData< EBCellFAB > &a_x, Real a_scale)
LayoutData< RefCountedPtr< EBSTENCIL_T > > m_colorEBStencilDomLo[EBAMRPO_NUMSTEN][SpaceDim]
Definition: EBAMRPoissonOp.H:930
int m_numPreCondIters
Definition: EBAMRPoissonOp.H:913
static void getAggregatedLayout(DisjointBoxLayout &a_dblCoar, const ProblemDomain &a_domainCoar, const EBIndexSpace *const a_ebisPtr, const int &a_maxBoxSize)
LayoutData< BaseIVFAB< Real > > m_one
Definition: EBAMRPoissonOp.H:938
virtual void getOpMatrix(const LevelData< EBCellFAB > &a_phi, const LevelData< EBCellFAB > &a_rhs)
void getVoFStencil(LayoutData< BaseIVFAB< VoFStencil > > const *&a_vofStencil)
Definition: EBAMRPoissonOp.H:791
EBLevelGrid m_eblgFine
Definition: EBAMRPoissonOp.H:882
void setRhsSetList(const LayoutData< Vector< VolIndex > > &a_list)
Definition: EBAMRPoissonOp.H:171
virtual void divideByIdentityCoef(LevelData< EBCellFAB > &a_rhs)
Definition: EBAMRPoissonOp.H:144
static void getDivFStencil(VoFStencil &a_vofStencil, const VolIndex &a_vof, const EBISBox &a_ebisBox, const IntVectSet &a_cfivs, const RealVect &a_dx, bool doFaceInterp=true)
DisjointBoxLayout m_dblCoarMG
Definition: EBAMRPoissonOp.H:979
RealVect m_dxCoar
Definition: EBAMRPoissonOp.H:893
void applyHomogeneousCFBCs(LevelData< EBCellFAB > &a_phi)
int phase() const
Definition: EBAMRPoissonOp.H:198
virtual void createCoarser(LevelData< EBCellFAB > &a_coarse, const LevelData< EBCellFAB > &a_fine, bool a_ghosted)