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.queryarr("type.xloylozlo",str)) { // 3D Corner
455 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YLO_ZLO][i] = bcmap[str[i]]; }
456 if (pp.queryarr("type.xloylozhi",str)) { // 3D Corner
457 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YLO_ZHI][i] = bcmap[str[i]]; }
458 if (pp.queryarr("type.xloyhizlo",str)) { // 3D Corner
459 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YHI_ZLO][i] = bcmap[str[i]]; }
460 if (pp.queryarr("type.xloyhizhi",str)) { // 3D Corner
461 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YHI_ZHI][i] = bcmap[str[i]]; }
462 if (pp.queryarr("type.xhiylozlo",str)) { // 3D Corner
463 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YLO_ZLO][i] = bcmap[str[i]]; }
464 if (pp.queryarr("type.xhiylozhi",str)) { // 3D Corner
465 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YLO_ZHI][i] = bcmap[str[i]]; }
466 if (pp.queryarr("type.xhiyhizlo",str)) { // 3D Corner
467 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YHI_ZLO][i] = bcmap[str[i]]; }
468 if (pp.queryarr("type.xhiyhizhi",str)) { // 3D Corner
469 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YHI_ZHI][i] = bcmap[str[i]]; }
470 #else // ignore unused inputs
471 pp.ignore("type.xloylozlo"); pp.ignore("type.xloylozhi"); pp.ignore("type.xloyhizlo"); pp.ignore("type.xloyhizhi");
472 pp.ignore("type.xhiylozlo"); pp.ignore("type.xhiylozhi"); pp.ignore("type.xhiyhizlo"); pp.ignore("type.xhiyhizhi");
473 #endif
474
475 #if AMREX_SPACEDIM==3
476 if (pp.queryarr("type.ylozlo",str)) { // 3D Edge
477 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YLO_ZLO][i] = bcmap[str[i]]; }
478 if (pp.queryarr("type.ylozhi",str)) { // 3D Edge
479 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YLO_ZHI][i] = bcmap[str[i]]; }
480 if (pp.queryarr("type.yhizlo",str)) { // 3D Edge
481 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YHI_ZLO][i] = bcmap[str[i]]; }
482 if (pp.queryarr("type.yhizhi",str)) { // 3D Edge
483 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YHI_ZHI][i] = bcmap[str[i]]; }
484 if (pp.queryarr("type.zloxlo",str)) { // 3D Edge
485 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZLO_XLO][i] = bcmap[str[i]]; }
486 if (pp.queryarr("type.zloxhi",str)) { // 3D Edge
487 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZLO_XHI][i] = bcmap[str[i]]; }
488 if (pp.queryarr("type.zhixlo",str)) { // 3D Edge
489 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZHI_XLO][i] = bcmap[str[i]]; }
490 if (pp.queryarr("type.zhixhi",str)) { // 3D Edge
491 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZHI_XHI][i] = bcmap[str[i]]; }
492 #else
493 pp.ignore("type.ylozlo"); pp.ignore("type.ylozhi"); pp.ignore("type.yhizlo"); pp.ignore("type.yhizhi");
494 pp.ignore("type.zloxlo"); pp.ignore("type.zloxhi"); pp.ignore("type.zhixlo"); pp.ignore("type.zhixhi");
495 #endif
496
497 if (pp.queryarr("type.xloylo",str)) { // 3D Edge / 2D Corner
498 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YLO][i] = bcmap[str[i]]; }
499 if (pp.queryarr("type.xloyhi",str)) { // 3D Edge / 2D Corner
500 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YHI][i] = bcmap[str[i]]; }
501 if (pp.queryarr("type.xhiylo",str)) { // 3D Edge / 2D Corner
502 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YLO][i] = bcmap[str[i]]; }
503 if (pp.queryarr("type.xhiyhi",str)) { // 3D Edge / 2D Corner
504 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YHI][i] = bcmap[str[i]]; }
505
506 if (pp.queryarr("type.xlo",str)) { // 3D Face / 2D Edge
507 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO][i] = bcmap[str[i]]; }
508 if (pp.queryarr("type.xhi",str)) { // 3D Face / 2D Edge
509 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI][i] = bcmap[str[i]]; }
510 if (pp.queryarr("type.ylo",str)) { // 3D Face / 2D Edge
511 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YLO][i] = bcmap[str[i]]; }
512 if (pp.queryarr("type.yhi",str)) { // 3D Face / 2D Edge
513 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YHI][i] = bcmap[str[i]]; }
514 #if AMREX_SPACEDIM==3
515 if (pp.queryarr("type.zlo",str)) { // 3D Face
516 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZLO][i] = bcmap[str[i]]; }
517 if (pp.queryarr("type.zhi",str)) { // 3D Face
518 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZHI][i] = bcmap[str[i]]; }
519 #else
520 pp.ignore("type.zlo"); pp.ignore("type.zhi");
521 #endif
522
523 // VALS
524 //std::vector<Set::Scalar> val;
525 std::vector<std::string> val;
526
527 #if AMREX_SPACEDIM==3
528 if (pp.queryarr("val.xloylozlo",val)) { // 3D Corner
529 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YLO_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
530 if (pp.queryarr("val.xloylozhi",val)) { // 3D Corner
531 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YLO_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
532 if (pp.queryarr("val.xloyhizlo",val)) { // 3D Corner
533 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YHI_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
534 if (pp.queryarr("val.xloyhizhi",val)) { // 3D Corner
535 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YHI_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
536 if (pp.queryarr("val.xhiylozlo",val)) { // 3D Corner
537 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YLO_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
538 if (pp.queryarr("val.xhiylozhi",val)) { // 3D Corner
539 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YLO_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
540 if (pp.queryarr("val.xhiyhizlo",val)) { // 3D Corner
541 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YHI_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
542 if (pp.queryarr("val.xhiyhizhi",val)) { // 3D Corner
543 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YHI_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
544 #else
545 pp.ignore("val.xloylozlo"); pp.ignore("val.xloylozhi"); pp.ignore("val.xloyhizlo"); pp.ignore("val.xloyhizhi");
546 pp.ignore("val.xhiylozlo"); pp.ignore("val.xhiylozhi"); pp.ignore("val.xhiyhizlo"); pp.ignore("val.xhiyhizhi");
547 #endif
548
549 #if AMREX_SPACEDIM==3
550 if (pp.queryarr("val.ylozlo",val)) { // 3D Edge
551 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YLO_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
552 if (pp.queryarr("val.ylozhi",val)) { // 3D Edge
553 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YLO_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
554 if (pp.queryarr("val.yhizlo",val)) { // 3D Edge
555 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YHI_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
556 if (pp.queryarr("val.yhizhi",val)) { // 3D Edge
557 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YHI_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
558 if (pp.queryarr("val.zloxlo",val)) { // 3D Edge
559 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZLO_XLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
560 if (pp.queryarr("val.zloxhi",val)) { // 3D Edge
561 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZLO_XHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
562 if (pp.queryarr("val.zhixlo",val)) { // 3D Edge
563 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZHI_XLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
564 if (pp.queryarr("val.zhixhi",val)) { // 3D Edge
565 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZHI_XHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
566 #else
567 pp.ignore("val.ylozlo"); pp.ignore("val.ylozhi"); pp.ignore("val.yhizlo"); pp.ignore("val.yhizhi");
568 pp.ignore("val.zloxlo"); pp.ignore("val.zloxhi"); pp.ignore("val.zhixlo"); pp.ignore("val.zhixhi");
569 #endif
570 if (pp.queryarr("val.xloylo",val)) { // 3D Edge / 2D Corner
571 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
572 if (pp.queryarr("val.xloyhi",val)) { // 3D Edge / 2D Corner
573 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
574 if (pp.queryarr("val.xhiylo",val)) { // 3D Edge / 2D Corner
575 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
576 if (pp.queryarr("val.xhiyhi",val)) { // 3D Edge / 2D Corner
577 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
578
579 if (pp.queryarr("val.xlo",val)) { // 3D Face / 2D Edge
580 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
581 if (pp.queryarr("val.xhi",val)) { // 3D Face / 2D Edge
582 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
583 if (pp.queryarr("val.ylo",val)) { // 3D Face / 2D Edge
584 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
585 if (pp.queryarr("val.yhi",val)) { // 3D Face / 2D Edge
586 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
587 #if AMREX_SPACEDIM==3
588 if (pp.queryarr("val.zlo",val)) { // 3D Face
589 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
590 if (pp.queryarr("val.zhi",val)) { // 3D Face
591 for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
592 #else
593 pp.ignore("val.zlo"); pp.ignore("val.zhi");
594 #endif
595
596 }
597};
598}
599}
600}
601#endif
SQRT2INV
#define SQRT2INV
Definition Elastic.H:80
SQRT3INV
#define SQRT3INV
Definition Elastic.H:79
pp_queryclass
#define pp_queryclass(...)
Definition ParmParse.H:130
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::ignore
void ignore(std::string name, const std::source_location &location=std::source_location::current())
Definition ParmParse.H:276
IO::ParmParse::queryarr
int queryarr(std::string name, std::vector< T > &value, const std::source_location &location=std::source_location::current())
Definition ParmParse.H:1024
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:18
Set::Scalar
amrex::Real Scalar
Definition Base.H:19
Set::Vector
Eigen::Matrix< amrex::Real, AMREX_SPACEDIM, 1 > Vector
Definition Base.H:21
Set::Matrix
Eigen::Matrix< amrex::Real, AMREX_SPACEDIM, AMREX_SPACEDIM > Matrix
Definition Base.H:24
Util::Abort
AMREX_GPU_HOST_DEVICE void Abort(const char *msg)
Definition Util.cpp:406
Generated by doxygen 1.9.8