Block-Structured AMR Software Framework
 
Loading...
Searching...
No Matches
AMReX_MLEBNodeFDLaplacian.H
Go to the documentation of this file.
1#ifndef AMREX_MLEBNODEFDLAPLACIAN_H_
2#define AMREX_MLEBNODEFDLAPLACIAN_H_
3#include <AMReX_Config.H>
4
5#include <AMReX_Array.H>
6#ifdef AMREX_USE_EB
8#endif
9#include <AMReX_MLNodeLinOp.H>
10
11#include <limits>
12
13namespace amrex {
14
15// Although the class has EB in the name, it works for non-EB build too.
16//
17// del dot (sigma grad phi) = rhs, for non-RZ
18// where phi and rhs are nodal multifab, and sigma is a tensor constant
19// with only diagonal components. The EB is assumed to be Dirichlet.
20//
21// del dot (sigma grad phi) - alpha/r^2 phi = rhs, for RZ where alpha is a
22// scalar constant that is zero by default.
23//
24// New feature: for non-RZ, sigma can also be a single-component
25// cell-centered multifab.
26
28 : public MLNodeLinOp
29{
30public:
31
32 MLEBNodeFDLaplacian () = default;
33
34#ifdef AMREX_USE_EB
36 const Vector<BoxArray>& a_grids,
37 const Vector<DistributionMapping>& a_dmap,
38 const LPInfo& a_info,
39 const Vector<EBFArrayBoxFactory const*>& a_factory);
40#endif
41
43 const Vector<BoxArray>& a_grids,
44 const Vector<DistributionMapping>& a_dmap,
45 const LPInfo& a_info);
46
47 ~MLEBNodeFDLaplacian () override = default;
48
53
54 void setSigma (Array<Real,AMREX_SPACEDIM> const& a_sigma) noexcept;
55
56 void setSigma (int amrlev, MultiFab const& a_sigma);
57
58 void setRZ (bool flag);
59
60 void setAlpha (Real a_alpha);
61
62#ifdef AMREX_USE_EB
63
64 // Phi on EB
65 void setEBDirichlet (Real a_phi_eb);
66 //
67 template <typename F>
68 std::enable_if_t<IsCallableR<Real,F,AMREX_D_DECL(Real,Real,Real)>::value>
69 setEBDirichlet (F const& f);
70
71 void define (const Vector<Geometry>& a_geom,
72 const Vector<BoxArray>& a_grids,
73 const Vector<DistributionMapping>& a_dmap,
74 const LPInfo& a_info,
75 const Vector<EBFArrayBoxFactory const*>& a_factory);
76
77 [[nodiscard]] std::unique_ptr<FabFactory<FArrayBox> > makeFactory (int amrlev, int mglev) const final;
78
79 [[nodiscard]] bool scaleRHS (int amrlev, MultiFab* rhs) const final;
80
81#endif
82
83 void define (const Vector<Geometry>& a_geom,
84 const Vector<BoxArray>& a_grids,
85 const Vector<DistributionMapping>& a_dmap,
86 const LPInfo& a_info);
87
88 [[nodiscard]] std::string name () const override { return std::string("MLEBNodeFDLaplacian"); }
89
90 void restriction (int amrlev, int cmglev, MultiFab& crse, MultiFab& fine) const final;
91 void interpolation (int amrlev, int fmglev, MultiFab& fine, const MultiFab& crse) const final;
92
93 void prepareForSolve () final;
94 void Fapply (int amrlev, int mglev, MultiFab& out, const MultiFab& in) const final;
95 void Fsmooth (int amrlev, int mglev, MultiFab& sol, const MultiFab& rhs) const final;
96 void normalize (int amrlev, int mglev, MultiFab& mf) const final;
97
98 void fixUpResidualMask (int amrlev, iMultiFab& resmsk) final;
99
100 [[nodiscard]] bool isSingular (int) const final { return false; }
101 [[nodiscard]] bool isBottomSingular () const final { return false; }
102
103 void compGrad (int amrlev, const Array<MultiFab*,AMREX_SPACEDIM>& grad,
104 MultiFab& sol, Location /*loc*/) const override;
105
106#if defined(AMREX_USE_HYPRE) && (AMREX_SPACEDIM > 1)
107 void fillIJMatrix (MFIter const& mfi,
109 Array4<int const> const& lid,
110 HypreNodeLap::Int* ncols,
111 HypreNodeLap::Int* cols,
112 Real* mat) const override;
113
114 void fillRHS (MFIter const& mfi,
115 Array4<int const> const& lid,
116 Real* rhs,
117 Array4<Real const> const& bfab) const override;
118#endif
119
120 void postSolve (Vector<MultiFab*> const& sol) const override;
121
122private:
125 bool m_has_sigma_mf = false;
126 Real m_s_phi_eb = std::numeric_limits<Real>::lowest();
128 int m_rz = false;
129 Real m_rz_alpha = 0._rt;
130};
131
132#ifdef AMREX_USE_EB
133
134template <typename F>
135std::enable_if_t<IsCallableR<Real,F,AMREX_D_DECL(Real,Real,Real)>::value>
136MLEBNodeFDLaplacian::setEBDirichlet (F const& f)
137{
139 for (int amrlev = 0; amrlev < m_num_amr_levels; ++amrlev) {
140 auto const* factory = dynamic_cast<EBFArrayBoxFactory const*>(m_factory[amrlev][0].get());
141 if (factory) {
142 Geometry const& geom = m_geom[amrlev][0];
143 auto const problo = geom.ProbLoArray();
144 auto const cellsize = geom.CellSizeArray();
145 m_phi_eb[amrlev].define(amrex::convert(m_grids[amrlev][0],IntVect(1)),
146 m_dmap[amrlev][0], 1, 1);
147 m_phi_eb[amrlev].setVal(0.0);
148 auto const& flags = factory->getMultiEBCellFlagFab();
149 auto const& levset = factory->getLevelSet();
150#ifdef AMREX_USE_OMP
151#pragma omp parallel if (Gpu::notInLaunchRegion())
152#endif
153 for (MFIter mfi(m_phi_eb[amrlev],TilingIfNotGPU()); mfi.isValid(); ++mfi)
154 {
155 const Box& ndbx = mfi.growntilebox();
156 const auto& flag = flags[mfi];
157 if (flag.getType() != FabType::regular) {
158 Array4<Real const> const lstarr = levset.const_array(mfi);
159 Array4<Real> const& phi = m_phi_eb[amrlev].array(mfi);
160 AMREX_HOST_DEVICE_FOR_3D(ndbx, i, j, k,
161 {
162 if (lstarr(i,j,k) >= Real(0.0)) {
163 phi(i,j,k) = f(AMREX_D_DECL(problo[0]+i*cellsize[0],
164 problo[1]+j*cellsize[1],
165 problo[2]+k*cellsize[2]));
166 }
167 });
168 }
169 }
170 }
171 }
172}
173
174#endif
175
176}
177
178#endif
#define AMREX_HOST_DEVICE_FOR_3D(...)
Definition AMReX_GpuLaunchMacrosC.nolint.H:106
Array4< Real > fine
Definition AMReX_InterpFaceRegister.cpp:90
Array4< Real const > crse
Definition AMReX_InterpFaceRegister.cpp:92
GpuArray< Real, AMREX_SPACEDIM > CellSizeArray() const noexcept
Definition AMReX_CoordSys.H:76
Definition AMReX_EBFabFactory.H:24
Rectangular problem domain geometry.
Definition AMReX_Geometry.H:73
GpuArray< Real, AMREX_SPACEDIM > ProbLoArray() const noexcept
Definition AMReX_Geometry.H:186
HYPRE_Int Int
Definition AMReX_HypreNodeLap.H:36
Definition AMReX_MFIter.H:57
bool isValid() const noexcept
Is the iterator valid i.e. is it associated with a FAB?
Definition AMReX_MFIter.H:141
Definition AMReX_MLEBNodeFDLaplacian.H:29
void restriction(int amrlev, int cmglev, MultiFab &crse, MultiFab &fine) const final
Definition AMReX_MLEBNodeFDLaplacian.cpp:167
std::string name() const override
Definition AMReX_MLEBNodeFDLaplacian.H:88
bool isBottomSingular() const final
Is the bottom of MG singular?
Definition AMReX_MLEBNodeFDLaplacian.H:101
MLEBNodeFDLaplacian & operator=(const MLEBNodeFDLaplacian &)=delete
int m_rz
Definition AMReX_MLEBNodeFDLaplacian.H:128
void Fsmooth(int amrlev, int mglev, MultiFab &sol, const MultiFab &rhs) const final
Definition AMReX_MLEBNodeFDLaplacian.cpp:526
void interpolation(int amrlev, int fmglev, MultiFab &fine, const MultiFab &crse) const final
Definition AMReX_MLEBNodeFDLaplacian.cpp:214
void setAlpha(Real a_alpha)
Definition AMReX_MLEBNodeFDLaplacian.cpp:67
MLEBNodeFDLaplacian(const MLEBNodeFDLaplacian &)=delete
Real m_rz_alpha
Definition AMReX_MLEBNodeFDLaplacian.H:129
void setRZ(bool flag)
Definition AMReX_MLEBNodeFDLaplacian.cpp:57
GpuArray< Real, AMREX_SPACEDIM > m_sigma
Definition AMReX_MLEBNodeFDLaplacian.H:123
void Fapply(int amrlev, int mglev, MultiFab &out, const MultiFab &in) const final
Definition AMReX_MLEBNodeFDLaplacian.cpp:403
void define(const Vector< Geometry > &a_geom, const Vector< BoxArray > &a_grids, const Vector< DistributionMapping > &a_dmap, const LPInfo &a_info)
Definition AMReX_MLEBNodeFDLaplacian.cpp:123
Vector< MultiFab > m_phi_eb
Definition AMReX_MLEBNodeFDLaplacian.H:127
void setSigma(Array< Real, AMREX_SPACEDIM > const &a_sigma) noexcept
Definition AMReX_MLEBNodeFDLaplacian.cpp:36
void normalize(int amrlev, int mglev, MultiFab &mf) const final
Definition AMReX_MLEBNodeFDLaplacian.cpp:631
MLEBNodeFDLaplacian(MLEBNodeFDLaplacian &&)=delete
Real m_s_phi_eb
Definition AMReX_MLEBNodeFDLaplacian.H:126
void prepareForSolve() final
Definition AMReX_MLEBNodeFDLaplacian.cpp:258
void compGrad(int amrlev, const Array< MultiFab *, AMREX_SPACEDIM > &grad, MultiFab &sol, Location) const override
Definition AMReX_MLEBNodeFDLaplacian.cpp:643
~MLEBNodeFDLaplacian() override=default
bool m_has_sigma_mf
Definition AMReX_MLEBNodeFDLaplacian.H:125
Vector< Vector< std::unique_ptr< MultiFab > > > m_sigma_mf
Definition AMReX_MLEBNodeFDLaplacian.H:124
void fixUpResidualMask(int amrlev, iMultiFab &resmsk) final
Definition AMReX_MLEBNodeFDLaplacian.cpp:637
void postSolve(Vector< MultiFab * > const &sol) const override
Definition AMReX_MLEBNodeFDLaplacian.cpp:759
bool isSingular(int) const final
Is it singular on given AMR level?
Definition AMReX_MLEBNodeFDLaplacian.H:100
Vector< Vector< std::unique_ptr< FabFactory< FAB > > > > m_factory
Definition AMReX_MLLinOp.H:610
Vector< Vector< BoxArray > > m_grids
Definition AMReX_MLLinOp.H:608
Vector< Vector< DistributionMapping > > m_dmap
Definition AMReX_MLLinOp.H:609
virtual std::unique_ptr< FabFactory< FAB > > makeFactory(int, int) const
Definition AMReX_MLLinOp.H:700
Vector< Vector< Geometry > > m_geom
first Vector is for amr level and second is mg level
Definition AMReX_MLLinOp.H:607
virtual bool scaleRHS(int, MF *) const
scale RHS to fix solvability
Definition AMReX_MLLinOp.H:458
LinOpEnumType::Location Location
Definition AMReX_MLLinOp.H:114
int m_num_amr_levels
Definition AMReX_MLLinOp.H:592
Definition AMReX_MLNodeLinOp.H:16
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
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > convert(const BoxND< dim > &b, const IntVectND< dim > &typ) noexcept
Returns a BoxND with different type.
Definition AMReX_Box.H:1435
IntVectND< AMREX_SPACEDIM > IntVect
Definition AMReX_BaseFwd.H:30
bool TilingIfNotGPU() noexcept
Definition AMReX_MFIter.H:12
std::array< T, N > Array
Definition AMReX_Array.H:24
Definition AMReX_Array4.H:61
Definition AMReX_Array.H:34
Test if a given type T is callable with arguments of type Args...
Definition AMReX_TypeTraits.H:207
Definition AMReX_MLLinOp.H:35