Block-Structured AMR Software Framework
AMReX_MLCurlCurl.H
Go to the documentation of this file.
1 #ifndef AMREX_ML_CURL_CURL_H_
2 #define AMREX_ML_CURL_CURL_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_MLLinOp.H>
6 #include <AMReX_MLCurlCurl_K.H>
7 
8 namespace amrex {
9 
25  : public MLLinOpT<Array<MultiFab,3> >
26 {
27 public:
29  using RT = typename MLLinOpT<MF>::RT;
30  using BCType = typename MLLinOpT<MF>::BCType;
31  using BCMode = typename MLLinOpT<MF>::BCMode;
33  using Location = typename MLLinOpT<MF>::Location;
34 
35  MLCurlCurl () = default;
36  MLCurlCurl (const Vector<Geometry>& a_geom,
37  const Vector<BoxArray>& a_grids,
38  const Vector<DistributionMapping>& a_dmap,
39  const LPInfo& a_info = LPInfo());
40 
41  void define (const Vector<Geometry>& a_geom,
42  const Vector<BoxArray>& a_grids,
43  const Vector<DistributionMapping>& a_dmap,
44  const LPInfo& a_info = LPInfo());
45 
46  void setScalars (RT a_alpha, RT a_beta) noexcept;
47 
49  void setBeta (const Vector<Array<MultiFab const*,3>>& a_bcoefs);
50 
53  void prepareRHS (Vector<MF*> const& rhs) const;
54 
55  void setDirichletNodesToZero (int amrlev, int mglev, MF& a_mf) const override;
56 
57  [[nodiscard]] std::string name () const override {
58  return std::string("curl of curl");
59  }
60 
61  bool setUsePCG (bool flag) { return std::exchange(m_use_pcg, flag); }
62 
63  void setLevelBC (int amrlev, const MF* levelbcdata,
64  const MF* robinbc_a = nullptr,
65  const MF* robinbc_b = nullptr,
66  const MF* robinbc_f = nullptr) override;
67 
68  void restriction (int amrlev, int cmglev, MF& crse, MF& fine) const override;
69 
70  void interpolation (int amrlev, int fmglev, MF& fine, const MF& crse) const override;
71 
72  void apply (int amrlev, int mglev, MF& out, MF& in, BCMode bc_mode,
73  StateMode s_mode, const MLMGBndryT<MF>* bndry=nullptr) const override;
74 
75  void smooth (int amrlev, int mglev, MF& sol, const MF& rhs,
76  bool skip_fillboundary=false) const override;
77 
78  void solutionResidual (int amrlev, MF& resid, MF& x, const MF& b,
79  const MF* crse_bcdata=nullptr) override;
80 
81  void correctionResidual (int amrlev, int mglev, MF& resid, MF& x,
82  const MF& b, BCMode bc_mode,
83  const MF* crse_bcdata=nullptr) override;
84 
85  void prepareForSolve () override;
86 
87  [[nodiscard]] bool isSingular (int /*amrlev*/) const override { return false; }
88  [[nodiscard]] bool isBottomSingular () const override { return false; }
89 
90  RT xdoty (int amrlev, int mglev, const MF& x, const MF& y, bool local) const override;
91 
92  [[nodiscard]] RT normInf (int amrlev, MF const& mf, bool local) const override;
93 
94  void averageDownAndSync (Vector<MF>& sol) const override;
95 
96  [[nodiscard]] IntVect getNGrowVectRestriction () const override {
97  return IntVect(1);
98  }
99 
100  void make (Vector<Vector<MF> >& mf, IntVect const& ng) const override;
101 
102  [[nodiscard]] MF make (int amrlev, int mglev, IntVect const& ng) const override;
103 
104  [[nodiscard]] MF makeAlias (MF const& mf) const override;
105 
106  [[nodiscard]] MF makeCoarseMG (int amrlev, int mglev, IntVect const& ng) const override;
107 
108  [[nodiscard]] MF makeCoarseAmr (int famrlev, IntVect const& ng) const override;
109 
110 // public for cuda
111 
112  void smooth4 (int amrlev, int mglev, MF& sol, MF const& rhs, int color) const;
113 
114  void compresid (int amrlev, int mglev, MF& resid, MF const& b) const;
115 
116  void applyPhysBC (int amrlev, int mglev, MultiFab& mf, CurlCurlStateType type) const;
117 
118 private:
119 
120  void applyBC (int amrlev, int mglev, MF& in, CurlCurlStateType type) const;
121 
122  [[nodiscard]] iMultiFab const& getDotMask (int amrlev, int mglev, int idim) const;
123 
124  [[nodiscard]] CurlCurlDirichletInfo getDirichletInfo (int amrlev, int mglev) const;
125  [[nodiscard]] CurlCurlSymmetryInfo getSymmetryInfo (int amrlev, int mglev) const;
126 
127  RT m_alpha = std::numeric_limits<RT>::lowest();
128  RT m_beta = std::numeric_limits<RT>::lowest();
129 
131 #if (AMREX_SPACEDIM == 3)
132  {IntVect(0,1,1), IntVect(1,0,1), IntVect(1,1,0)};
133 #else
134  {IntVect(0,1), IntVect(1,0), IntVect(1,1)};
135 #endif
136 
138  static constexpr int m_ncomp = 1;
139  Vector<Vector<std::unique_ptr<Gpu::DeviceScalar
142  bool m_use_pcg = false;
143 };
144 
145 }
146 
147 #endif
Array4< Real > fine
Definition: AMReX_InterpFaceRegister.cpp:90
Array4< Real const > crse
Definition: AMReX_InterpFaceRegister.cpp:92
Definition: AMReX_LUSolver.H:16
curl (alpha curl E) + beta E = rhs
Definition: AMReX_MLCurlCurl.H:26
void setLevelBC(int amrlev, const MF *levelbcdata, const MF *robinbc_a=nullptr, const MF *robinbc_b=nullptr, const MF *robinbc_f=nullptr) override
Definition: AMReX_MLCurlCurl.cpp:155
void applyBC(int amrlev, int mglev, MF &in, CurlCurlStateType type) const
Definition: AMReX_MLCurlCurl.cpp:630
void smooth(int amrlev, int mglev, MF &sol, const MF &rhs, bool skip_fillboundary=false) const override
Definition: AMReX_MLCurlCurl.cpp:312
void restriction(int amrlev, int cmglev, MF &crse, MF &fine) const override
Definition: AMReX_MLCurlCurl.cpp:162
Vector< Vector< std::unique_ptr< Gpu::DeviceScalar< LUSolver< AMREX_SPACEDIM *2, RT > > > > > m_lusolver
Definition: AMReX_MLCurlCurl.H:140
void setDirichletNodesToZero(int amrlev, int mglev, MF &a_mf) const override
Definition: AMReX_MLCurlCurl.cpp:103
void interpolation(int amrlev, int fmglev, MF &fine, const MF &crse) const override
Definition: AMReX_MLCurlCurl.cpp:195
bool isSingular(int) const override
Is it singular on given AMR level?
Definition: AMReX_MLCurlCurl.H:87
CurlCurlSymmetryInfo getSymmetryInfo(int amrlev, int mglev) const
Definition: AMReX_MLCurlCurl.cpp:826
CurlCurlDirichletInfo getDirichletInfo(int amrlev, int mglev) const
Definition: AMReX_MLCurlCurl.cpp:792
typename MLLinOpT< MF >::BCType BCType
Definition: AMReX_MLCurlCurl.H:30
typename MLLinOpT< MF >::RT RT
Definition: AMReX_MLCurlCurl.H:29
RT normInf(int amrlev, MF const &mf, bool local) const override
Definition: AMReX_MLCurlCurl.cpp:556
void prepareForSolve() override
Definition: AMReX_MLCurlCurl.cpp:454
Vector< Vector< Array< std::unique_ptr< iMultiFab >, 3 > > > m_dotmask
Definition: AMReX_MLCurlCurl.H:137
Vector< Vector< Array< std::unique_ptr< MultiFab >, 3 > > > m_bcoefs
Definition: AMReX_MLCurlCurl.H:141
void averageDownAndSync(Vector< MF > &sol) const override
Definition: AMReX_MLCurlCurl.cpp:561
void solutionResidual(int amrlev, MF &resid, MF &x, const MF &b, const MF *crse_bcdata=nullptr) override
Definition: AMReX_MLCurlCurl.cpp:390
void setScalars(RT a_alpha, RT a_beta) noexcept
Definition: AMReX_MLCurlCurl.cpp:36
void make(Vector< Vector< MF > > &mf, IntVect const &ng) const override
Definition: AMReX_MLCurlCurl.cpp:573
RT m_beta
Definition: AMReX_MLCurlCurl.H:128
void smooth4(int amrlev, int mglev, MF &sol, MF const &rhs, int color) const
Definition: AMReX_MLCurlCurl.cpp:328
IntVect getNGrowVectRestriction() const override
Definition: AMReX_MLCurlCurl.H:96
static constexpr int m_ncomp
Definition: AMReX_MLCurlCurl.H:138
MLCurlCurl()=default
RT m_alpha
Definition: AMReX_MLCurlCurl.H:127
bool m_use_pcg
Definition: AMReX_MLCurlCurl.H:142
std::string name() const override
Definition: AMReX_MLCurlCurl.H:57
iMultiFab const & getDotMask(int amrlev, int mglev, int idim) const
Definition: AMReX_MLCurlCurl.cpp:781
void define(const Vector< Geometry > &a_geom, const Vector< BoxArray > &a_grids, const Vector< DistributionMapping > &a_dmap, const LPInfo &a_info=LPInfo())
Definition: AMReX_MLCurlCurl.cpp:13
typename MLLinOpT< MF >::Location Location
Definition: AMReX_MLCurlCurl.H:33
void correctionResidual(int amrlev, int mglev, MF &resid, MF &x, const MF &b, BCMode bc_mode, const MF *crse_bcdata=nullptr) override
Definition: AMReX_MLCurlCurl.cpp:399
typename MLLinOpT< MF >::StateMode StateMode
Definition: AMReX_MLCurlCurl.H:32
void compresid(int amrlev, int mglev, MF &resid, MF const &b) const
Definition: AMReX_MLCurlCurl.cpp:408
void apply(int amrlev, int mglev, MF &out, MF &in, BCMode bc_mode, StateMode s_mode, const MLMGBndryT< MF > *bndry=nullptr) const override
Definition: AMReX_MLCurlCurl.cpp:226
MF makeAlias(MF const &mf) const override
Definition: AMReX_MLCurlCurl.cpp:591
RT xdoty(int amrlev, int mglev, const MF &x, const MF &y, bool local) const override
Definition: AMReX_MLCurlCurl.cpp:541
bool setUsePCG(bool flag)
Definition: AMReX_MLCurlCurl.H:61
MF makeCoarseMG(int amrlev, int mglev, IntVect const &ng) const override
Definition: AMReX_MLCurlCurl.cpp:601
typename MLLinOpT< MF >::BCMode BCMode
Definition: AMReX_MLCurlCurl.H:31
void prepareRHS(Vector< MF * > const &rhs) const
Definition: AMReX_MLCurlCurl.cpp:94
Array< MultiFab, 3 > MF
Definition: AMReX_MLCurlCurl.H:28
void applyPhysBC(int amrlev, int mglev, MultiFab &mf, CurlCurlStateType type) const
Definition: AMReX_MLCurlCurl.cpp:648
void setBeta(const Vector< Array< MultiFab const *, 3 >> &a_bcoefs)
This is needed only if there is variable beta coefficient.
Definition: AMReX_MLCurlCurl.cpp:43
Array< IntVect, 3 > m_etype
Definition: AMReX_MLCurlCurl.H:134
MF makeCoarseAmr(int famrlev, IntVect const &ng) const override
Definition: AMReX_MLCurlCurl.cpp:616
bool isBottomSingular() const override
Is the bottom of MG singular?
Definition: AMReX_MLCurlCurl.H:88
Definition: AMReX_MLLinOp.H:98
LinOpBCType BCType
Definition: AMReX_MLLinOp.H:111
typename FabDataType< MF >::value_type RT
Definition: AMReX_MLLinOp.H:109
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_iMultiFab.H:32
Definition: AMReX_Amr.cpp:49
CurlCurlStateType
Definition: AMReX_MLCurlCurl_K.H:310
IntVectND< AMREX_SPACEDIM > IntVect
Definition: AMReX_BaseFwd.H:30
std::array< T, N > Array
Definition: AMReX_Array.H:24
Definition: AMReX_MLCurlCurl_K.H:185
Definition: AMReX_MLCurlCurl_K.H:248
Definition: AMReX_GpuMemory.H:56
Definition: AMReX_MLLinOp.H:35
StateMode
Definition: AMReX_MLLinOp.H:86
BCMode
Definition: AMReX_MLLinOp.H:85
Location
Definition: AMReX_MLLinOp.H:87