Line data Source code
1 : //
2 : // This is the most commonly used standard boundary condition implementation.
3 : // The name "Constant" refers to the invariance of the BC value or character along each face
4 : // of the simulation domain; however, you may cause a change in the value with time by using a
5 : // :ref:`Numeric::Interpolator::Linear` string.
6 : //
7 : // The BC on each face is specified with a :code:`type` string that specifies the nature of
8 : // the BC, and a :code:`val` that specifies the value.
9 : // By default, all types are Dirichlet and all values are 0.0.
10 : // The complete list of BC types are below:
11 : //
12 : // `Dirichlet <https://en.wikipedia.org/wiki/Dirichlet_boundary_condition>`_ boundary condition
13 : // types can be specified with :code:`Dirichlet`, :code:`dirichlet`, :code:`EXT_DIR`.
14 : // `Neumann <https://en.wikipedia.org/wiki/Neumann_boundary_condition>`_ boundary condition
15 : // types are specified with :code:`Neumann`, :code:`neumann`.
16 : // `Periodic <https://en.wikipedia.org/wiki/Periodic_boundary_conditions>`_ boundary conditions
17 : // can be specified with :code:`Periodic`, :code:`periodic`, :code:`INT_DIR`.
18 : // **Important**: Ensure that your geometry is specified to be periodic using :code:`geometry.is_periodic`.
19 : // For instance, if periodic in x, set to :code:`1 0 0`
20 : //
21 : // The BC values can be specified either as a number (e.g. :code:`1.0`) or using a linear interpolator
22 : // string (e.g. :code:`(0,1:2.3,-0.1)`).
23 : //
24 : // The number of values and types must be either 0 (for defaults), 1 (to set the same for all field components)
25 : // or N (where N=number of field components).
26 : //
27 :
28 : #ifndef BC_CONSTANT_H_
29 : #define BC_CONSTANT_H_
30 :
31 : #include <AMReX_ParallelDescriptor.H>
32 : #include <AMReX_ParmParse.H>
33 : #include <AMReX_BCRec.H>
34 : #include <AMReX_PhysBCFunct.H>
35 : #include <AMReX_Array.H>
36 :
37 : #include "Set/Set.H"
38 : #include "BC/BC.H"
39 : #include "Numeric/Interpolator/Linear.H"
40 :
41 : namespace BC
42 : {
43 : class Constant
44 : : public BC<Set::Scalar>
45 : {
46 :
47 : #if AMREX_SPACEDIM==2
48 : enum Face {
49 : XLO, YLO, XHI, YHI,
50 : INT
51 : };
52 : #elif AMREX_SPACEDIM==3
53 : enum Face {
54 : XLO, YLO, ZLO, XHI, YHI, ZHI, // 6
55 : INT
56 : };
57 : #endif
58 :
59 : public:
60 : static constexpr const char* name = "constant";
61 :
62 : //Constant (amrex::Vector<amrex::Geometry> &_geom) : geom(_geom) {};
63 12 : Constant(int a_ncomp, Unit a_unit = Unit::Less()) : m_ncomp(a_ncomp), unit(a_unit) {};
64 36 : Constant(int a_ncomp, IO::ParmParse& pp, std::string name) : m_ncomp(a_ncomp)
65 : {
66 36 : pp.queryclass(name, *this);
67 36 : };
68 3 : Constant(int a_ncomp, Unit a_unit, IO::ParmParse& pp, std::string name) :
69 3 : m_ncomp(a_ncomp), unit(a_unit)
70 : {
71 3 : pp.queryclass(name, *this);
72 3 : };
73 : Constant(int ncomp, amrex::Vector<std::string> bc_hi_str,
74 : amrex::Vector<std::string> bc_lo_str,
75 : AMREX_D_DECL(amrex::Vector<amrex::Real> _bc_lo_1,
76 : amrex::Vector<amrex::Real> _bc_lo_2,
77 : amrex::Vector<amrex::Real> _bc_lo_3),
78 : AMREX_D_DECL(amrex::Vector<amrex::Real> _bc_hi_1,
79 : amrex::Vector<amrex::Real> _bc_hi_2,
80 : amrex::Vector<amrex::Real> _bc_hi_3));
81 :
82 12 : static Constant ZeroNeumann(int ncomp = 1, Unit unit = Unit::Less())
83 : {
84 12 : Constant ret = Constant(ncomp, unit);
85 :
86 60 : for (int d = 0; d < m_nfaces; d++)
87 : {
88 48 : ret.m_bc_type[d].clear();
89 120 : for (int n = 0; n < ncomp; n++)
90 : {
91 72 : ret.m_bc_type[d].push_back((int)amrex::LinOpBCType::Neumann);
92 72 : ret.m_bc_val[d].push_back (0.0);;
93 : }
94 : }
95 :
96 12 : return ret;
97 0 : }
98 :
99 : static Constant ZeroDirichlet(int ncomp = 1, Unit unit = Unit::Less())
100 : {
101 : Constant ret = Constant(ncomp, unit);
102 :
103 : for (int d = 0; d < m_nfaces; d++)
104 : for (int n = 0; n < ncomp; n++)
105 : {
106 : ret.m_bc_type[n][d] = (int)amrex::LinOpBCType::Dirichlet;
107 : ret.m_bc_val [n][d] = 0.0;
108 : }
109 :
110 : return ret;
111 : }
112 :
113 :
114 90 : virtual ~Constant() {};
115 :
116 : virtual void FillBoundary(amrex::BaseFab<Set::Scalar>& in, const amrex::Box& box,
117 : int ngrow, int dcomp, int ncomp, amrex::Real time,
118 : Orientation face = Orientation::All,
119 : const amrex::Mask* mask = nullptr) override;
120 :
121 : using BC::FillBoundary;
122 :
123 : amrex::BCRec GetBCRec() override;
124 : virtual amrex::Array<int, AMREX_SPACEDIM> IsPeriodic() override;
125 : virtual amrex::Periodicity Periodicity() const override;
126 : virtual amrex::Periodicity Periodicity(const amrex::Box& b) override;
127 :
128 :
129 :
130 : template<class T>
131 : const amrex::Array<amrex::Array<T, AMREX_SPACEDIM>, 2> GetBCTypes()
132 : {
133 : return { {{AMREX_D_DECL((T)m_bc_type[Face::XLO][0],(T)m_bc_type[Face::YLO][0],(T)m_bc_type[Face::ZLO][0])},
134 : {AMREX_D_DECL((T)m_bc_type[Face::XLO][0],(T)m_bc_type[Face::YLO][0],(T)m_bc_type[Face::ZLO][0])}} };
135 : }
136 :
137 :
138 : private:
139 : #if AMREX_SPACEDIM==2
140 : static const int m_nfaces = 4;
141 : #elif AMREX_SPACEDIM==3
142 : static const int m_nfaces = 6;
143 : #endif
144 :
145 : unsigned int m_ncomp = 0;
146 : Unit unit = Unit::Less();
147 :
148 : //int bc_lo[BL_SPACEDIM];
149 : //int bc_hi[BL_SPACEDIM];
150 : //amrex::Vector<amrex::Real> AMREX_D_DECL(bc_lo_1, bc_lo_2, bc_lo_3);
151 : //amrex::Vector<amrex::Real> AMREX_D_DECL(bc_hi_1, bc_hi_2, bc_hi_3);
152 :
153 : std::array<std::vector<int>, m_nfaces> m_bc_type;
154 : std::array<std::vector<Numeric::Interpolator::Linear<Set::Scalar>>, m_nfaces> m_bc_val;
155 :
156 : public:
157 39 : static void Parse(Constant& value, IO::ParmParse& pp)
158 : {
159 39 : std::map<std::string, int> bcmap;
160 78 : bcmap["BOGUS_BC"] = amrex::BCType::mathematicalBndryTypes::bogus;
161 78 : bcmap["INT_DIR"] = amrex::BCType::mathematicalBndryTypes::int_dir;
162 78 : bcmap["REFLECT_ODD"] = amrex::BCType::mathematicalBndryTypes::reflect_odd;
163 78 : bcmap["INT_DIR"] = amrex::BCType::mathematicalBndryTypes::int_dir;
164 78 : bcmap["REFLECT_EVEN"] = amrex::BCType::mathematicalBndryTypes::reflect_even;
165 78 : bcmap["FOEXTRAP"] = amrex::BCType::mathematicalBndryTypes::foextrap;
166 78 : bcmap["EXT_DIR"] = amrex::BCType::mathematicalBndryTypes::ext_dir;
167 78 : bcmap["HOEXTRAP"] = amrex::BCType::mathematicalBndryTypes::hoextrap;
168 78 : bcmap["Interior"] = amrex::BCType::mathematicalBndryTypes::int_dir;
169 78 : bcmap["Inflow"] = amrex::BCType::mathematicalBndryTypes::ext_dir;
170 78 : bcmap["Outflow"] = amrex::BCType::mathematicalBndryTypes::foextrap;
171 78 : bcmap["Symmetry"] = amrex::BCType::mathematicalBndryTypes::reflect_even;
172 78 : bcmap["SlipWall"] = amrex::BCType::mathematicalBndryTypes::ext_dir;
173 78 : bcmap["NoSlipWall"] = amrex::BCType::mathematicalBndryTypes::ext_dir;
174 : // From <AMReX_LO_BCTYPES.H>
175 78 : bcmap["interior"] = (int)amrex::LinOpBCType::interior;
176 78 : bcmap["Dirichlet"] = (int)amrex::LinOpBCType::Dirichlet;
177 78 : bcmap["dirichlet"] = (int)amrex::LinOpBCType::Dirichlet;
178 78 : bcmap["Neumann"] = (int)amrex::LinOpBCType::Neumann;
179 78 : bcmap["NEUMANN"] = (int)amrex::LinOpBCType::Neumann;
180 78 : bcmap["neumann"] = (int)amrex::LinOpBCType::Neumann;
181 78 : bcmap["reflect_odd"] = (int)amrex::LinOpBCType::reflect_odd;
182 78 : bcmap["Marshak"] = (int)amrex::LinOpBCType::Marshak;
183 78 : bcmap["SanchezPomraning"] = (int)amrex::LinOpBCType::SanchezPomraning;
184 78 : bcmap["inflow"] = (int)amrex::LinOpBCType::inflow;
185 78 : bcmap["Periodic"] = (int)amrex::LinOpBCType::Periodic;
186 39 : bcmap["periodic"] = (int)amrex::LinOpBCType::Periodic;
187 :
188 :
189 :
190 39 : value.m_bc_type[Face::XLO].clear(); value.m_bc_val[Face::XLO].clear();
191 39 : value.m_bc_type[Face::XHI].clear(); value.m_bc_val[Face::XHI].clear();
192 39 : value.m_bc_type[Face::YLO].clear(); value.m_bc_val[Face::YLO].clear();
193 39 : value.m_bc_type[Face::YHI].clear(); value.m_bc_val[Face::YHI].clear();
194 : #if AMREX_SPACEDIM == 3
195 0 : value.m_bc_type[Face::ZLO].clear(); value.m_bc_val[Face::ZLO].clear();
196 0 : value.m_bc_type[Face::ZHI].clear(); value.m_bc_val[Face::ZHI].clear();
197 : #endif
198 :
199 : // TYPES
200 :
201 39 : std::vector<std::string> str;
202 156 : pp.queryarr_default("type.xlo", str,{"dirichlet"}); // BC type on the lower x edge (2d) face (3d)
203 79 : for (unsigned int i = 0; i < str.size(); i++) if (!bcmap.count(str[i])) Util::Exception(INFO, "Invalid BC: ", str[i]);
204 71 : if (str.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_type[Face::XLO].push_back(bcmap[str[i]]);
205 8 : else if (str.size() == 1) value.m_bc_type[Face::XLO].resize(value.m_ncomp, bcmap[str[0]]);
206 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC type args: expected ", value.m_ncomp, " or 1 but got ", str.size());
207 156 : pp.queryarr_default("type.xhi", str,{"dirichlet"}); // BC type on the upper x edge (2d) face (3d)
208 79 : for (unsigned int i = 0; i < str.size(); i++) if (!bcmap.count(str[i])) Util::Exception(INFO, "Invalid BC: ", str[i]);
209 71 : if (str.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_type[Face::XHI].push_back(bcmap[str[i]]);
210 8 : else if (str.size() == 1) value.m_bc_type[Face::XHI].resize(value.m_ncomp, bcmap[str[0]]);
211 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC type args: expected ", value.m_ncomp, " or 1 but got ", str.size());
212 156 : pp.queryarr_default("type.ylo", str,{"dirichlet"}); // BC type on the lower y edge (2d) face (3d)
213 79 : for (unsigned int i = 0; i < str.size(); i++) if (!bcmap.count(str[i])) Util::Exception(INFO, "Invalid BC: ", str[i]);
214 71 : if (str.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_type[Face::YLO].push_back(bcmap[str[i]]);
215 8 : else if (str.size() == 1) value.m_bc_type[Face::YLO].resize(value.m_ncomp, bcmap[str[0]]);
216 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC type args: expected ", value.m_ncomp, " or 1 but got ", str.size());
217 156 : pp.queryarr_default("type.yhi", str,{"dirichlet"}); // BC type on the upper y edge (2d) face (3d)
218 79 : for (unsigned int i = 0; i < str.size(); i++) if (!bcmap.count(str[i])) Util::Exception(INFO, "Invalid BC: ", str[i]);
219 71 : if (str.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_type[Face::YHI].push_back(bcmap[str[i]]);
220 8 : else if (str.size() == 1) value.m_bc_type[Face::YHI].resize(value.m_ncomp, bcmap[str[0]]);
221 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC type args: expected ", value.m_ncomp, " or 1 but got ", str.size());
222 156 : pp.queryarr_default("type.zlo", str,{"dirichlet"}); // BC type on the lower z face (processed but ignored in 2d to prevent unused input errors)
223 : #if AMREX_SPACEDIM==3
224 0 : for (unsigned int i = 0; i < str.size(); i++) if (!bcmap.count(str[i])) Util::Exception(INFO, "Invalid BC: ", str[i]);
225 0 : if (str.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_type[Face::ZLO].push_back(bcmap[str[i]]);
226 0 : else if (str.size() == 1) value.m_bc_type[Face::ZLO].resize(value.m_ncomp, bcmap[str[0]]);
227 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC type args: expected ", value.m_ncomp, " or 1 but got ", str.size());
228 : #endif
229 156 : pp.queryarr_default("type.zhi", str,{"dirichlet"}); // BC type on the upper z face (processed but ignored in 2d to prevent unused input errors)
230 : #if AMREX_SPACEDIM==3
231 0 : for (unsigned int i = 0; i < str.size(); i++) if (!bcmap.count(str[i])) Util::Exception(INFO, "Invalid BC: ", str[i]);
232 0 : if (str.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_type[Face::ZHI].push_back(bcmap[str[i]]);
233 0 : else if (str.size() == 1) value.m_bc_type[Face::ZHI].resize(value.m_ncomp, bcmap[str[0]]);
234 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC type args: expected ", value.m_ncomp, " or 1 but got ", str.size());
235 : #endif
236 :
237 : // VALS
238 39 : std::vector<std::string> val;
239 39 : value.m_bc_val[Face::XLO].clear();
240 156 : pp.queryarr_default("val.xlo", val,{"0.0"}); // BC value on the lower x edge (2d) face (3d)
241 39 : if (val.size() == value.m_ncomp)
242 63 : for (unsigned int i = 0; i < value.m_ncomp; i++)
243 32 : value.m_bc_val[Face::XLO].push_back(Numeric::Interpolator::Linear<Set::Scalar>(val[i],Unit::Time(),value.unit));
244 8 : else if (val.size() == 1)
245 8 : value.m_bc_val[Face::XLO].resize(value.m_ncomp, Numeric::Interpolator::Linear<Set::Scalar>(val[0]));
246 0 : else if (val.size() == 0)
247 0 : value.m_bc_val[Face::XLO].resize(value.m_ncomp, 0.0);
248 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC value args: expected ", value.m_ncomp, " or 0 or 1 but got ", val.size());
249 39 : value.m_bc_val[Face::XHI].clear();
250 156 : pp.queryarr_default("val.xhi", val,{"0.0"}); // BC value on the upper x edge (2d) face (3d)
251 71 : if (val.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_val[Face::XHI].push_back(Numeric::Interpolator::Linear<Set::Scalar>(val[i],Unit::Time(),value.unit));
252 8 : else if (val.size() == 1) value.m_bc_val[Face::XHI].resize(value.m_ncomp, Numeric::Interpolator::Linear<Set::Scalar>(val[0]));
253 0 : else if (val.size() == 0) value.m_bc_val[Face::XHI].resize(value.m_ncomp, 0.0);
254 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC value args: expected ", value.m_ncomp, " or 0 or 1 but got ", val.size());
255 39 : value.m_bc_val[Face::YLO].clear();
256 156 : pp.queryarr_default("val.ylo", val,{"0.0"}); // BC value on the lower y edge (2d) face (3d)
257 73 : if (val.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_val[Face::YLO].push_back(Numeric::Interpolator::Linear<Set::Scalar>(val[i],Unit::Time(),value.unit));
258 7 : else if (val.size() == 1) value.m_bc_val[Face::YLO].resize(value.m_ncomp, Numeric::Interpolator::Linear<Set::Scalar>(val[0]));
259 0 : else if (val.size() == 0) value.m_bc_val[Face::YLO].resize(value.m_ncomp, 0.0);
260 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC value args: expected ", value.m_ncomp, " or 0 or 1 but got ", val.size());
261 39 : value.m_bc_val[Face::YHI].clear();
262 156 : pp.queryarr_default("val.yhi", val,{"0.0"}); // BC value on the upper y edge (2d) face (3d)
263 73 : if (val.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_val[Face::YHI].push_back(Numeric::Interpolator::Linear<Set::Scalar>(val[i],Unit::Time(),value.unit));
264 7 : else if (val.size() == 1) value.m_bc_val[Face::YHI].resize(value.m_ncomp, Numeric::Interpolator::Linear<Set::Scalar>(val[0]));
265 0 : else if (val.size() == 0) value.m_bc_val[Face::YHI].resize(value.m_ncomp, 0.0);
266 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC value args: expected ", value.m_ncomp, " or 0 or 1 but got ", val.size());
267 156 : pp.queryarr_default("val.zlo", val,{"0.0"}); // BC value on the lower z face (processed but ignored in 2d to prevent unused input errors)
268 : #if AMREX_SPACEDIM==3
269 0 : value.m_bc_val[Face::ZLO].clear();
270 0 : if (val.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_val[Face::ZLO].push_back(Numeric::Interpolator::Linear<Set::Scalar>(val[i],Unit::Time(),value.unit));
271 0 : else if (val.size() == 1) value.m_bc_val[Face::ZLO].resize(value.m_ncomp, Numeric::Interpolator::Linear<Set::Scalar>(val[0]));
272 0 : else if (val.size() == 0) value.m_bc_val[Face::ZLO].resize(value.m_ncomp, 0.0);
273 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC value args: expected ", value.m_ncomp, " or 0 or 1 but got ", val.size());
274 : #endif
275 156 : pp.queryarr_default("val.zhi", val,{"0.0"}); // BC value on the upper z face (processed but ignored in 2d to prevent unused input errors)
276 : #if AMREX_SPACEDIM==3
277 0 : value.m_bc_val[Face::ZHI].clear();
278 0 : if (val.size() == value.m_ncomp) for (unsigned int i = 0; i < value.m_ncomp; i++) value.m_bc_val[Face::ZHI].push_back(Numeric::Interpolator::Linear<Set::Scalar>(val[i],Unit::Time(),value.unit));
279 0 : else if (val.size() == 1) value.m_bc_val[Face::ZHI].resize(value.m_ncomp, Numeric::Interpolator::Linear<Set::Scalar>(val[0]));
280 0 : else if (val.size() == 0) value.m_bc_val[Face::ZHI].resize(value.m_ncomp, 0.0);
281 0 : else Util::Exception(INFO, "Incorrect number of ", pp.prefix(), " BC value args: expected ", value.m_ncomp, " or 0 or 1 but got ", val.size());
282 : #endif
283 39 : }
284 :
285 : };
286 : }
287 : #endif
|