Chombo + EB + MF  3.2
ViscousTensorOp.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 _VISCOUSTENSOROP_H_
12 #define _VISCOUSTENSOROP_H_
13 
14 #include "AMRMultiGrid.H"
15 #include "REAL.H"
16 #include "Box.H"
17 #include "LevelDataOps.H"
18 #include "BCFunc.H"
19 #include "FArrayBox.H"
20 #include "FluxBox.H"
21 #include "CFIVS.H"
22 #include "AMRTGA.H"
23 #include "TensorCFInterp.H"
24 #include "CoarseAverage.H"
25 #include "LevelFluxRegister.H"
26 #include "AMRIO.H"
27 #include "NamespaceHeader.H"
28 
29 #define VTOP_DEFAULT_SAFETY 0.9
30 ///
31 /**
32  operator is alpha a I + (divF) = alpha*acoef I + beta*div(eta(grad B + grad BT) + lambda I div B )
33  beta, lambda, and eta are incorporated into the flux F
34 */
35 class ViscousTensorOp : public LevelTGAHelmOp<LevelData<FArrayBox>, FluxBox>
36 {
37 public:
38 
40  {
44  };
45 
46  virtual void diagonalScale(LevelData<FArrayBox>& a_rhs, bool a_kappaWeighted)
47  {
48  diagonalScale(a_rhs);
49  }
50 
51  virtual void setAlphaAndBeta(const Real& a_alpha,
52  const Real& a_beta)
53  {
54  m_alpha = a_alpha;
55  m_beta = a_beta;
56  defineRelCoef();
57  }
58 
59  virtual void diagonalScale(LevelData<FArrayBox>& a_rhs)
60  {
61  DisjointBoxLayout grids = a_rhs.disjointBoxLayout();
62  for (DataIterator dit = grids.dataIterator(); dit.ok(); ++dit)
63  {
64  for (int idir = 0; idir < SpaceDim; idir++)
65  {
66  int isrc = 0; int idst = idir; int inco = 1;
67  a_rhs[dit()].mult((*m_acoef)[dit()], isrc, idst, inco);
68  }
69  }
70  }
71 
73  {
74  DisjointBoxLayout grids = a_rhs.disjointBoxLayout();
75  for (DataIterator dit = grids.dataIterator(); dit.ok(); ++dit)
76  {
77  for (int idir = 0; idir < SpaceDim; idir++)
78  {
79  int isrc = 0; int idst = idir; int inco = 1;
80  a_rhs[dit()].divide((*m_acoef)[dit()], isrc, idst, inco);
81  }
82  }
83  }
84 
85  ///
86  /**
87  */
88  virtual ~ViscousTensorOp()
89  {
90 #if 0
91  // these explicitly de-reference the storage
96 #endif
97  }
98 
99  ///
100  /**
101  */
102  ViscousTensorOp(const DisjointBoxLayout& a_grids,
103  const DisjointBoxLayout& a_gridsFine,
104  const DisjointBoxLayout& a_gridsCoar,
105  const RefCountedPtr<LevelData<FluxBox> >& a_eta,
106  const RefCountedPtr<LevelData<FluxBox> >& a_lambda,
107  const RefCountedPtr<LevelData<FArrayBox> >& a_acoef,
108  Real a_alpha,
109  Real a_beta,
110  int a_refToFine,
111  int a_refToCoar,
112  const ProblemDomain& a_domain,
113  const Real& a_dxLevel,
114  const Real& a_dxCoar,
115  BCFunc a_bc,
116  Real a_safety = VTOP_DEFAULT_SAFETY,
117  Real a_relaxTolerance = 0.0,
118  int a_relaxMinIter = 1000
119  );
120 
121  ///
122  /**
123  */
124  ViscousTensorOp(const DisjointBoxLayout& a_grids,
125  const DisjointBoxLayout& a_gridsFine,
126  const DisjointBoxLayout& a_gridsCoar,
127  const RefCountedPtr<LevelData<FluxBox> >& a_eta,
128  const RefCountedPtr<LevelData<FluxBox> >& a_lambda,
129  const RefCountedPtr<LevelData<FArrayBox> >& a_acoef,
130  Real a_alpha,
131  Real a_beta,
132  int a_refToFine,
133  int a_refToCoar,
134  const ProblemDomain& a_domain,
135  const Real& a_dxLevel,
136  const Real& a_dxCoar,
137  BCHolder& a_bc,
138  Real a_safety = VTOP_DEFAULT_SAFETY,
139  Real a_relaxTolerance = 0.0,
140  int a_relaxMinIter = 1000
141  );
142 
143  ///
144  /**
145  */
146  void define(const DisjointBoxLayout& a_grids,
147  const DisjointBoxLayout& a_gridsFine,
148  const DisjointBoxLayout& a_gridsCoar,
149  const RefCountedPtr<LevelData<FluxBox> >& a_eta,
150  const RefCountedPtr<LevelData<FluxBox> >& a_lambda,
151  const RefCountedPtr<LevelData<FArrayBox> >& a_acoef,
152  Real a_alpha,
153  Real a_beta,
154  int a_refToFine,
155  int a_refToCoar,
156  const ProblemDomain& a_domain,
157  const Real& a_dxLevel,
158  const Real& a_dxCoar,
159  BCHolder& a_bc,
160  Real a_safety = VTOP_DEFAULT_SAFETY,
161  Real a_relaxTolerance = 0.0,
162  int a_relaxMinIter = 1000
163  );
164 
165  virtual void residual( LevelData<FArrayBox>& a_lhs,
166  const LevelData<FArrayBox>& a_phi,
167  const LevelData<FArrayBox>& a_rhs,
168  bool a_homogeneous = false);
169 
170 
171  virtual void residualNoExchange( LevelData<FArrayBox>& a_lhs,
172  const LevelData<FArrayBox>& a_phi,
173  const LevelData<FArrayBox>& a_rhs,
174  bool a_homogeneous = false);
175 
176  virtual void preCond( LevelData<FArrayBox>& a_correction,
177  const LevelData<FArrayBox>& a_residual);
178 
179  virtual void applyOp( LevelData<FArrayBox>& a_lhs,
180  const LevelData<FArrayBox>& a_phi,
181  bool a_homogeneous = false);
182 
183  // version of applyOp which doesn't call exchange
184  virtual void applyOpNoExchange( LevelData<FArrayBox>& a_lhs,
185  const LevelData<FArrayBox>& a_phi,
186  bool a_homogeneous = false);
187 
188  virtual void create( LevelData<FArrayBox>& a_lhs,
189  const LevelData<FArrayBox>& a_rhs);
190  virtual void createCoarsened( LevelData<FArrayBox>& a_lhs,
191  const LevelData<FArrayBox>& a_rhs,
192  const int& a_refRat);
193 
194  void reflux(const LevelData<FArrayBox>& a_phiFine,
195  const LevelData<FArrayBox>& a_phi,
196  LevelData<FArrayBox>& residual,
197  AMRLevelOp<LevelData<FArrayBox> >* a_finerOp);
198 
199  //use at thy peril
200  virtual void getFlux(FluxBox& a_flux,
201  const LevelData<FArrayBox>& a_data,
202  const Box& a_grid,
203  const DataIndex& a_dit,
204  Real a_scale)
205  {
206  const FArrayBox& data = a_data[a_dit];
207  a_flux.define(a_grid, SpaceDim);
208  for (int idir = 0; idir < SpaceDim; idir++)
209  {
210  const FArrayBox& etaFace = (*m_eta)[a_dit][idir];
211  const FArrayBox& lambdaFace = (*m_lambda)[a_dit][idir];
212  const FArrayBox& gradData = m_grad[a_dit];
213  Box faceBox = a_flux[idir].box();
214  Box domFaceBox = surroundingNodes(m_domain.domainBox(), idir);
215  faceBox &= domFaceBox;
216  getFlux(a_flux[idir], data, gradData, etaFace, lambdaFace, faceBox, idir, 1);
217  a_flux[idir] *= a_scale;
218  }
219  }
221  const LevelData<FArrayBox>& a_phi)
222  {
223  this->fillGrad(a_phi);
224  computeOperatorNoBCs(a_lhs, a_phi);
225  }
226 
227  void getFlux(FArrayBox& a_flux,
228  const FArrayBox& a_data,
229  const FArrayBox& a_gradData,
230  const FArrayBox& a_etaFace,
231  const FArrayBox& a_lambdaFace,
232  const Box& a_facebox,
233  int a_dir,
234  int ref = 1);
235 
236  /// utility function which computes operator after all bc's have been set
238  const LevelData<FArrayBox>& a_phi);
239 
240  virtual void assign( LevelData<FArrayBox>& a_lhs,
241  const LevelData<FArrayBox>& a_rhs) ;
242  virtual Real dotProduct(const LevelData<FArrayBox>& a_1,
243  const LevelData<FArrayBox>& a_2) ;
244  virtual void incr( LevelData<FArrayBox>& a_lhs,
245  const LevelData<FArrayBox>& a_x,
246  Real a_scale) ;
247  virtual void axby( LevelData<FArrayBox>& a_lhs,
248  const LevelData<FArrayBox>& a_x,
249  const LevelData<FArrayBox>& a_y,
250  Real a, Real b) ;
251 
252  virtual void scale(LevelData<FArrayBox>& a_lhs, const Real& a_scale) ;
253 
254  virtual Real norm(const LevelData<FArrayBox>& a_x, int a_ord);
255 
256  virtual Real dx() const
257  {
258  return m_dx;
259  }
260 
261  virtual Real dxCrse() const
262  {
263  return m_dxCrse;
264  }
265 
266  virtual void setToZero( LevelData<FArrayBox>& a_x);
267  /*@}*/
268 
269  /**
270  \name MGLevelOp functions */
271  /*@{*/
272 
273  virtual void relax(LevelData<FArrayBox>& a_e,
274  const LevelData<FArrayBox>& a_residual,
275  int iterations);
276 
277  virtual void createCoarser(LevelData<FArrayBox>& a_coarse,
278  const LevelData<FArrayBox>& a_fine,
279  bool ghosted);
280  /**
281  calculate restricted residual
282  a_resCoarse[2h] = I[h->2h] (rhsFine[h] - L[h](phiFine[h])
283  */
284  virtual void restrictResidual(LevelData<FArrayBox>& a_resCoarse,
285  LevelData<FArrayBox>& a_phiFine,
286  const LevelData<FArrayBox>& a_rhsFine);
287 
288  /**
289  correct the fine solution based on coarse correction
290  a_phiThisLevel += I[2h->h](a_correctCoarse)
291  */
292  virtual void prolongIncrement(LevelData<FArrayBox>& a_phiThisLevel,
293  const LevelData<FArrayBox>& a_correctCoarse);
294 
295  /*@}*/
296 
297  /**
298  \name AMRLevelOp functions */
299  /*@{*/
300 
301  /** returns 1 when there are no coarser AMRLevelOp objects */
302  virtual int refToCoarser()
303  {
304  return m_refToCoar;
305  }
306 
307  /** a_residual = a_rhs - L(a_phi, a_phiFine, a_phiCoarse) */
308  virtual void AMRResidual(LevelData<FArrayBox>& a_residual,
309  const LevelData<FArrayBox>& a_phiFine,
310  const LevelData<FArrayBox>& a_phi,
311  const LevelData<FArrayBox>& a_phiCoarse,
312  const LevelData<FArrayBox>& a_rhs,
313  bool a_homogeneousPhysBC,
314  AMRLevelOp<LevelData<FArrayBox> >* a_finerOp);
315 
316  /** residual assuming no more coarser AMR levels */
317 
318  virtual void AMRResidualNC(LevelData<FArrayBox>& a_residual,
319  const LevelData<FArrayBox>& a_phiFine,
320  const LevelData<FArrayBox>& a_phi,
321  const LevelData<FArrayBox>& a_rhs,
322  bool a_homogeneousPhysBC,
323  AMRLevelOp<LevelData<FArrayBox> >* a_finerOp);
324 
325  /** a_residual = a_rhs - L(a_phi, a_phiCoarse) */
326  virtual void AMRResidualNF(LevelData<FArrayBox>& a_residual,
327  const LevelData<FArrayBox>& a_phi,
328  const LevelData<FArrayBox>& a_phiCoarse,
329  const LevelData<FArrayBox>& a_rhs,
330  bool a_homogeneousPhysBC);
331 
332  /** a_resCoarse = I[h-2h]( a_residual - L(a_correction, a_coarseCorrection))
333  it is assumed that a_resCoarse has already been filled in with the coarse
334  version of AMRResidualNF and that this operation is free to overwrite
335  in the overlap regions.
336  */
337 
338  virtual void AMRRestrict(LevelData<FArrayBox>& a_resCoarse,
339  const LevelData<FArrayBox>& a_residual,
340  const LevelData<FArrayBox>& a_correction,
341  const LevelData<FArrayBox>& a_coarseCorrection,
342  bool a_skip_res = false );
343 
344  /** a_correction += I[h->h](a_coarseCorrection) */
345  virtual void AMRProlong(LevelData<FArrayBox>& a_correction,
346  const LevelData<FArrayBox>& a_coarseCorrection);
347 
348  /** a_residual = a_residual - L(a_correction, a_coarseCorrection) */
349  virtual void AMRUpdateResidual(LevelData<FArrayBox>& a_residual,
350  const LevelData<FArrayBox>& a_correction,
351  const LevelData<FArrayBox>& a_coarseCorrection);
352 
353  ///
354  /**
355  compute norm over all cells on coarse not covered by finer
356  */
357  virtual Real AMRNorm(const LevelData<FArrayBox>& a_coarseResid,
358  const LevelData<FArrayBox>& a_fineResid,
359  const int& a_refRat,
360  const int& a_ord);
361 
362  virtual void outputLevel(LevelData<FArrayBox>& a_rhs,
363  string& a_name)
364  {
365  string fname(a_name);
366  fname.append(".hdf5");
367 #ifdef CH_HDF5
368  writeLevelname(&a_rhs, fname.c_str());
369 #endif
370  }
371 
372  //debugging hook
373  virtual void outputAMR(Vector<LevelData<FArrayBox>* >& a_rhs, string& a_name)
374  {
375  // for now, klugs is to use outputLevel on level 0
376  outputLevel(*a_rhs[0], a_name);
377  }
378 
380  ///
381  /**
382  does homogeneous coarse/fine interpolation for phi
383  */
384  void
386  const DataIndex& a_datInd,
387  int a_idir,
388  Side::LoHiSide a_hiorlo);
389 
390  ///
391  /** does homogeneous coarse/fine interpolation for tangential gradient
392  (needs phi ghost cells to be filled in first, so should call
393  homogeneousCFInterpPhi first)
394  */
396  const LevelData<FArrayBox>& a_phi,
397  const DataIndex& a_datInd,
398  int a_idir,
399  Side::LoHiSide a_hiorlo);
400 
402  const DataIndex& a_datInd,
403  const int a_idir,
404  const Side::LoHiSide a_hiorlo,
405  const IntVectSet& a_interpIVS);
406  ///
407  /**
408  Apply the AMR operator, including coarse-fine matching
409  */
410  void AMROperator( LevelData<FArrayBox>& a_LofPhi,
411  const LevelData<FArrayBox>& a_phiFine,
412  const LevelData<FArrayBox>& a_phi,
413  const LevelData<FArrayBox>& a_phiCoarse,
414  bool a_homogeneousDomBC,
415  AMRLevelOp<LevelData<FArrayBox> >* a_finerOp);
416 
417  void AMROperatorNF( LevelData<FArrayBox>& a_LofPhi,
418  const LevelData<FArrayBox>& a_phi,
419  const LevelData<FArrayBox>& a_phiCoarse,
420  bool a_homogeneousBC);
421 
422  virtual void AMROperatorNC( LevelData<FArrayBox>& a_LofPhi,
423  const LevelData<FArrayBox>& a_phiFine,
424  const LevelData<FArrayBox>& a_phi,
425  bool a_homogeneousBC,
426  AMRLevelOp<LevelData<FArrayBox> >* a_finerOp);
427 
428  void cellGrad(FArrayBox& a_gradPhi,
429  const FArrayBox& a_phi,
430  const Box& a_grid);
431 
432  void cfinterp(const LevelData<FArrayBox>& a_phiFine,
433  const LevelData<FArrayBox>& a_phiCoarse);
434 
435  void fillGrad(const LevelData<FArrayBox>& a_phiFine);
436 
437  /// optimized version of fillGrad which doesn't call exchange on phi
438  void fillGradNoExchange(const LevelData<FArrayBox>& a_phiFine);
439 
440  static void loHiCenterFace(Box& a_loBox,
441  int& a_hasLo,
442  Box& a_hiBox,
443  int& a_hasHi,
444  Box& a_centerBox,
445  const ProblemDomain& a_eblg,
446  const Box& a_inBox,
447  const int& a_dir);
448 
449  static void getFaceDivAndGrad(FArrayBox& a_faceDiv,
450  FArrayBox& a_faceGrad,
451  const FArrayBox& a_data,
452  const FArrayBox& a_gradData,
453  const ProblemDomain& a_domain,
454  const Box& a_faceBox,
455  const int& a_faceDir,
456  const Real a_dx);
457 
458  static void getFluxFromDivAndGrad(FArrayBox& a_flux,
459  const FArrayBox& a_faceDiv,
460  const FArrayBox& a_faceGrad,
461  const FArrayBox& a_etaFace,
462  const FArrayBox& a_lambdaFace,
463  const Box& a_faceBox,
464  int a_dir);
465 
466  ///take cell centered divergence of the inputs.
467  /**
468  not part of the operator evaluation. this is to
469  test whether the divergence of b converges at h^2
470  as b does if b is analytically divergence-free.
471  */
473  const LevelData<FArrayBox>& a_phi,
474  const LevelData<FArrayBox>* a_phiC);
475 
476  /// prolongation type
477  /** make this public/static to make it easy to change, while
478  also ensuring that it remains constent across all instances
479  */
480  static int s_prolongType;
481 
482 
483  /// if true, only do exchange once per GSRB pass, rather than for each color
484  /** make this public/static to make it easy to change, while
485  also ensuring that it remains constent across all instances
486  */
487  static bool s_lazy_gsrb;
488 
489  /// access function
490 
494  Real getAlpha() const {return m_alpha;}
495  Real getBeta() const {return m_beta;}
496 
497 
498 protected:
499  void defineRelCoef();
508 
509  //b is a scalar in 2d, vector in 3d
510  int m_ncomp;
514  //relaxation coef
515 public: // I want to get this from an owner of this -mfa
517 protected:
518  // underrelaxation parameter
520  // stopping tolerance for relaxation
522  // minimum number of smooths before tolerance check
524 
529 
532  // will need these to do tangential gradient computations
536 
537 private:
538  ///weak construction is bad
540  {
541  MayDay::Error("invalid operator");
542  }
543  //copy constructor and operator= disallowed for all the usual reasons
545  {
546  MayDay::Error("invalid operator");
547  }
548  void operator=(const ViscousTensorOp& a_opin)
549  {
550  MayDay::Error("invalid operator");
551  }
552 
553 };
554 
555 ///
556 /**
557  Factory to create ViscousTensorOps
558 */
559 class ViscousTensorOpFactory: public AMRLevelOpFactory<LevelData<FArrayBox> >
560 {
561 public:
563  {
564 #if 0
565  for (int ivec = 0; ivec < m_eta.size(); ivec++)
566  {
570  }
571 #endif
572 
573  }
574 
576  const Vector<RefCountedPtr<LevelData<FluxBox> > >& a_eta,
577  const Vector<RefCountedPtr<LevelData<FluxBox> > >& a_lambda,
578  const Vector<RefCountedPtr<LevelData<FArrayBox> > >& a_acoef,
579  Real a_alpha,
580  Real a_beta,
581  const Vector<int>& a_refRatios,
582  const ProblemDomain& a_domainCoar,
583  const Real& a_dxCoar,
584  BCFunc a_bc,
585  Real a_safety = VTOP_DEFAULT_SAFETY,
586  Real a_relaxTolerance = 0.0,
587  int a_relaxMinIter = 1000
588  );
589 
590  ViscousTensorOpFactory(const Vector<DisjointBoxLayout>& a_grids,
591  const Vector<RefCountedPtr<LevelData<FluxBox> > >& a_eta,
592  const Vector<RefCountedPtr<LevelData<FluxBox> > >& a_lambda,
593  const Vector<RefCountedPtr<LevelData<FArrayBox> > >& a_acoef,
594  Real a_alpha,
595  Real a_beta,
596  const Vector<int>& a_refRatios,
597  const ProblemDomain& a_domainCoar,
598  const Real& a_dxCoar,
599  BCHolder a_bc,
600  Real a_safety = VTOP_DEFAULT_SAFETY,
601  Real a_relaxTolerance = 0.0,
602  int a_relaxMinIter = 1000
603  );
604 
605  virtual void define(const Vector<DisjointBoxLayout>& a_grids,
606  const Vector<RefCountedPtr<LevelData<FluxBox> > >& a_eta,
607  const Vector<RefCountedPtr<LevelData<FluxBox> > >& a_lambda,
608  const Vector<RefCountedPtr<LevelData<FArrayBox> > >& a_acoef,
609  Real a_alpha,
610  Real a_beta,
611  const Vector<int>& a_refRatios,
612  const ProblemDomain& a_domainCoar,
613  const Real& a_dxCoar,
614  BCHolder a_bc,
615  Real a_safety = VTOP_DEFAULT_SAFETY,
616  Real a_relaxTolerance = 0.0,
617  int a_relaxMinIter = 1000
618  );
619 
620  ///
621  virtual ViscousTensorOp*
622  MGnewOp(const ProblemDomain& a_FineindexSpace,
623  int depth,
624  bool homoOnly = true);
625 
626  ///
627  virtual ViscousTensorOp* AMRnewOp(const ProblemDomain& a_indexSpace);
628 
629  ///
630  virtual int refToFiner(const ProblemDomain&) const;
631 
632  /// set any relevant default values
633  void setDefaults();
634 
635  // 0 = arithmetic, 1 = harmonic
637 
638 
639 private:
647  // refinement to next coarser level
655 
656  ///weak construction is bad
658  {
659  MayDay::Error("invalid operator");
660  }
661  //copy constructor and operator= disallowed for all the usual reasons
662  ViscousTensorOpFactory(const ViscousTensorOpFactory& a_opin)
663  {
664  MayDay::Error("invalid operator");
665  }
666 
667  void operator=(const ViscousTensorOpFactory& a_opin)
668  {
669  MayDay::Error("invalid operator");
670  }
671 };
672 
673 extern void
674 coarsenStuff(LevelData<FluxBox> & a_etaCoar,
675  LevelData<FluxBox> & a_lambdaCoar,
676  const LevelData<FluxBox> & a_etaFine,
677  const LevelData<FluxBox> & a_lambdaFine,
678  const int & a_refToDepth);
679 
680 #include "NamespaceFooter.H"
681 #endif
int m_ncomp
Definition: ViscousTensorOp.H:510
virtual void create(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
RefCountedPtr< LevelData< FArrayBox > > getACoef() const
Definition: ViscousTensorOp.H:493
ProblemDomain m_domain
Definition: ViscousTensorOp.H:513
Real m_alpha
Definition: ViscousTensorOp.H:503
void reflux(const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, LevelData< FArrayBox > &residual, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
void defineRelCoef()
virtual void residual(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs, bool a_homogeneous=false)
Real m_dxCrse
Definition: ViscousTensorOp.H:512
Definition: ViscousTensorOp.H:42
void computeOperatorNoBCs(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi)
utility function which computes operator after all bc&#39;s have been set
A reference-counting handle class.
Definition: RefCountedPtr.H:173
LevelData< FArrayBox > m_relaxCoef
Definition: ViscousTensorOp.H:516
void fillGradNoExchange(const LevelData< FArrayBox > &a_phiFine)
optimized version of fillGrad which doesn&#39;t call exchange on phi
An irregular domain on an integer lattice.
Definition: IntVectSet.H:44
virtual ~ViscousTensorOp()
Definition: ViscousTensorOp.H:88
virtual Real dxCrse() const
Definition: ViscousTensorOp.H:261
Definition: BCFunc.H:136
ViscousTensorOp()
weak construction is bad
Definition: ViscousTensorOp.H:539
Real m_beta
Definition: ViscousTensorOp.H:504
A class to facilitate interaction with physical boundary conditions.
Definition: ProblemDomain.H:141
void homogeneousCFInterpTanGrad(LevelData< FArrayBox > &a_tanGrad, const LevelData< FArrayBox > &a_phi, const DataIndex &a_datInd, int a_idir, Side::LoHiSide a_hiorlo)
LayoutData< TensorFineStencilSet > m_hiTanStencilSets[SpaceDim]
Definition: ViscousTensorOp.H:533
virtual void applyOp(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, bool a_homogeneous=false)
virtual void incr(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_x, Real a_scale)
void AMROperatorNF(LevelData< FArrayBox > &a_LofPhi, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, bool a_homogeneousBC)
virtual void relax(LevelData< FArrayBox > &a_e, const LevelData< FArrayBox > &a_residual, int iterations)
LevelData< FArrayBox > m_grad
Definition: ViscousTensorOp.H:525
one dimensional dynamic array
Definition: Vector.H:53
void coarsenStuff(LevelData< FluxBox > &a_etaCoar, LevelData< FluxBox > &a_lambdaCoar, const LevelData< FluxBox > &a_etaFine, const LevelData< FluxBox > &a_lambdaFine, const int &a_refToDepth)
BCHolder m_bc
Definition: ViscousTensorOp.H:507
Vector< DisjointBoxLayout > m_boxes
Definition: ViscousTensorOp.H:645
A strange but true thing to make copying from one boxlayoutdata to another fast.
Definition: Copier.H:152
void define(const DisjointBoxLayout &a_grids, const DisjointBoxLayout &a_gridsFine, const DisjointBoxLayout &a_gridsCoar, const RefCountedPtr< LevelData< FluxBox > > &a_eta, const RefCountedPtr< LevelData< FluxBox > > &a_lambda, const RefCountedPtr< LevelData< FArrayBox > > &a_acoef, Real a_alpha, Real a_beta, int a_refToFine, int a_refToCoar, const ProblemDomain &a_domain, const Real &a_dxLevel, const Real &a_dxCoar, BCHolder &a_bc, Real a_safety=VTOP_DEFAULT_SAFETY, Real a_relaxTolerance=0.0, int a_relaxMinIter=1000)
Definition: ViscousTensorOp.H:35
static int s_coefficientAverageType
Definition: ViscousTensorOp.H:636
RefCountedPtr< LevelData< FluxBox > > m_lambda
Definition: ViscousTensorOp.H:501
static void getFluxFromDivAndGrad(FArrayBox &a_flux, const FArrayBox &a_faceDiv, const FArrayBox &a_faceGrad, const FArrayBox &a_etaFace, const FArrayBox &a_lambdaFace, const Box &a_faceBox, int a_dir)
Real getBeta() const
Definition: ViscousTensorOp.H:495
Quadratic coarse-fine interpolation utility for tensors.
Definition: TensorCFInterp.H:36
virtual bool ok() const
return true if this iterator is still in its Layout
Definition: LayoutIterator.H:117
virtual void axby(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_x, const LevelData< FArrayBox > &a_y, Real a, Real b)
Definition: DataIterator.H:190
virtual void getFlux(FluxBox &a_flux, const LevelData< FArrayBox > &a_data, const Box &a_grid, const DataIndex &a_dit, Real a_scale)
Definition: ViscousTensorOp.H:200
void fillGrad(const LevelData< FArrayBox > &a_phiFine)
These functions are part of the LevelTGA interface......
LayoutData< CFIVS > m_loCFIVS[SpaceDim]
Definition: ViscousTensorOp.H:530
LevelDataOps< FArrayBox > m_levelOps
Definition: ViscousTensorOp.H:526
virtual void diagonalScale(LevelData< FArrayBox > &a_rhs, bool a_kappaWeighted)
Definition: ViscousTensorOp.H:46
Real m_beta
Definition: ViscousTensorOp.H:651
BCHolder m_bc
Definition: ViscousTensorOp.H:649
void(* BCFunc)(FArrayBox &a_state, const Box &a_valid, const ProblemDomain &a_domain, Real a_dx, bool a_homogeneous)
Definition: BCFunc.H:30
void divergenceCC(LevelData< FArrayBox > &a_div, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > *a_phiC)
take cell centered divergence of the inputs.
virtual Real dx() const
Definition: ViscousTensorOp.H:256
const int SpaceDim
Definition: SPACE.H:38
Definition: AMRMultiGrid.H:39
virtual void outputLevel(LevelData< FArrayBox > &a_rhs, string &a_name)
Definition: ViscousTensorOp.H:362
virtual void setToZero(LevelData< FArrayBox > &a_x)
virtual void createCoarser(LevelData< FArrayBox > &a_coarse, const LevelData< FArrayBox > &a_fine, bool ghosted)
void operator=(const ViscousTensorOp &a_opin)
Definition: ViscousTensorOp.H:548
Definition: ViscousTensorOp.H:41
static void loHiCenterFace(Box &a_loBox, int &a_hasLo, Box &a_hiBox, int &a_hasHi, Box &a_centerBox, const ProblemDomain &a_eblg, const Box &a_inBox, const int &a_dir)
Vector< RefCountedPtr< LevelData< FArrayBox > > > m_acoef
Definition: ViscousTensorOp.H:642
A FArrayBox-like container for face-centered fluxes.
Definition: FluxBox.H:22
void cellGrad(FArrayBox &a_gradPhi, const FArrayBox &a_phi, const Box &a_grid)
void applyOpNoBoundary(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi)
Definition: ViscousTensorOp.H:220
virtual void AMRResidualNF(LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, const LevelData< FArrayBox > &a_rhs, bool a_homogeneousPhysBC)
virtual void restrictResidual(LevelData< FArrayBox > &a_resCoarse, LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_rhsFine)
Vector< ProblemDomain > m_domains
Definition: ViscousTensorOp.H:644
#define VTOP_DEFAULT_SAFETY
Definition: ViscousTensorOp.H:29
virtual void diagonalScale(LevelData< FArrayBox > &a_rhs)
Definition: ViscousTensorOp.H:59
Real m_alpha
Definition: ViscousTensorOp.H:650
prolongationType
Definition: ViscousTensorOp.H:39
virtual int refToCoarser()
Definition: ViscousTensorOp.H:302
virtual void divideByIdentityCoef(LevelData< FArrayBox > &a_rhs)
Definition: ViscousTensorOp.H:72
virtual void outputAMR(Vector< LevelData< FArrayBox > * > &a_rhs, string &a_name)
Definition: ViscousTensorOp.H:373
int m_refToCoar
Definition: ViscousTensorOp.H:505
Real m_dx
Definition: ViscousTensorOp.H:511
RefCountedPtr< LevelData< FluxBox > > getLambda() const
Definition: ViscousTensorOp.H:492
double Real
Definition: REAL.H:33
Box surroundingNodes(const Box &b, int dir)
Definition: Box.H:2161
LayoutData< TensorFineStencilSet > m_loTanStencilSets[SpaceDim]
Definition: ViscousTensorOp.H:534
Definition: ViscousTensorOp.H:559
Real m_relaxTolerance
Definition: ViscousTensorOp.H:653
ViscousTensorOp(const ViscousTensorOp &a_opin)
Definition: ViscousTensorOp.H:544
virtual void setAlphaAndBeta(const Real &a_alpha, const Real &a_beta)
Definition: ViscousTensorOp.H:51
static bool s_lazy_gsrb
if true, only do exchange once per GSRB pass, rather than for each color
Definition: ViscousTensorOp.H:487
A BoxLayout that has a concept of disjointedness.
Definition: DisjointBoxLayout.H:30
virtual void assign(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
LoHiSide
Definition: LoHiSide.H:27
void interpOnIVSHomo(LevelData< FArrayBox > &a_phif, const DataIndex &a_datInd, const int a_idir, const Side::LoHiSide a_hiorlo, const IntVectSet &a_interpIVS)
Copier m_exchangeCopier
Definition: ViscousTensorOp.H:527
virtual void residualNoExchange(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs, bool a_homogeneous=false)
void homogeneousCFInterpPhi(LevelData< FArrayBox > &a_phif, const DataIndex &a_datInd, int a_idir, Side::LoHiSide a_hiorlo)
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.
void homogeneousCFInterp(LevelData< FArrayBox > &a_phif)
virtual void AMRResidualNC(LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs, bool a_homogeneousPhysBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
void cfinterp(const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phiCoarse)
virtual void AMRRestrict(LevelData< FArrayBox > &a_resCoarse, const LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection, bool a_skip_res=false)
virtual void applyOpNoExchange(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, bool a_homogeneous=false)
RefCountedPtr< LevelData< FluxBox > > getEta() const
access function
Definition: ViscousTensorOp.H:491
LayoutData< CFIVS > m_hiCFIVS[SpaceDim]
Definition: ViscousTensorOp.H:531
A Rectangular Domain on an Integer Lattice.
Definition: Box.H:469
int m_relaxMinIter
Definition: ViscousTensorOp.H:654
Vector< Real > m_dx
Definition: ViscousTensorOp.H:646
int m_relaxMinIter
Definition: ViscousTensorOp.H:523
void operator=(const ViscousTensorOpFactory &a_opin)
Definition: ViscousTensorOp.H:667
Definition: DataIndex.H:114
Vector< RefCountedPtr< LevelData< FluxBox > > > m_eta
Definition: ViscousTensorOp.H:640
const DisjointBoxLayout & disjointBoxLayout() const
Definition: LevelData.H:225
virtual void createCoarsened(LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs, const int &a_refRat)
virtual void AMRProlong(LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection)
const Box & box() const
Returns cell-centered box which defines fluxBox.
Definition: FArrayBox.H:45
Definition: AMRTGA.H:159
virtual void preCond(LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_residual)
ViscousTensorOpFactory()
weak construction is bad
Definition: ViscousTensorOp.H:657
void AMROperator(LevelData< FArrayBox > &a_LofPhi, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, bool a_homogeneousDomBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
virtual void AMRResidual(LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, const LevelData< FArrayBox > &a_rhs, bool a_homogeneousPhysBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
Real m_safety
Definition: ViscousTensorOp.H:519
RefCountedPtr< LevelData< FluxBox > > m_eta
Definition: ViscousTensorOp.H:500
virtual Real norm(const LevelData< FArrayBox > &a_x, int a_ord)
Real m_relaxTolerance
Definition: ViscousTensorOp.H:521
Vector< int > m_refRatios
Definition: ViscousTensorOp.H:648
Real m_safety
Definition: ViscousTensorOp.H:652
virtual void prolongIncrement(LevelData< FArrayBox > &a_phiThisLevel, const LevelData< FArrayBox > &a_correctCoarse)
void writeLevelname(const LevelData< FArrayBox > *a_dataPtr, const char *a_filename)
void define(const Box &bx, int n=1)
Resize FluxBox similar to BaseFab::resize()
static int s_prolongType
prolongation type
Definition: ViscousTensorOp.H:480
Real getAlpha() const
Definition: ViscousTensorOp.H:494
const Box & domainBox() const
Returns the logical computational domain.
Definition: ProblemDomain.H:887
virtual Real AMRNorm(const LevelData< FArrayBox > &a_coarseResid, const LevelData< FArrayBox > &a_fineResid, const int &a_refRat, const int &a_ord)
Definition: AMRMultiGrid.H:233
DataIterator dataIterator() const
Parallel iterator.
virtual void scale(LevelData< FArrayBox > &a_lhs, const Real &a_scale)
virtual ~ViscousTensorOpFactory()
Definition: ViscousTensorOp.H:562
virtual Real dotProduct(const LevelData< FArrayBox > &a_1, const LevelData< FArrayBox > &a_2)
Vector< IntVect > m_colors
Definition: ViscousTensorOp.H:535
LevelData< FArrayBox > m_phic
Definition: ViscousTensorOp.H:643
virtual void AMRUpdateResidual(LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection)
virtual void AMROperatorNC(LevelData< FArrayBox > &a_LofPhi, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, bool a_homogeneousBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
TensorCFInterp m_interpWithCoarser
Definition: ViscousTensorOp.H:528
static void getFaceDivAndGrad(FArrayBox &a_faceDiv, FArrayBox &a_faceGrad, const FArrayBox &a_data, const FArrayBox &a_gradData, const ProblemDomain &a_domain, const Box &a_faceBox, const int &a_faceDir, const Real a_dx)
int m_refToFine
Definition: ViscousTensorOp.H:506
ViscousTensorOpFactory(const ViscousTensorOpFactory &a_opin)
Definition: ViscousTensorOp.H:662
Vector< RefCountedPtr< LevelData< FluxBox > > > m_lambda
Definition: ViscousTensorOp.H:641
Definition: ViscousTensorOp.H:43
RefCountedPtr< LevelData< FArrayBox > > m_acoef
Definition: ViscousTensorOp.H:502