Block-Structured AMR Software Framework
AMReX_MLEBABecLap.H
Go to the documentation of this file.
1 #ifndef AMREX_MLEBABECLAP_H_
2 #define AMREX_MLEBABECLAP_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_EBFabFactory.H>
6 #include <AMReX_MLCellABecLap.H>
7 #include <AMReX_Array.H>
8 #include <limits>
9 
10 namespace amrex {
11 
12 // (alpha * a - beta * (del dot b grad)) phi
13 
15  : public MLCellABecLap
16 {
17 public:
18 
19  MLEBABecLap () = default;
20  MLEBABecLap (const Vector<Geometry>& a_geom,
21  const Vector<BoxArray>& a_grids,
22  const Vector<DistributionMapping>& a_dmap,
23  const LPInfo& a_info,
24  const Vector<EBFArrayBoxFactory const*>& a_factory,
25  int a_ncomp = 1);
26 
27  ~MLEBABecLap () override;
28 
29  MLEBABecLap (const MLEBABecLap&) = delete;
30  MLEBABecLap (MLEBABecLap&&) = delete;
31  MLEBABecLap& operator= (const MLEBABecLap&) = delete;
33 
34  void define (const Vector<Geometry>& a_geom,
35  const Vector<BoxArray>& a_grids,
36  const Vector<DistributionMapping>& a_dmap,
37  const LPInfo& a_info,
38  const Vector<EBFArrayBoxFactory const*>& a_factory,
39  int a_ncomp = 1);
40 
41  void setPhiOnCentroid ();
42 
43  void setScalars (Real a, Real b);
44  void setACoeffs (int amrlev, const MultiFab& alpha);
45  void setACoeffs (int amrlev, Real alpha);
46 
47  void setBCoeffs (int amrlev, const Array<MultiFab const*,AMREX_SPACEDIM>& beta,
48  Location a_beta_loc);
49  void setBCoeffs (int amrlev, const Array<MultiFab const*,AMREX_SPACEDIM>& beta)
50  {setBCoeffs (amrlev, beta, Location::FaceCenter);}
51 
52  void setBCoeffs (int amrlev, Real beta);
53  void setBCoeffs (int amrlev, Vector<Real> const& beta);
54 
55  // Tells the solver that EB boundaries have Dirichlet bc's specified by "phi"
56  void setEBDirichlet (int amrlev, const MultiFab& phi, const MultiFab& beta);
57  void setEBDirichlet (int amrlev, const MultiFab& phi, Real beta);
58  void setEBDirichlet (int amrlev, const MultiFab& phi, Vector<Real> const& beta);
59 
60  // Tells the solver that EB boundaries have homogeneous Dirichlet bc's
61  void setEBHomogDirichlet (int amrlev, const MultiFab& beta);
62  void setEBHomogDirichlet (int amrlev, Real beta);
63  void setEBHomogDirichlet (int amrlev, Vector<Real> const& beta);
64 
65  int getNComp () const override { return m_ncomp; }
66 
67  bool needsUpdate () const override {
69  }
70  void update () override;
71 
72  std::unique_ptr<FabFactory<FArrayBox> > makeFactory (int amrlev, int mglev) const final;
73 
74  bool isCrossStencil () const override { return false; }
75 
76  void applyBC (int amrlev, int mglev, MultiFab& in, BCMode bc_mode, StateMode s_mode,
77  const MLMGBndry* bndry=nullptr, bool skip_fillboundary=false) const final;
78  void apply (int amrlev, int mglev, MultiFab& out, MultiFab& in, BCMode bc_mode,
79  StateMode s_mode, const MLMGBndry* bndry=nullptr) const override;
80  void compGrad (int amrlev, const Array<MultiFab*,AMREX_SPACEDIM>& grad,
81  MultiFab& sol, Location loc) const final;
82 
83  void prepareForSolve () override;
84  bool isSingular (int amrlev) const override { return m_is_singular[amrlev]; }
85  bool isBottomSingular () const override { return m_is_singular[0]; }
86  void Fapply (int amrlev, int mglev, MultiFab& out, const MultiFab& in) const final;
87  void Fsmooth (int amrlev, int mglev, MultiFab& sol, const MultiFab& rhs, int redblack) const final;
88  void FFlux (int amrlev, const MFIter& mfi,
90  const FArrayBox& sol, Location loc,
91  int face_only=0) const final;
92 
93  void normalize (int amrlev, int mglev, MultiFab& mf) const final;
94 
95  Real getAScalar () const final { return m_a_scalar; }
96  Real getBScalar () const final { return m_b_scalar; }
97  MultiFab const* getACoeffs (int amrlev, int mglev) const final
98  { return &(m_a_coeffs[amrlev][mglev]); }
99  Array<MultiFab const*,AMREX_SPACEDIM> getBCoeffs (int amrlev, int mglev) const final
100  { return amrex::GetArrOfConstPtrs(m_b_coeffs[amrlev][mglev]); }
101 
102  std::unique_ptr<MLLinOp> makeNLinOp (int /*grid_size*/) const final {
103  amrex::Abort("MLABecLaplacian::makeNLinOp: Not implemented");
104  return std::unique_ptr<MLLinOp>{};
105  }
106 
107  void restriction (int, int, MultiFab& crse, MultiFab& fine) const final;
108 
109  void interpolation (int amrlev, int fmglev, MultiFab& fine, const MultiFab& crse) const final;
110 
111  void averageDownSolutionRHS (int camrlev, MultiFab& crse_sol, MultiFab& crse_rhs,
112  const MultiFab& fine_sol, const MultiFab& fine_rhs) final;
113 
114  void getEBFluxes (const Vector<MultiFab*>& a_flux,
115  const Vector<MultiFab*>& a_sol) const override;
116 
117  void applyRobinBCTermsCoeffs ();
118 
119 #if defined(AMREX_USE_HYPRE) && (AMREX_SPACEDIM > 1)
120  [[nodiscard]] std::unique_ptr<Hypre> makeHypre (Hypre::Interface hypre_interface) const override;
121 #endif
122 
123 #ifdef AMREX_USE_PETSC
124  [[nodiscard]] std::unique_ptr<PETScABecLap> makePETSc () const override;
125 #endif
126 
127  Real m_a_scalar = std::numeric_limits<Real>::quiet_NaN();
128  Real m_b_scalar = std::numeric_limits<Real>::quiet_NaN();
131 
132  bool m_scalars_set = false;
133  bool m_acoef_set = false;
134 
135 protected:
136 
137  int m_ncomp = 1;
138 
139  bool m_needs_update = true;
140 
141  Location m_beta_loc; // Location of coefficients: face centers or face centroids
142  Location m_phi_loc; // Location of solution variable: cell centers or cell centroids
143 
145 
148 
150 
151  mutable int m_is_eb_inhomog;
152 
153  //
154  // functions
155  //
156  bool isEBDirichlet () const noexcept { return m_eb_phi[0] != nullptr; }
157 
158  void averageDownCoeffsSameAmrLevel (int amrlev, Vector<MultiFab>& a,
160  const Vector<MultiFab*>& b_eb);
161  void averageDownCoeffs ();
162  void averageDownCoeffsToCoarseAmrLevel (int flev);
163 
164  [[nodiscard]] bool supportRobinBC () const noexcept override { return true; }
165 };
166 
167 }
168 
169 #endif
Array4< Real > fine
Definition: AMReX_InterpFaceRegister.cpp:90
Array4< Real const > crse
Definition: AMReX_InterpFaceRegister.cpp:92
A Fortran Array of REALs.
Definition: AMReX_FArrayBox.H:229
Interface
Definition: AMReX_Hypre.H:21
Definition: AMReX_MFIter.H:57
Definition: AMReX_MLCellABecLap.H:13
bool needsUpdate() const override
Does it need update if it's reused?
Definition: AMReX_MLCellABecLap.H:46
typename MLLinOpT< MF >::Location Location
Definition: AMReX_MLCellABecLap.H:19
typename MLLinOpT< MF >::BCMode BCMode
Definition: AMReX_MLCellLinOp.H:28
typename MLLinOpT< MF >::StateMode StateMode
Definition: AMReX_MLCellLinOp.H:29
Definition: AMReX_MLEBABecLap.H:16
~MLEBABecLap() override
bool isCrossStencil() const override
Definition: AMReX_MLEBABecLap.H:74
MLEBABecLap(MLEBABecLap &&)=delete
void prepareForSolve() override
Definition: AMReX_MLEBABecLap.cpp:686
int m_is_eb_inhomog
Definition: AMReX_MLEBABecLap.H:151
MLEBABecLap(const MLEBABecLap &)=delete
bool m_acoef_set
Definition: AMReX_MLEBABecLap.H:133
Vector< Vector< std::unique_ptr< MultiFab > > > m_eb_b_coeffs
Definition: AMReX_MLEBABecLap.H:147
void setPhiOnCentroid()
Definition: AMReX_MLEBABecLap.cpp:103
Real getBScalar() const final
Definition: AMReX_MLEBABecLap.H:96
Real m_b_scalar
Definition: AMReX_MLEBABecLap.H:128
MultiFab const * getACoeffs(int amrlev, int mglev) const final
Definition: AMReX_MLEBABecLap.H:97
Vector< std::unique_ptr< MultiFab > > m_eb_phi
Definition: AMReX_MLEBABecLap.H:146
bool supportRobinBC() const noexcept override
Definition: AMReX_MLEBABecLap.H:164
void Fapply(int amrlev, int mglev, MultiFab &out, const MultiFab &in) const final
Definition: AMReX_MLEBABecLap_F.cpp:12
void restriction(int, int, MultiFab &crse, MultiFab &fine) const final
Definition: AMReX_MLEBABecLap.cpp:950
void applyRobinBCTermsCoeffs()
Definition: AMReX_MLEBABecLap.cpp:1297
void setACoeffs(int amrlev, const MultiFab &alpha)
Definition: AMReX_MLEBABecLap.cpp:125
void normalize(int amrlev, int mglev, MultiFab &mf) const final
Definition: AMReX_MLEBABecLap.cpp:857
Vector< Vector< MultiFab > > m_a_coeffs
Definition: AMReX_MLEBABecLap.H:129
void averageDownSolutionRHS(int camrlev, MultiFab &crse_sol, MultiFab &crse_rhs, const MultiFab &fine_sol, const MultiFab &fine_rhs) final
Definition: AMReX_MLEBABecLap.cpp:1005
void interpolation(int amrlev, int fmglev, MultiFab &fine, const MultiFab &crse) const final
Definition: AMReX_MLEBABecLap.cpp:958
MLEBABecLap & operator=(const MLEBABecLap &)=delete
void applyBC(int amrlev, int mglev, MultiFab &in, BCMode bc_mode, StateMode s_mode, const MLMGBndry *bndry=nullptr, bool skip_fillboundary=false) const final
Definition: AMReX_MLEBABecLap.cpp:1015
void averageDownCoeffs()
Definition: AMReX_MLEBABecLap.cpp:607
bool isEBDirichlet() const noexcept
Definition: AMReX_MLEBABecLap.H:156
void averageDownCoeffsSameAmrLevel(int amrlev, Vector< MultiFab > &a, Vector< Array< MultiFab, AMREX_SPACEDIM > > &b, const Vector< MultiFab * > &b_eb)
Definition: AMReX_MLEBABecLap.cpp:632
void Fsmooth(int amrlev, int mglev, MultiFab &sol, const MultiFab &rhs, int redblack) const final
Definition: AMReX_MLEBABecLap_F.cpp:157
void getEBFluxes(const Vector< MultiFab * > &a_flux, const Vector< MultiFab * > &a_sol) const override
Definition: AMReX_MLEBABecLap.cpp:1224
Vector< int > m_is_singular
Definition: AMReX_MLEBABecLap.H:149
std::unique_ptr< FabFactory< FArrayBox > > makeFactory(int amrlev, int mglev) const final
Definition: AMReX_MLEBABecLap.cpp:34
void setEBDirichlet(int amrlev, const MultiFab &phi, const MultiFab &beta)
Definition: AMReX_MLEBABecLap.cpp:190
Vector< Vector< iMultiFab > > m_cc_mask
Definition: AMReX_MLEBABecLap.H:144
int m_ncomp
Definition: AMReX_MLEBABecLap.H:137
void apply(int amrlev, int mglev, MultiFab &out, MultiFab &in, BCMode bc_mode, StateMode s_mode, const MLMGBndry *bndry=nullptr) const override
Definition: AMReX_MLEBABecLap.cpp:1179
void averageDownCoeffsToCoarseAmrLevel(int flev)
Definition: AMReX_MLEBABecLap.cpp:663
bool needsUpdate() const override
Does it need update if it's reused?
Definition: AMReX_MLEBABecLap.H:67
void setEBHomogDirichlet(int amrlev, const MultiFab &beta)
Definition: AMReX_MLEBABecLap.cpp:400
void setBCoeffs(int amrlev, const Array< MultiFab const *, AMREX_SPACEDIM > &beta)
Definition: AMReX_MLEBABecLap.H:49
Array< MultiFab const *, AMREX_SPACEDIM > getBCoeffs(int amrlev, int mglev) const final
Definition: AMReX_MLEBABecLap.H:99
bool m_needs_update
Definition: AMReX_MLEBABecLap.H:139
void FFlux(int amrlev, const MFIter &mfi, const Array< FArrayBox *, AMREX_SPACEDIM > &flux, const FArrayBox &sol, Location loc, int face_only=0) const final
Definition: AMReX_MLEBABecLap_F.cpp:322
Real getAScalar() const final
Definition: AMReX_MLEBABecLap.H:95
MLEBABecLap()=default
Real m_a_scalar
Definition: AMReX_MLEBABecLap.H:127
int getNComp() const override
Return number of components.
Definition: AMReX_MLEBABecLap.H:65
void compGrad(int amrlev, const Array< MultiFab *, AMREX_SPACEDIM > &grad, MultiFab &sol, Location loc) const final
Definition: AMReX_MLEBABecLap.cpp:731
std::unique_ptr< MLLinOp > makeNLinOp(int) const final
Definition: AMReX_MLEBABecLap.H:102
bool isSingular(int amrlev) const override
Is it singular on given AMR level?
Definition: AMReX_MLEBABecLap.H:84
void setScalars(Real a, Real b)
Definition: AMReX_MLEBABecLap.cpp:109
Vector< Vector< Array< MultiFab, AMREX_SPACEDIM > > > m_b_coeffs
Definition: AMReX_MLEBABecLap.H:130
void update() override
Update for reuse.
Definition: AMReX_MLEBABecLap.cpp:1188
void define(const Vector< Geometry > &a_geom, const Vector< BoxArray > &a_grids, const Vector< DistributionMapping > &a_dmap, const LPInfo &a_info, const Vector< EBFArrayBoxFactory const * > &a_factory, int a_ncomp=1)
Definition: AMReX_MLEBABecLap.cpp:44
bool isBottomSingular() const override
Is the bottom of MG singular?
Definition: AMReX_MLEBABecLap.H:85
void setBCoeffs(int amrlev, const Array< MultiFab const *, AMREX_SPACEDIM > &beta, Location a_beta_loc)
Definition: AMReX_MLEBABecLap.cpp:141
Location m_phi_loc
Definition: AMReX_MLEBABecLap.H:142
bool m_scalars_set
Definition: AMReX_MLEBABecLap.H:132
Location m_beta_loc
Definition: AMReX_MLEBABecLap.H:141
Definition: AMReX_MLMGBndry.H:12
A collection (stored as an array) of FArrayBox objects.
Definition: AMReX_MultiFab.H:38
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition: AMReX_Vector.H:27
Definition: AMReX_Amr.cpp:49
std::array< T const *, AMREX_SPACEDIM > GetArrOfConstPtrs(const std::array< T, AMREX_SPACEDIM > &a) noexcept
Definition: AMReX_Array.H:864
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
std::unique_ptr< Hypre > makeHypre(const BoxArray &grids, const DistributionMapping &dmap, const Geometry &geom, MPI_Comm comm_, Hypre::Interface interface, const iMultiFab *overset_mask)
Definition: AMReX_Hypre.cpp:12
std::array< T, N > Array
Definition: AMReX_Array.H:24
Definition: AMReX_MLLinOp.H:35