Chombo + EB: ResistivityOp.H Source File

00001 #ifdef CH_LANG_CC
00002 /*
00003 *      _______              __
00004 *     / ___/ /  ___  __ _  / /  ___
00005 *    / /__/ _ \/ _ \/  V \/ _ \/ _ \
00006 *    \___/_//_/\___/_/_/_/_.__/\___/
00007 *    Please refer to Copyright.txt, in Chombo's root directory.
00008 */
00009 #endif
00010 
00011 // dtg 2007
00012 
00013 #ifndef _RESISTIVITYOP_H_
00014 #define _RESISTIVITYOP_H_
00015 
00016 #include "AMRMultiGrid.H"
00017 #include "REAL.H"
00018 #include "Box.H"
00019 #include "LevelDataOps.H"
00020 #include "BCFunc.H"
00021 #include "FArrayBox.H"
00022 #include "FluxBox.H"
00023 #include "CFIVS.H"
00024 #include "TensorCFInterp.H"
00025 #include "CoarseAverage.H"
00026 #include "AMRTGA.H"
00027 #include "NamespaceHeader.H"
00028 
00030 
00040 class ResistivityOp : public LevelTGAHelmOp<LevelData<FArrayBox>, FluxBox>
00041 {
00042 public:
00043 
00044   virtual void setAlphaAndBeta(const Real& a_alpha,
00045                                const Real& a_beta)
00046   {
00047     m_alpha = a_alpha;
00048     m_beta = a_beta;
00049     setLambda();
00050   }
00051 
00052   virtual void diagonalScale(LevelData<FArrayBox>& a_rhs)
00053   {
00054     //does not apply here
00055   }
00056 
00058 
00060   virtual ~ResistivityOp() { ; }
00061 
00062   //use at thy peril
00063   virtual void getFlux(FluxBox&             a_flux,
00064                const LevelData<FArrayBox>&  a_data,
00065                const Box&       a_grid,
00066                const DataIndex& a_dit,
00067                Real             a_scale);
00068 
00070 
00072   ResistivityOp(const DisjointBoxLayout&                    a_grids,
00073                 const DisjointBoxLayout&                    a_gridsFine,
00074                 const DisjointBoxLayout&                    a_gridsCoar,
00075                 const RefCountedPtr<LevelData<FluxBox> >&   a_eta,
00076                 Real                                        a_alpha,
00077                 Real                                        a_beta,
00078                 int                                         a_refToFine,
00079                 int                                         a_refToCoar,
00080                 const ProblemDomain&                        a_domain,
00081                 const Real&                                 a_dxLevel,
00082                 const Real&                                 a_dxCoar,
00083                 BCFunc                                      a_bc);
00084 
00085   virtual void residual(  LevelData<FArrayBox>& a_lhs,
00086                           const LevelData<FArrayBox>& a_phi,
00087                           const LevelData<FArrayBox>& a_rhs,
00088                           bool a_homogeneous = false);
00089 
00090   virtual void preCond(   LevelData<FArrayBox>& a_correction,
00091                           const LevelData<FArrayBox>& a_residual);
00092 
00093   virtual void applyOp(   LevelData<FArrayBox>& a_lhs,
00094                           const LevelData<FArrayBox>& a_phi,
00095                           bool a_homogeneous = false);
00096 
00097   virtual void applyOpNoBoundary( LevelData<FArrayBox>& a_lhs,
00098                                   const LevelData<FArrayBox>& a_phi);
00099 
00100   virtual void create(    LevelData<FArrayBox>& a_lhs,
00101                           const LevelData<FArrayBox>& a_rhs);
00102   virtual void createCoarsened(    LevelData<FArrayBox>& a_lhs,
00103                                    const LevelData<FArrayBox>& a_rhs,
00104                                    const int& a_refRat);
00105 
00106   void reflux(const LevelData<FArrayBox>& a_phiFine,
00107               const LevelData<FArrayBox>& a_phi,
00108               LevelData<FArrayBox>& residual,
00109               AMRLevelOp<LevelData<FArrayBox> >* a_finerOp);
00110 
00111   void getFlux(FArrayBox&       a_flux,
00112                const FArrayBox& a_data,
00113                const FArrayBox& a_gradData,
00114                const FArrayBox& a_etaFace,
00115                const Box&       a_facebox,
00116                int              a_dir,
00117                int ref = 1);
00118 
00120   void computeOperatorNoBCs(LevelData<FArrayBox>& a_lhs,
00121                             const LevelData<FArrayBox>& a_phi);
00122 
00123   virtual void assign(    LevelData<FArrayBox>& a_lhs,
00124                           const LevelData<FArrayBox>& a_rhs) ;
00125   virtual Real dotProduct(const LevelData<FArrayBox>& a_1,
00126                           const LevelData<FArrayBox>& a_2) ;
00127   virtual void incr( LevelData<FArrayBox>& a_lhs,
00128                      const LevelData<FArrayBox>& a_x,
00129                      Real a_scale) ;
00130   virtual void axby( LevelData<FArrayBox>& a_lhs,
00131                      const LevelData<FArrayBox>& a_x,
00132                      const LevelData<FArrayBox>& a_y,
00133                      Real a, Real b) ;
00134 
00135   virtual void scale(LevelData<FArrayBox>& a_lhs, const Real& a_scale) ;
00136 
00137   virtual Real norm(const LevelData<FArrayBox>& a_x, int a_ord);
00138 
00139   virtual void setToZero( LevelData<FArrayBox>& a_x);
00145 
00146   virtual void relax(LevelData<FArrayBox>& a_e,
00147                      const LevelData<FArrayBox>& a_residual,
00148                      int iterations);
00149 
00150   virtual void createCoarser(LevelData<FArrayBox>& a_coarse,
00151                              const LevelData<FArrayBox>& a_fine,
00152                              bool ghosted);
00157   virtual void restrictResidual(LevelData<FArrayBox>& a_resCoarse,
00158                                 LevelData<FArrayBox>& a_phiFine,
00159                                 const LevelData<FArrayBox>& a_rhsFine);
00160 
00165   virtual void prolongIncrement(LevelData<FArrayBox>& a_phiThisLevel,
00166                                 const LevelData<FArrayBox>& a_correctCoarse);
00167 
00173 
00175   virtual int refToCoarser(){ return m_refToCoar; }
00176 
00178   virtual void AMRResidual(LevelData<FArrayBox>& a_residual,
00179                            const LevelData<FArrayBox>& a_phiFine,
00180                            const LevelData<FArrayBox>& a_phi,
00181                            const LevelData<FArrayBox>& a_phiCoarse,
00182                            const LevelData<FArrayBox>& a_rhs,
00183                            bool a_homogeneousPhysBC,
00184                            AMRLevelOp<LevelData<FArrayBox> >* a_finerOp);
00185 
00188   virtual void AMRResidualNC(LevelData<FArrayBox>& a_residual,
00189                              const LevelData<FArrayBox>& a_phiFine,
00190                              const LevelData<FArrayBox>& a_phi,
00191                              const LevelData<FArrayBox>& a_rhs,
00192                              bool a_homogeneousPhysBC,
00193                              AMRLevelOp<LevelData<FArrayBox> >* a_finerOp);
00194 
00196   virtual void AMRResidualNF(LevelData<FArrayBox>& a_residual,
00197                              const LevelData<FArrayBox>& a_phi,
00198                              const LevelData<FArrayBox>& a_phiCoarse,
00199                              const LevelData<FArrayBox>& a_rhs,
00200                              bool a_homogeneousPhysBC);
00201 
00208   virtual void AMRRestrict(LevelData<FArrayBox>& a_resCoarse,
00209                            const LevelData<FArrayBox>& a_residual,
00210                            const LevelData<FArrayBox>& a_correction,
00211                            const LevelData<FArrayBox>& a_coarseCorrection);
00212 
00214   virtual void AMRProlong(LevelData<FArrayBox>& a_correction,
00215                           const LevelData<FArrayBox>& a_coarseCorrection);
00216 
00218   virtual void AMRUpdateResidual(LevelData<FArrayBox>& a_residual,
00219                                  const LevelData<FArrayBox>& a_correction,
00220                                  const LevelData<FArrayBox>& a_coarseCorrection);
00221 
00223 
00226   virtual Real AMRNorm(const LevelData<FArrayBox>& a_coarseResid,
00227                        const LevelData<FArrayBox>& a_fineResid,
00228                        const int&                  a_refRat,
00229                        const int&                  a_ord);
00230 
00231   void homogeneousCFInterp(LevelData<FArrayBox>& a_phif);
00233 
00236   void
00237   homogeneousCFInterpPhi(LevelData<FArrayBox>& a_phif,
00238                          const DataIndex& a_datInd,
00239                          int a_idir,
00240                          Side::LoHiSide a_hiorlo);
00241 
00243 
00247   void homogeneousCFInterpTanGrad(LevelData<FArrayBox>& a_tanGrad,
00248                                   const LevelData<FArrayBox>& a_phi,
00249                                   const DataIndex& a_datInd,
00250                                   int a_idir,
00251                                   Side::LoHiSide a_hiorlo);
00252 
00253   void interpOnIVSHomo(LevelData<FArrayBox>& a_phif,
00254                        const DataIndex& a_datInd,
00255                        const int a_idir,
00256                        const Side::LoHiSide a_hiorlo,
00257                        const IntVectSet& a_interpIVS);
00259 
00262   void AMROperator(      LevelData<FArrayBox>& a_LofPhi,
00263                    const LevelData<FArrayBox>& a_phiFine,
00264                    const LevelData<FArrayBox>& a_phi,
00265                    const LevelData<FArrayBox>& a_phiCoarse,
00266                    bool a_homogeneousDomBC,
00267                    AMRLevelOp<LevelData<FArrayBox> >*  a_finerOp);
00268 
00269   void AMROperatorNF(      LevelData<FArrayBox>& a_LofPhi,
00270                      const LevelData<FArrayBox>& a_phi,
00271                      const LevelData<FArrayBox>& a_phiCoarse,
00272                      bool a_homogeneousBC);
00273 
00274   virtual void AMROperatorNC(      LevelData<FArrayBox>& a_LofPhi,
00275                              const LevelData<FArrayBox>& a_phiFine,
00276                              const LevelData<FArrayBox>& a_phi,
00277                              bool a_homogeneousBC,
00278                              AMRLevelOp<LevelData<FArrayBox> >* a_finerOp);
00279 
00280   void cellGrad(FArrayBox&             a_gradPhi,
00281                 const  FArrayBox&      a_phi,
00282                 const Box&             a_grid);
00283 
00284   void cfinterp(const LevelData<FArrayBox>& a_phiFine,
00285                 const LevelData<FArrayBox>& a_phiCoarse);
00286 
00287   virtual void fillGrad(const LevelData<FArrayBox>& a_phiFine);
00288 
00289   void loHiCenterFace(Box&                 a_loBox,
00290                       int&                 a_hasLo,
00291                       Box&                 a_hiBox,
00292                       int&                 a_hasHi,
00293                       Box&                 a_centerBox,
00294                       const ProblemDomain& a_eblg,
00295                       const Box&           a_inBox,
00296                       const int&           a_dir);
00297 
00298   void getFaceDivAndGrad(FArrayBox&             a_faceDiv,
00299                          FArrayBox&             a_faceGrad,
00300                          const FArrayBox&       a_data,
00301                          const FArrayBox&       a_gradData,
00302                          const ProblemDomain&   a_domain,
00303                          const Box&             a_faceBox,
00304                          const int&             a_faceDir,
00305                          const Real             a_dx);
00306 
00308 
00313   void divergenceCC(LevelData<FArrayBox>&       a_div,
00314                     const LevelData<FArrayBox>& a_phi,
00315                     const LevelData<FArrayBox>* a_phiC);
00316 protected:
00317   void setLambda();
00318   RefCountedPtr<LevelData<FluxBox> >         m_eta;
00319   DisjointBoxLayout                          m_grids;
00320   Real                                       m_alpha;
00321   Real                                       m_beta;
00322   int                                        m_refToCoar;
00323   int                                        m_refToFine;
00324   BCFunc                                     m_bc;
00325 
00326   //b is a scalar in 2d, vector in 3d
00327   static const int        s_nComp; //==3
00328   static const int        s_nGradComp; //==3*SpaceDim
00329   Real                    m_dx;
00330   Real                    m_dxCrse;
00331   ProblemDomain           m_domain;
00332   //relaxation coef
00333   LevelData<FArrayBox>    m_lambda;
00334   LevelData<FArrayBox>    m_grad;
00335   LevelDataOps<FArrayBox> m_levelOps;
00336   Copier                  m_exchangeCopier;
00337   TensorCFInterp          m_interpWithCoarser;
00338 
00339   LayoutData<CFIVS> m_loCFIVS[SpaceDim];
00340   LayoutData<CFIVS> m_hiCFIVS[SpaceDim];
00341   // will need these to do tangential gradient computations
00342   LayoutData<TensorFineStencilSet> m_hiTanStencilSets[SpaceDim];
00343   LayoutData<TensorFineStencilSet> m_loTanStencilSets[SpaceDim];
00344   Vector<IntVect> m_colors;
00345 private:
00347   ResistivityOp()
00348   {
00349     MayDay::Error("invalid operator");
00350   }
00351   //copy constructor and operator= disallowed for all the usual reasons
00352   ResistivityOp(const ResistivityOp& a_opin)
00353   {
00354     MayDay::Error("invalid operator");
00355   }
00356   void operator=(const ResistivityOp& a_opin)
00357   {
00358     MayDay::Error("invalid operator");
00359   }
00360 };
00361 
00363 
00366 class ResistivityOpFactory: public AMRLevelOpFactory<LevelData< FArrayBox> >
00367 {
00368 public:
00369   virtual ~ResistivityOpFactory(){;}
00370 
00371   ResistivityOpFactory(const Vector<DisjointBoxLayout>&                     a_grids,
00372                        const Vector<RefCountedPtr<LevelData<FluxBox> > >&   a_eta,
00373                        Real                                                 a_alpha,
00374                        Real                                                 a_beta,
00375                        const Vector<int>&                                   a_refRatios,
00376                        const ProblemDomain&                                 a_domainCoar,
00377                        const Real&                                          a_dxCoar,
00378                        BCFunc                                               a_bc);
00379 
00381   virtual MGLevelOp< LevelData<FArrayBox> >*
00382   MGnewOp(const ProblemDomain& a_FineindexSpace,
00383           int depth,
00384           bool homoOnly = true);
00385 
00387   virtual AMRLevelOp< LevelData<FArrayBox> >*
00388   AMRnewOp(const ProblemDomain& a_indexSpace);
00389 
00391   virtual int refToFiner(const ProblemDomain&) const;
00392 
00393 private:
00394   Vector<RefCountedPtr<LevelData<FluxBox> > >    m_eta;
00395   LevelData<FArrayBox>                           m_phic;
00396   Vector<ProblemDomain>                          m_domains;
00397   Vector<DisjointBoxLayout>                      m_boxes;
00398   Vector<Real>                                   m_dx;
00399   // refinement to next coarser level
00400   Vector<int>                                    m_refRatios;
00401   BCFunc  m_bc;
00402   Real m_alpha, m_beta;
00403 
00405   ResistivityOpFactory()
00406   {
00407     MayDay::Error("invalid operator");
00408   }
00409   //copy constructor and operator= disallowed for all the usual reasons
00410   ResistivityOpFactory(const ResistivityOpFactory& a_opin)
00411   {
00412     MayDay::Error("invalid operator");
00413   }
00414 
00415   void operator=(const ResistivityOpFactory& a_opin)
00416   {
00417     MayDay::Error("invalid operator");
00418   }
00419 };
00420 
00421 #include "NamespaceFooter.H"
00422 #endif