Alamo
Constant.H
Go to the documentation of this file.
1//
2// This BC defines boundary conditions that are constant with respect to space.
3// However they may change in time using the :ref:`Numeric::Interpolator::Linear` method.
4//
5// In 2D, BCs are prescribed for each edge (4) and corner (4) on the boundary, for a total of 8 possible regions.
6// In 3D, BCs are prescribed for each face (6), each edge (12), and each corner (8), for a total of 26 possible regions.
7// Care must be taken to define BCs for edges/corners consistently with the faces/edges, or you will get poor
8// convergence and inaccurate behavior.
9// (See :ref:`BC::Operator::Elastic::TensionTest` for a reduced BC with streamlined load options.)
10//
11// To define a boundary condition, you must define both a "type" and a "value" in each direction.
12// The type can be "displacement", "disp", "neumann", "traction", "trac".
13// The value is the corresponding value.
14// All BCs are, by default, dirichlet (displacement) with a value of zero.
15//
16
17#ifndef BC_OPERATOR_ELASTIC_CONSTANT_H
18#define BC_OPERATOR_ELASTIC_CONSTANT_H
19
20// #include "Operator/Elastic.H"
21#include "IO/ParmParse.H"
22#include "Numeric/Interpolator/Linear.H"
23#include "BC/Operator/Elastic/Elastic.H"
24
25namespace BC
26{
27namespace Operator
28{
29namespace Elastic
30{
31class Constant : public Elastic
32{
33public:
34 static constexpr const char* name = "constant";
35
36 Constant()
37 {
38 // By default, all boundary conditions are displacement
39 for (int face = 0; face < m_nfaces; face++)
40 for (int direction = 0; direction < AMREX_SPACEDIM; direction++)
41 {
42 m_bc_type[face][direction] = Type::Displacement;
43 m_bc_val [face][direction] = 0.0;
44 }
45 };
46
47 Constant(IO::ParmParse &pp, std::string name): Constant()
48 { pp_queryclass(name,*this); }
49
50 ~Constant() {};
51
52
53 void
54 Set(const Face face,
55 const Direction direction,
56 const Type type,
57 const Set::Scalar value)
58 {
59 m_bc_type[face][direction] = type;
60 m_bc_val [face][direction] = value;
61 }
62
63 void
64 Set(const Face /*face*/,
65 const Direction /*direction*/,
66 const Type /*type*/,
67 const Set::Scalar /*value*/,
68 amrex::Vector<amrex::MultiFab *> &/*a_rhs*/,
69 const amrex::Vector<amrex::Geometry> &/*a_geom*/)
70 {
71 Util::Abort(INFO,"This has been depricated");
72 }
73
74 void
75 Set(const Face face,
76 const Direction direction,
77 const Type type,
78 const Set::Scalar value,
79 amrex::Vector<amrex::MultiFab> &a_rhs,
80 const amrex::Vector<amrex::Geometry> &a_geom)
81 {
82 amrex::Vector<amrex::MultiFab *> pa_rhs = amrex::GetVecOfPtrs(a_rhs);
83 Set(face,direction,type,value,pa_rhs,a_geom);
84 }
85
86 void
87 Set(const Face face,
88 const Direction direction,
89 const Type type,
90 const Set::Scalar value,
91 amrex::Vector<std::unique_ptr<amrex::MultiFab> > &a_rhs,
92 const amrex::Vector<amrex::Geometry> &a_geom)
93 {
94 amrex::Vector<amrex::MultiFab *> pa_rhs = amrex::GetVecOfPtrs(a_rhs);
95 Set(face,direction,type,value,pa_rhs,a_geom);
96 }
97
98 //using Elastic::Init;
99
100 virtual void
101 Init(amrex::FabArray<amrex::BaseFab<Set::Vector>> *a_rhs,
102 const amrex::Geometry &a_geom,
103 bool a_homogeneous = false) const override
104 {
105 amrex::Box domain(a_geom.growPeriodicDomain(2));
106 domain.convert(amrex::IntVect::TheNodeVector());
107 const amrex::Dim3 lo= amrex::lbound(domain), hi = amrex::ubound(domain);
108 for (amrex::MFIter mfi(*a_rhs, amrex::TilingIfNotGPU()); mfi.isValid(); ++mfi)
109 {
110 amrex::Box bx = mfi.grownnodaltilebox();
111 bx = bx & domain;
112 amrex::Array4<Set::Vector> const& rhs = a_rhs->array(mfi);
113 amrex::ParallelFor (bx,[=] AMREX_GPU_DEVICE(int i, int j, int k) {
114
115 for (int dir = 0; dir < AMREX_SPACEDIM; dir++)
116 {
117 Face face = Face::INT;
118
119 #if AMREX_SPACEDIM == 2
120
121 if (i==lo.x && j==lo.y) face = Face::XLO_YLO;
122 else if (i==lo.x && j==hi.y) face = Face::XLO_YHI;
123 else if (i==hi.x && j==lo.y) face = Face::XHI_YLO;
124 else if (i==hi.x && j==hi.y) face = Face::XHI_YHI;
125
126 else if (i==lo.x) face = Face::XLO;
127 else if (i==hi.x) face = Face::XHI;
128 else if (j==lo.y) face = Face::YLO;
129 else if (j==hi.y) face = Face::YHI;
130
131 #elif AMREX_SPACEDIM == 3
132
133 if (i==lo.x && j==lo.y && k==lo.z) face = Face::XLO_YLO_ZLO;
134 else if (i==lo.x && j==lo.y && k==hi.z) face = Face::XLO_YLO_ZHI;
135 else if (i==lo.x && j==hi.y && k==lo.z) face = Face::XLO_YHI_ZLO;
136 else if (i==lo.x && j==hi.y && k==hi.z) face = Face::XLO_YHI_ZHI;
137 else if (i==hi.x && j==lo.y && k==lo.z) face = Face::XHI_YLO_ZLO;
138 else if (i==hi.x && j==lo.y && k==hi.z) face = Face::XHI_YLO_ZHI;
139 else if (i==hi.x && j==hi.y && k==lo.z) face = Face::XHI_YHI_ZLO;
140 else if (i==hi.x && j==hi.y && k==hi.z) face = Face::XHI_YHI_ZHI;
141
142 else if (j==lo.y && k==lo.z) face = Face::YLO_ZLO;
143 else if (j==lo.y && k==hi.z) face = Face::YLO_ZHI;
144 else if (j==hi.y && k==lo.z) face = Face::YHI_ZLO;
145 else if (j==hi.y && k==hi.z) face = Face::YHI_ZHI;
146 else if (k==lo.z && i==lo.x) face = Face::ZLO_XLO;
147 else if (k==lo.z && i==hi.x) face = Face::ZLO_XHI;
148 else if (k==hi.z && i==lo.x) face = Face::ZHI_XLO;
149 else if (k==hi.z && i==hi.x) face = Face::ZHI_XHI;
150 else if (i==lo.x && j==lo.y) face = Face::XLO_YLO;
151 else if (i==lo.x && j==hi.y) face = Face::XLO_YHI;
152 else if (i==hi.x && j==lo.y) face = Face::XHI_YLO;
153 else if (i==hi.x && j==hi.y) face = Face::XHI_YHI;
154
155 else if (i==lo.x) face = Face::XLO;
156 else if (i==hi.x) face = Face::XHI;
157 else if (j==lo.y) face = Face::YLO;
158 else if (j==hi.y) face = Face::YHI;
159 else if (k==lo.z) face = Face::ZLO;
160 else if (k==hi.z) face = Face::ZHI;
161
162 #endif
163
164 if (a_geom.isPeriodic(0) && a_geom.isPeriodic(1))
165 {
166 if (face == Face::XLO || face == Face::XHI ||
167 face == Face::YLO || face == Face::YHI) continue;
168 }
169 if (a_geom.isPeriodic(0))
170 {
171 if (face == Face::XLO || face==Face::XHI) continue;
172 }
173 if (a_geom.isPeriodic(1))
174 {
175 if (face == Face::YLO || face==Face::YHI) continue;
176 }
177 #if AMREX_SPACEDIM == 3
178 if (a_geom.isPeriodic(2))
179 {
180 if (face == Face::ZLO || face==Face::ZHI) continue;
181 }
182 #endif
183
184 if (!(face == Face::INT))
185 {
186 if (a_homogeneous && m_bc_type[face][dir] == Type::Displacement)
187 rhs(i,j,k)(dir) = 0.0;
188 else
189 rhs(i,j,k)(dir) = m_bc_val[face][dir](m_time);
190 }
191 }
192 });
193 }
194 }
195
196 using Elastic::Init;
197 virtual void
198 Init(amrex::MultiFab * a_rhs,
199 const amrex::Geometry &a_geom,
200 bool a_homogeneous = false) const override
201 {
202 amrex::Box domain(a_geom.Domain());
203 domain.convert(amrex::IntVect::TheNodeVector());
204 const amrex::Dim3 lo= amrex::lbound(domain), hi = amrex::ubound(domain);
205 for (amrex::MFIter mfi(*a_rhs, amrex::TilingIfNotGPU()); mfi.isValid(); ++mfi)
206 {
207 amrex::Box bx = mfi.tilebox();
208 bx.grow(2);
209 bx = bx & domain;
210 amrex::Array4<amrex::Real> const& rhs = a_rhs->array(mfi);
211 amrex::ParallelFor (bx,[=] AMREX_GPU_DEVICE(int i, int j, int k) {
212
213 for (int dir = 0; dir < AMREX_SPACEDIM; dir++)
214 {
215 Face face = Face::INT;
216
217 #if AMREX_SPACEDIM == 2
218
219 if (i==lo.x && j==lo.y) face = Face::XLO_YLO;
220 else if (i==lo.x && j==hi.y) face = Face::XLO_YHI;
221 else if (i==hi.x && j==lo.y) face = Face::XHI_YLO;
222 else if (i==hi.x && j==hi.y) face = Face::XHI_YHI;
223
224 else if (i==lo.x) face = Face::XLO;
225 else if (i==hi.x) face = Face::XHI;
226 else if (j==lo.y) face = Face::YLO;
227 else if (j==hi.y) face = Face::YHI;
228
229 #elif AMREX_SPACEDIM == 3
230
231 if (i==lo.x && j==lo.y && k==lo.z) face = Face::XLO_YLO_ZLO;
232 else if (i==lo.x && j==lo.y && k==hi.z) face = Face::XLO_YLO_ZHI;
233 else if (i==lo.x && j==hi.y && k==lo.z) face = Face::XLO_YHI_ZLO;
234 else if (i==lo.x && j==hi.y && k==hi.z) face = Face::XLO_YHI_ZHI;
235 else if (i==hi.x && j==lo.y && k==lo.z) face = Face::XHI_YLO_ZLO;
236 else if (i==hi.x && j==lo.y && k==hi.z) face = Face::XHI_YLO_ZHI;
237 else if (i==hi.x && j==hi.y && k==lo.z) face = Face::XHI_YHI_ZLO;
238 else if (i==hi.x && j==hi.y && k==hi.z) face = Face::XHI_YHI_ZHI;
239
240 else if (j==lo.y && k==lo.z) face = Face::YLO_ZLO;
241 else if (j==lo.y && k==hi.z) face = Face::YLO_ZHI;
242 else if (j==hi.y && k==lo.z) face = Face::YHI_ZLO;
243 else if (j==hi.y && k==hi.z) face = Face::YHI_ZHI;
244 else if (k==lo.z && i==lo.x) face = Face::ZLO_XLO;
245 else if (k==lo.z && i==hi.x) face = Face::ZLO_XHI;
246 else if (k==hi.z && i==lo.x) face = Face::ZHI_XLO;
247 else if (k==hi.z && i==hi.x) face = Face::ZHI_XHI;
248 else if (i==lo.x && j==lo.y) face = Face::XLO_YLO;
249 else if (i==lo.x && j==hi.y) face = Face::XLO_YHI;
250 else if (i==hi.x && j==lo.y) face = Face::XHI_YLO;
251 else if (i==hi.x && j==hi.y) face = Face::XHI_YHI;
252
253 else if (i==lo.x) face = Face::XLO;
254 else if (i==hi.x) face = Face::XHI;
255 else if (j==lo.y) face = Face::YLO;
256 else if (j==hi.y) face = Face::YHI;
257 else if (k==lo.z) face = Face::ZLO;
258 else if (k==hi.z) face = Face::ZHI;
259
260 #endif
261
262 if (a_geom.isPeriodic(0) && a_geom.isPeriodic(1))
263 {
264 if (face == Face::XLO || face == Face::XHI ||
265 face == Face::YLO || face == Face::YHI) continue;
266 }
267 if (a_geom.isPeriodic(0))
268 {
269 if (face == Face::XLO || face==Face::XHI) continue;
270 }
271 if (a_geom.isPeriodic(1))
272 {
273 if (face == Face::YLO || face==Face::YHI) continue;
274 }
275
276 if (!(face == Face::INT))
277 {
278 if (a_homogeneous && m_bc_type[face][dir] == Type::Displacement)
279 rhs(i,j,k,dir) = 0.0;
280 else
281 rhs(i,j,k,dir) = m_bc_val[face][dir](m_time);
282 }
283 }
284 });
285 }
286 }
287
288 virtual
289 std::array<Type,AMREX_SPACEDIM> getType (
290 const int &i, const int &j, [[maybe_unused]] const int &k,
291 const amrex::Box &domain) override
292 {
293 amrex::IntVect m(AMREX_D_DECL(i,j,k));
294 const amrex::Dim3 lo = amrex::lbound(domain), hi = amrex::ubound(domain);
295
296 // Corners
297 #if AMREX_SPACEDIM == 2
298 if (m[0] == lo.x && m[1] == lo.y) return m_bc_type[Face::XLO_YLO];
299 if (m[0] == lo.x && m[1] == hi.y) return m_bc_type[Face::XLO_YHI];
300 if (m[0] == hi.x && m[1] == lo.y) return m_bc_type[Face::XHI_YLO];
301 if (m[0] == hi.x && m[1] == hi.y) return m_bc_type[Face::XHI_YHI];
302 if (m[0] == lo.x) return m_bc_type[Face::XLO];
303 if (m[0] == hi.x) return m_bc_type[Face::XHI];
304 if (m[1] == lo.y) return m_bc_type[Face::YLO];
305 if (m[1] == hi.y) return m_bc_type[Face::YHI];
306
307 #elif AMREX_SPACEDIM == 3
308
309 if (m[0] == lo.x && m[1] == lo.y && m[2] == lo.z) return m_bc_type[Face::XLO_YLO_ZLO];
310 if (m[0] == lo.x && m[1] == lo.y && m[2] == hi.z) return m_bc_type[Face::XLO_YLO_ZHI];
311 if (m[0] == lo.x && m[1] == hi.y && m[2] == lo.z) return m_bc_type[Face::XLO_YHI_ZLO];
312 if (m[0] == lo.x && m[1] == hi.y && m[2] == hi.z) return m_bc_type[Face::XLO_YHI_ZHI];
313 if (m[0] == hi.x && m[1] == lo.y && m[2] == lo.z) return m_bc_type[Face::XHI_YLO_ZLO];
314 if (m[0] == hi.x && m[1] == lo.y && m[2] == hi.z) return m_bc_type[Face::XHI_YLO_ZHI];
315 if (m[0] == hi.x && m[1] == hi.y && m[2] == lo.z) return m_bc_type[Face::XHI_YHI_ZLO];
316 if (m[0] == hi.x && m[1] == hi.y && m[2] == hi.z) return m_bc_type[Face::XHI_YHI_ZHI];
317
318 if (m[1] == lo.y && m[2] == lo.z) return m_bc_type[Face::YLO_ZLO];
319 if (m[1] == lo.y && m[2] == hi.z) return m_bc_type[Face::YLO_ZHI];
320 if (m[1] == hi.y && m[2] == lo.z) return m_bc_type[Face::YHI_ZLO];
321 if (m[1] == hi.y && m[2] == hi.z) return m_bc_type[Face::YHI_ZHI];
322 if (m[2] == lo.z && m[0] == lo.x) return m_bc_type[Face::ZLO_XLO];
323 if (m[2] == lo.z && m[0] == hi.x) return m_bc_type[Face::ZLO_XHI];
324 if (m[2] == hi.z && m[0] == lo.x) return m_bc_type[Face::ZHI_XLO];
325 if (m[2] == hi.z && m[0] == hi.x) return m_bc_type[Face::ZHI_XHI];
326 if (m[0] == lo.x && m[1] == lo.y) return m_bc_type[Face::XLO_YLO];
327 if (m[0] == lo.x && m[1] == hi.y) return m_bc_type[Face::XLO_YHI];
328 if (m[0] == hi.x && m[1] == lo.y) return m_bc_type[Face::XHI_YLO];
329 if (m[0] == hi.x && m[1] == hi.y) return m_bc_type[Face::XHI_YHI];
330
331 if (m[0] == lo.x) return m_bc_type[Face::XLO];
332 if (m[0] == hi.x) return m_bc_type[Face::XHI];
333 if (m[1] == lo.y) return m_bc_type[Face::YLO];
334 if (m[1] == hi.y) return m_bc_type[Face::YHI];
335 if (m[2] == lo.z) return m_bc_type[Face::ZLO];
336 if (m[2] == hi.z) return m_bc_type[Face::ZHI];
337
338 #endif
339
340 return {AMREX_D_DECL(Type::None,Type::None,Type::None)};
341 }
342
343
344 AMREX_FORCE_INLINE
345 Set::Vector operator () (const Set::Vector &u,
346 const Set::Matrix &gradu,
347 const Set::Matrix &sigma,
348 const int &i, const int &j, const int &k,
349 const amrex::Box &domain) override
350 {
351 (void)i; (void)j; (void)k; // Suppress "unused variable" warnings
352 //Set::Vector f;
353
354 const amrex::Dim3 lo= amrex::lbound(domain), hi = amrex::ubound(domain);
355
356 amrex::IntVect m(AMREX_D_DECL(i,j,k));
357
358 // Corners
359 #if AMREX_SPACEDIM == 2
360
361 if (m[0] == lo.x && m[1] == lo.y) return set(m_bc_type[Face::XLO_YLO], u, gradu, sigma, Set::Vector(-SQRT2INV, -SQRT2INV));
362 if (m[0] == lo.x && m[1] == hi.y) return set(m_bc_type[Face::XLO_YHI], u, gradu, sigma, Set::Vector(-SQRT2INV, +SQRT2INV));
363 if (m[0] == hi.x && m[1] == lo.y) return set(m_bc_type[Face::XHI_YLO], u, gradu, sigma, Set::Vector(+SQRT2INV, -SQRT2INV));
364 if (m[0] == hi.x && m[1] == hi.y) return set(m_bc_type[Face::XHI_YHI], u, gradu, sigma, Set::Vector(+SQRT2INV, +SQRT2INV));
365
366 if (m[0] == lo.x) return set(m_bc_type[Face::XLO], u, gradu, sigma, Set::Vector(-1, 0));
367 if (m[0] == hi.x) return set(m_bc_type[Face::XHI], u, gradu, sigma, Set::Vector(+1, 0));
368 if (m[1] == lo.y) return set(m_bc_type[Face::YLO], u, gradu, sigma, Set::Vector( 0,-1));
369 if (m[1] == hi.y) return set(m_bc_type[Face::YHI], u, gradu, sigma, Set::Vector( 0,+1));
370
371 #elif AMREX_SPACEDIM == 3
372
373 if (m[0] == lo.x && m[1] == lo.y && m[2] == lo.z) return set(m_bc_type[Face::XLO_YLO_ZLO], u, gradu, sigma, Set::Vector(-SQRT3INV,-SQRT3INV,-SQRT3INV));
374 if (m[0] == lo.x && m[1] == lo.y && m[2] == hi.z) return set(m_bc_type[Face::XLO_YLO_ZHI], u, gradu, sigma, Set::Vector(-SQRT3INV,-SQRT3INV,+SQRT3INV));
375 if (m[0] == lo.x && m[1] == hi.y && m[2] == lo.z) return set(m_bc_type[Face::XLO_YHI_ZLO], u, gradu, sigma, Set::Vector(-SQRT3INV,+SQRT3INV,-SQRT3INV));
376 if (m[0] == lo.x && m[1] == hi.y && m[2] == hi.z) return set(m_bc_type[Face::XLO_YHI_ZHI], u, gradu, sigma, Set::Vector(-SQRT3INV,+SQRT3INV,+SQRT3INV));
377 if (m[0] == hi.x && m[1] == lo.y && m[2] == lo.z) return set(m_bc_type[Face::XHI_YLO_ZLO], u, gradu, sigma, Set::Vector(+SQRT3INV,-SQRT3INV,-SQRT3INV));
378 if (m[0] == hi.x && m[1] == lo.y && m[2] == hi.z) return set(m_bc_type[Face::XHI_YLO_ZHI], u, gradu, sigma, Set::Vector(+SQRT3INV,-SQRT3INV,+SQRT3INV));
379 if (m[0] == hi.x && m[1] == hi.y && m[2] == lo.z) return set(m_bc_type[Face::XHI_YHI_ZLO], u, gradu, sigma, Set::Vector(+SQRT3INV,+SQRT3INV,-SQRT3INV));
380 if (m[0] == hi.x && m[1] == hi.y && m[2] == hi.z) return set(m_bc_type[Face::XHI_YHI_ZHI], u, gradu, sigma, Set::Vector(+SQRT3INV,+SQRT3INV,+SQRT3INV));
381
382 if (m[1] == lo.y && m[2] == lo.z) return set(m_bc_type[Face::YLO_ZLO], u, gradu, sigma, Set::Vector(0, -SQRT2INV,-SQRT2INV));
383 if (m[1] == lo.y && m[2] == hi.z) return set(m_bc_type[Face::YLO_ZHI], u, gradu, sigma, Set::Vector(0, -SQRT2INV,+SQRT2INV));
384 if (m[1] == hi.y && m[2] == lo.z) return set(m_bc_type[Face::YHI_ZLO], u, gradu, sigma, Set::Vector(0, +SQRT2INV,-SQRT2INV));
385 if (m[1] == hi.y && m[2] == hi.z) return set(m_bc_type[Face::YHI_ZHI], u, gradu, sigma, Set::Vector(0, +SQRT2INV,+SQRT2INV));
386 if (m[2] == lo.z && m[0] == lo.x) return set(m_bc_type[Face::ZLO_XLO], u, gradu, sigma, Set::Vector(-SQRT2INV, 0, -SQRT2INV));
387 if (m[2] == lo.z && m[0] == hi.x) return set(m_bc_type[Face::ZLO_XHI], u, gradu, sigma, Set::Vector(+SQRT2INV, 0, -SQRT2INV));
388 if (m[2] == hi.z && m[0] == lo.x) return set(m_bc_type[Face::ZHI_XLO], u, gradu, sigma, Set::Vector(-SQRT2INV, 0, +SQRT2INV));
389 if (m[2] == hi.z && m[0] == hi.x) return set(m_bc_type[Face::ZHI_XHI], u, gradu, sigma, Set::Vector(+SQRT2INV, 0, +SQRT2INV));
390 if (m[0] == lo.x && m[1] == lo.y) return set(m_bc_type[Face::XLO_YLO], u, gradu, sigma, Set::Vector(-SQRT2INV, -SQRT2INV,0 ));
391 if (m[0] == lo.x && m[1] == hi.y) return set(m_bc_type[Face::XLO_YHI], u, gradu, sigma, Set::Vector(-SQRT2INV, +SQRT2INV,0 ));
392 if (m[0] == hi.x && m[1] == lo.y) return set(m_bc_type[Face::XHI_YLO], u, gradu, sigma, Set::Vector(+SQRT2INV, -SQRT2INV,0 ));
393 if (m[0] == hi.x && m[1] == hi.y) return set(m_bc_type[Face::XHI_YHI], u, gradu, sigma, Set::Vector(+SQRT2INV, +SQRT2INV,0 ));
394
395 if (m[0] == lo.x) return set(m_bc_type[Face::XLO], u, gradu, sigma, Set::Vector(-1, 0, 0));
396 if (m[0] == hi.x) return set(m_bc_type[Face::XHI], u, gradu, sigma, Set::Vector(+1, 0, 0));
397 if (m[1] == lo.y) return set(m_bc_type[Face::YLO], u, gradu, sigma, Set::Vector( 0,-1, 0));
398 if (m[1] == hi.y) return set(m_bc_type[Face::YHI], u, gradu, sigma, Set::Vector( 0,+1, 0));
399 if (m[2] == lo.z) return set(m_bc_type[Face::ZLO], u, gradu, sigma, Set::Vector( 0, 0,-1));
400 if (m[2] == hi.z) return set(m_bc_type[Face::ZHI], u, gradu, sigma, Set::Vector( 0, 0,+1));
401
402 #endif
403
404 Util::Abort(INFO,"Boundary condition error");
405 return Set::Vector::Zero();
406 }
407
408 AMREX_FORCE_INLINE
409 Set::Vector set(std::array<Type,AMREX_SPACEDIM> &bc_type,
410 const Set::Vector &u, const Set::Matrix &gradu, const Set::Matrix &sigma, Set::Vector n) const
411 {
412 Set::Vector f = Set::Vector::Zero();
413 for (int i = 0; i < AMREX_SPACEDIM; i++)
414 {
415 if (bc_type[i] == Type::Displacement)
416 f(i) = u(i);
417 else if (bc_type[i] == Type::Traction)
418 f(i) = (sigma*n)(i);
419 else if (bc_type[i] == Type::Neumann)
420 f(i) = (gradu*n)(i);
421 else
422 continue;
423 }
424 return f;
425 }
426
427protected:
428
429 #if AMREX_SPACEDIM==2
430 static const int m_nfaces = 8;
431 #elif AMREX_SPACEDIM==3
432 static const int m_nfaces = 26;
433 #endif
434
435 std::array<std::array<Type, AMREX_SPACEDIM>, m_nfaces> m_bc_type;
436 std::array<std::array<Numeric::Interpolator::Linear<Set::Scalar>,AMREX_SPACEDIM>, m_nfaces> m_bc_val;
437
438public:
439 static void Parse(Constant & value, IO::ParmParse & pp)
440 {
441 std::map<std::string, Type> bcmap;
442 bcmap["displacement"] = Type::Displacement;
443 bcmap["disp"] = Type::Displacement;
444 bcmap["neumann"] = Type::Neumann;
445 bcmap["traction"] = Type::Traction;
446 bcmap["trac"] = Type::Traction;
447 bcmap["periodic"] = Type::Periodic;
448
449 std::vector<std::string> str;
450
451 // TYPES
452
453 #if AMREX_SPACEDIM==3
454 if (pp.contains("type.xloylozlo")) {
455 pp_queryarr("type.xloylozlo",str); // 3D Corner
456 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YLO_ZLO][i] = bcmap[str[i]]; }
457 if (pp.contains("type.xloylozhi")) {
458 pp_queryarr("type.xloylozhi",str); // 3D Corner
459 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YLO_ZHI][i] = bcmap[str[i]]; }
460 if (pp.contains("type.xloyhizlo")) {
461 pp_queryarr("type.xloyhizlo",str); // 3D Corner
462 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YHI_ZLO][i] = bcmap[str[i]]; }
463 if (pp.contains("type.xloyhizhi")) {
464 pp_queryarr("type.xloyhizhi",str); // 3D Corner
465 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YHI_ZHI][i] = bcmap[str[i]]; }
466 if (pp.contains("type.xhiylozlo")) {
467 pp_queryarr("type.xhiylozlo",str); // 3D Corner
468 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YLO_ZLO][i] = bcmap[str[i]]; }
469 if (pp.contains("type.xhiylozhi")) {
470 pp_queryarr("type.xhiylozhi",str); // 3D Corner
471 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YLO_ZHI][i] = bcmap[str[i]]; }
472 if (pp.contains("type.xhiyhizlo")) {
473 pp_queryarr("type.xhiyhizlo",str); // 3D Corner
474 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YHI_ZLO][i] = bcmap[str[i]]; }
475 if (pp.contains("type.xhiyhizhi")) {
476 pp_queryarr("type.xhiyhizhi",str); // 3D Corner
477 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YHI_ZHI][i] = bcmap[str[i]]; }
478 #else // ignore unused inputs
479 pp.ignore("type.xloylozlo"); pp.ignore("type.xloylozhi"); pp.ignore("type.xloyhizlo"); pp.ignore("type.xloyhizhi");
480 pp.ignore("type.xhiylozlo"); pp.ignore("type.xhiylozhi"); pp.ignore("type.xhiyhizlo"); pp.ignore("type.xhiyhizhi");
481 #endif
482
483 #if AMREX_SPACEDIM==3
484 if (pp.contains("type.ylozlo")) {
485 pp_queryarr("type.ylozlo",str); // 3D Edge
486 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YLO_ZLO][i] = bcmap[str[i]]; }
487 if (pp.contains("type.ylozhi")) {
488 pp_queryarr("type.ylozhi",str); // 3D Edge
489 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YLO_ZHI][i] = bcmap[str[i]]; }
490 if (pp.contains("type.yhizlo")) {
491 pp_queryarr("type.yhizlo",str); // 3D Edge
492 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YHI_ZLO][i] = bcmap[str[i]]; }
493 if (pp.contains("type.yhizhi")) {
494 pp_queryarr("type.yhizhi",str); // 3D Edge
495 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YHI_ZHI][i] = bcmap[str[i]]; }
496 if (pp.contains("type.zloxlo")) {
497 pp_queryarr("type.zloxlo",str); // 3D Edge
498 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZLO_XLO][i] = bcmap[str[i]]; }
499 if (pp.contains("type.zloxhi")) {
500 pp_queryarr("type.zloxhi",str); // 3D Edge
501 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZLO_XHI][i] = bcmap[str[i]]; }
502 if (pp.contains("type.zhixlo")) {
503 pp_queryarr("type.zhixlo",str); // 3D Edge
504 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZHI_XLO][i] = bcmap[str[i]]; }
505 if (pp.contains("type.zhixhi")) {
506 pp_queryarr("type.zhixhi",str); // 3D Edge
507 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZHI_XHI][i] = bcmap[str[i]]; }
508 #else
509 pp.ignore("type.ylozlo"); pp.ignore("type.ylozhi"); pp.ignore("type.yhizlo"); pp.ignore("type.yhizhi");
510 pp.ignore("type.zloxlo"); pp.ignore("type.zloxhi"); pp.ignore("type.zhixlo"); pp.ignore("type.zhixhi");
511 #endif
512
513 if (pp.contains("type.xloylo")) {
514 pp_queryarr("type.xloylo",str); // 3D Edge / 2D Corner
515 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YLO][i] = bcmap[str[i]]; }
516 if (pp.contains("type.xloyhi")) {
517 pp_queryarr("type.xloyhi",str); // 3D Edge / 2D Corner
518 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YHI][i] = bcmap[str[i]]; }
519 if (pp.contains("type.xhiylo")) {
520 pp_queryarr("type.xhiylo",str); // 3D Edge / 2D Corner
521 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YLO][i] = bcmap[str[i]]; }
522 if (pp.contains("type.xhiyhi")) {
523 pp_queryarr("type.xhiyhi",str); // 3D Edge / 2D Corner
524 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YHI][i] = bcmap[str[i]]; }
525
526 if (pp.contains("type.xlo")) {
527 pp_queryarr("type.xlo",str); // 3D Face / 2D Edge
528 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO][i] = bcmap[str[i]]; }
529 if (pp.contains("type.xhi")) {
530 pp_queryarr("type.xhi",str); // 3D Face / 2D Edge
531 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI][i] = bcmap[str[i]]; }
532 if (pp.contains("type.ylo")) {
533 pp_queryarr("type.ylo",str); // 3D Face / 2D Edge
534 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YLO][i] = bcmap[str[i]]; }
535 if (pp.contains("type.yhi")) {
536 pp_queryarr("type.yhi",str); // 3D Face / 2D Edge
537 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YHI][i] = bcmap[str[i]]; }
538 #if AMREX_SPACEDIM==3
539 if (pp.contains("type.zlo")) {
540 pp_queryarr("type.zlo",str); // 3D Face
541 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZLO][i] = bcmap[str[i]]; }
542 if (pp.contains("type.zhi")) {
543 pp_queryarr("type.zhi",str); // 3D Face
544 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZHI][i] = bcmap[str[i]]; }
545 #else
546 pp.ignore("type.zlo"); pp.ignore("type.zhi");
547 #endif
548
549 // VALS
550 //std::vector<Set::Scalar> val;
551 std::vector<std::string> val;
552
553 #if AMREX_SPACEDIM==3
554 if (pp.contains("val.xloylozlo")) {
555 pp_queryarr("val.xloylozlo",val); // 3D Corner
556 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YLO_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
557 if (pp.contains("val.xloylozhi")) {
558 pp_queryarr("val.xloylozhi",val); // 3D Corner
559 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YLO_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
560 if (pp.contains("val.xloyhizlo")) {
561 pp_queryarr("val.xloyhizlo",val); // 3D Corner
562 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YHI_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
563 if (pp.contains("val.xloyhizhi")) {
564 pp_queryarr("val.xloyhizhi",val); // 3D Corner
565 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YHI_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
566 if (pp.contains("val.xhiylozlo")) {
567 pp_queryarr("val.xhiylozlo",val); // 3D Corner
568 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YLO_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
569 if (pp.contains("val.xhiylozhi")) {
570 pp_queryarr("val.xhiylozhi",val); // 3D Corner
571 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YLO_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
572 if (pp.contains("val.xhiyhizlo")) {
573 pp_queryarr("val.xhiyhizlo",val); // 3D Corner
574 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YHI_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
575 if (pp.contains("val.xhiyhizhi")) {
576 pp_queryarr("val.xhiyhizhi",val); // 3D Corner
577 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YHI_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
578 #else
579 pp.ignore("val.xloylozlo"); pp.ignore("val.xloylozhi"); pp.ignore("val.xloyhizlo"); pp.ignore("val.xloyhizhi");
580 pp.ignore("val.xhiylozlo"); pp.ignore("val.xhiylozhi"); pp.ignore("val.xhiyhizlo"); pp.ignore("val.xhiyhizhi");
581 #endif
582
583 #if AMREX_SPACEDIM==3
584 if (pp.contains("val.ylozlo")) {
585 pp_queryarr("val.ylozlo",val); // 3D Edge
586 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YLO_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
587 if (pp.contains("val.ylozhi")) {
588 pp_queryarr("val.ylozhi",val); // 3D Edge
589 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YLO_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
590 if (pp.contains("val.yhizlo")) {
591 pp_queryarr("val.yhizlo",val); // 3D Edge
592 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YHI_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
593 if (pp.contains("val.yhizhi")) {
594 pp_queryarr("val.yhizhi",val); // 3D Edge
595 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YHI_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
596 if (pp.contains("val.zloxlo")) {
597 pp_queryarr("val.zloxlo",val); // 3D Edge
598 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZLO_XLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
599 if (pp.contains("val.zloxhi")) {
600 pp_queryarr("val.zloxhi",val); // 3D Edge
601 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZLO_XHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
602 if (pp.contains("val.zhixlo")) {
603 pp_queryarr("val.zhixlo",val); // 3D Edge
604 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZHI_XLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
605 if (pp.contains("val.zhixhi")) {
606 pp_queryarr("val.zhixhi",val); // 3D Edge
607 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZHI_XHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
608 #else
609 pp.ignore("val.ylozlo"); pp.ignore("val.ylozhi"); pp.ignore("val.yhizlo"); pp.ignore("val.yhizhi");
610 pp.ignore("val.zloxlo"); pp.ignore("val.zloxhi"); pp.ignore("val.zhixlo"); pp.ignore("val.zhixhi");
611 #endif
612 if (pp.contains("val.xloylo")) {
613 pp_queryarr("val.xloylo",val); // 3D Edge / 2D Corner
614 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
615 if (pp.contains("val.xloyhi")) {
616 pp_queryarr("val.xloyhi",val); // 3D Edge / 2D Corner
617 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
618 if (pp.contains("val.xhiylo")) {
619 pp_queryarr("val.xhiylo",val); // 3D Edge / 2D Corner
620 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
621 if (pp.contains("val.xhiyhi")) {
622 pp_queryarr("val.xhiyhi",val); // 3D Edge / 2D Corner
623 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
624
625 if (pp.contains("val.xlo")) {
626 pp_queryarr("val.xlo",val); // 3D Face / 2D Edge
627 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
628 if (pp.contains("val.xhi")) {
629 pp_queryarr("val.xhi",val); // 3D Face / 2D Edge
630 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
631 if (pp.contains("val.ylo")) {
632 pp_queryarr("val.ylo",val); // 3D Face / 2D Edge
633 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
634 if (pp.contains("val.yhi")) {
635 pp_queryarr("val.yhi",val); // 3D Face / 2D Edge
636 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
637 #if AMREX_SPACEDIM==3
638 if (pp.contains("val.zlo")) {
639 pp_queryarr("val.zlo",val); // 3D Face
640 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
641 if (pp.contains("val.zhi")) {
642 pp_queryarr("val.zhi",val); // 3D Face
643 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
644 #else
645 pp.ignore("val.zlo"); pp.ignore("val.zhi");
646 #endif
647
648 }
649};
650}
651}
652}
653#endif
SQRT2INV
#define SQRT2INV
Definition Elastic.H:80
SQRT3INV
#define SQRT3INV
Definition Elastic.H:79
pp_queryarr
#define pp_queryarr(...)
Definition ParmParse.H:105
pp_queryclass
#define pp_queryclass(...)
Definition ParmParse.H:109
INFO
#define INFO
Definition Util.H:24
BC::Operator::Elastic::Constant::Set
void Set(const Face, const Direction, const Type, const Set::Scalar, amrex::Vector< amrex::MultiFab * > &, const amrex::Vector< amrex::Geometry > &)
Definition Constant.H:64
BC::Operator::Elastic::Constant::Set
void Set(const Face face, const Direction direction, const Type type, const Set::Scalar value, amrex::Vector< std::unique_ptr< amrex::MultiFab > > &a_rhs, const amrex::Vector< amrex::Geometry > &a_geom)
Definition Constant.H:87
BC::Operator::Elastic::Constant::m_bc_type
std::array< std::array< Type, AMREX_SPACEDIM >, m_nfaces > m_bc_type
Definition Constant.H:435
BC::Operator::Elastic::Constant::set
AMREX_FORCE_INLINE Set::Vector set(std::array< Type, AMREX_SPACEDIM > &bc_type, const Set::Vector &u, const Set::Matrix &gradu, const Set::Matrix &sigma, Set::Vector n) const
Definition Constant.H:409
BC::Operator::Elastic::Constant::Set
void Set(const Face face, const Direction direction, const Type type, const Set::Scalar value, amrex::Vector< amrex::MultiFab > &a_rhs, const amrex::Vector< amrex::Geometry > &a_geom)
Definition Constant.H:75
BC::Operator::Elastic::Constant::Constant
Constant(IO::ParmParse &pp, std::string name)
Definition Constant.H:47
BC::Operator::Elastic::Constant::name
static constexpr const char * name
Definition Constant.H:34
BC::Operator::Elastic::Constant::m_bc_val
std::array< std::array< Numeric::Interpolator::Linear< Set::Scalar >, AMREX_SPACEDIM >, m_nfaces > m_bc_val
Definition Constant.H:436
BC::Operator::Elastic::Constant::Parse
static void Parse(Constant &value, IO::ParmParse &pp)
Definition Constant.H:439
BC::Operator::Elastic::Constant::Init
virtual void Init(amrex::FabArray< amrex::BaseFab< Set::Vector > > *a_rhs, const amrex::Geometry &a_geom, bool a_homogeneous=false) const override
Definition Constant.H:101
BC::Operator::Elastic::Constant::getType
virtual std::array< Type, AMREX_SPACEDIM > getType(const int &i, const int &j, const int &k, const amrex::Box &domain) override
Definition Constant.H:289
BC::Operator::Elastic::Constant::Set
void Set(const Face face, const Direction direction, const Type type, const Set::Scalar value)
Definition Constant.H:54
BC::Operator::Elastic::Constant::operator()
AMREX_FORCE_INLINE Set::Vector operator()(const Set::Vector &u, const Set::Matrix &gradu, const Set::Matrix &sigma, const int &i, const int &j, const int &k, const amrex::Box &domain) override
Definition Constant.H:345
BC::Operator::Elastic::Constant::Init
virtual void Init(amrex::MultiFab *a_rhs, const amrex::Geometry &a_geom, bool a_homogeneous=false) const override
Definition Constant.H:198
BC::Operator::Elastic::Elastic::Init
virtual void Init(amrex::MultiFab *a_rhs, const amrex::Geometry &a_geom, bool a_homogeneous=false) const =0
IO::ParmParse::contains
bool contains(std::string name)
Definition ParmParse.H:173
IO::ParmParse::ignore
void ignore(std::string name)
Definition ParmParse.H:137
BC
Collection of boundary condition (BC) objects.
Definition BC.cpp:5
Operator
Documentation for operator namespace.
Definition Diagonal.cpp:14
Set
A collection of data types and symmetry-reduced data structures.
Definition Base.H:17
Set::Scalar
amrex::Real Scalar
Definition Base.H:18
Set::Vector
Eigen::Matrix< amrex::Real, AMREX_SPACEDIM, 1 > Vector
Definition Base.H:19
Set::Matrix
Eigen::Matrix< amrex::Real, AMREX_SPACEDIM, AMREX_SPACEDIM > Matrix
Definition Base.H:22
Util::Abort
void Abort(const char *msg)
Definition Util.cpp:268
Generated by doxygen 1.9.8