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