Line data Source code
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 :
25 : namespace BC
26 : {
27 : namespace Operator
28 : {
29 : namespace Elastic
30 : {
31 : class Constant : public Elastic
32 : {
33 : public:
34 : static constexpr const char* name = "constant";
35 :
36 28 : Constant()
37 28 : {
38 : // By default, all boundary conditions are displacement
39 360 : for (int face = 0; face < m_nfaces; face++)
40 1152 : for (int direction = 0; direction < AMREX_SPACEDIM; direction++)
41 : {
42 820 : m_bc_type[face][direction] = Type::Displacement;
43 820 : m_bc_val [face][direction] = 0.0;
44 : }
45 28 : };
46 :
47 17 : Constant(IO::ParmParse &pp, std::string name): Constant()
48 17 : { pp_queryclass(name,*this); }
49 :
50 45 : ~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 784 : Init(amrex::FabArray<amrex::BaseFab<Set::Vector>> *a_rhs,
102 : const amrex::Geometry &a_geom,
103 : bool a_homogeneous = false) const override
104 : {
105 784 : amrex::Box domain(a_geom.growPeriodicDomain(2));
106 1568 : domain.convert(amrex::IntVect::TheNodeVector());
107 784 : const amrex::Dim3 lo= amrex::lbound(domain), hi = amrex::ubound(domain);
108 1701 : for (amrex::MFIter mfi(*a_rhs, amrex::TilingIfNotGPU()); mfi.isValid(); ++mfi)
109 : {
110 917 : amrex::Box bx = mfi.grownnodaltilebox();
111 917 : bx = bx & domain;
112 917 : amrex::Array4<Set::Vector> const& rhs = a_rhs->array(mfi);
113 917 : amrex::ParallelFor (bx,[=] AMREX_GPU_DEVICE(int i, int j, int k) {
114 :
115 821805 : for (int dir = 0; dir < AMREX_SPACEDIM; dir++)
116 : {
117 572870 : Face face = Face::INT;
118 :
119 : #if AMREX_SPACEDIM == 2
120 :
121 347870 : if (i==lo.x && j==lo.y) face = Face::XLO_YLO;
122 347580 : else if (i==lo.x && j==hi.y) face = Face::XLO_YHI;
123 347300 : else if (i==hi.x && j==lo.y) face = Face::XHI_YLO;
124 347020 : else if (i==hi.x && j==hi.y) face = Face::XHI_YHI;
125 :
126 346738 : else if (i==lo.x) face = Face::XLO;
127 343426 : else if (i==hi.x) face = Face::XHI;
128 340574 : else if (j==lo.y) face = Face::YLO;
129 337346 : else if (j==hi.y) face = Face::YHI;
130 :
131 : #elif AMREX_SPACEDIM == 3
132 :
133 225000 : if (i==lo.x && j==lo.y && k==lo.z) face = Face::XLO_YLO_ZLO;
134 223200 : else if (i==lo.x && j==lo.y && k==hi.z) face = Face::XLO_YLO_ZHI;
135 221400 : else if (i==lo.x && j==hi.y && k==lo.z) face = Face::XLO_YHI_ZLO;
136 219600 : else if (i==lo.x && j==hi.y && k==hi.z) face = Face::XLO_YHI_ZHI;
137 217800 : else if (i==hi.x && j==lo.y && k==lo.z) face = Face::XHI_YLO_ZLO;
138 216000 : else if (i==hi.x && j==lo.y && k==hi.z) face = Face::XHI_YLO_ZHI;
139 214200 : else if (i==hi.x && j==hi.y && k==lo.z) face = Face::XHI_YHI_ZLO;
140 212400 : else if (i==hi.x && j==hi.y && k==hi.z) face = Face::XHI_YHI_ZHI;
141 :
142 210600 : else if (j==lo.y && k==lo.z) face = Face::YLO_ZLO;
143 205200 : else if (j==lo.y && k==hi.z) face = Face::YLO_ZHI;
144 199800 : else if (j==hi.y && k==lo.z) face = Face::YHI_ZLO;
145 194400 : else if (j==hi.y && k==hi.z) face = Face::YHI_ZHI;
146 189000 : else if (k==lo.z && i==lo.x) face = Face::ZLO_XLO;
147 183600 : else if (k==lo.z && i==hi.x) face = Face::ZLO_XHI;
148 178200 : else if (k==hi.z && i==lo.x) face = Face::ZHI_XLO;
149 172800 : else if (k==hi.z && i==hi.x) face = Face::ZHI_XHI;
150 167400 : else if (i==lo.x && j==lo.y) face = Face::XLO_YLO;
151 162000 : else if (i==lo.x && j==hi.y) face = Face::XLO_YHI;
152 156600 : else if (i==hi.x && j==lo.y) face = Face::XHI_YLO;
153 151200 : else if (i==hi.x && j==hi.y) face = Face::XHI_YHI;
154 :
155 145800 : else if (i==lo.x) face = Face::XLO;
156 129600 : else if (i==hi.x) face = Face::XHI;
157 113400 : else if (j==lo.y) face = Face::YLO;
158 97200 : else if (j==hi.y) face = Face::YHI;
159 81000 : else if (k==lo.z) face = Face::ZLO;
160 64800 : else if (k==hi.z) face = Face::ZHI;
161 :
162 : #endif
163 :
164 572870 : if (a_geom.isPeriodic(0) && a_geom.isPeriodic(1))
165 : {
166 0 : if (face == Face::XLO || face == Face::XHI ||
167 0 : face == Face::YLO || face == Face::YHI) continue;
168 : }
169 572870 : if (a_geom.isPeriodic(0))
170 : {
171 0 : if (face == Face::XLO || face==Face::XHI) continue;
172 : }
173 572870 : if (a_geom.isPeriodic(1))
174 : {
175 58092 : if (face == Face::YLO || face==Face::YHI) continue;
176 : }
177 : #if AMREX_SPACEDIM == 3
178 225000 : if (a_geom.isPeriodic(2))
179 : {
180 0 : if (face == Face::ZLO || face==Face::ZHI) continue;
181 : }
182 : #endif
183 :
184 571378 : if (!(face == Face::INT))
185 : {
186 188220 : if (a_homogeneous && m_bc_type[face][dir] == Type::Displacement)
187 0 : rhs(i,j,k)(dir) = 0.0;
188 : else
189 376440 : rhs(i,j,k)(dir) = m_bc_val[face][dir](m_time);
190 : }
191 : }
192 248935 : });
193 784 : }
194 784 : }
195 :
196 : using Elastic::Init;
197 : virtual void
198 0 : Init(amrex::MultiFab * a_rhs,
199 : const amrex::Geometry &a_geom,
200 : bool a_homogeneous = false) const override
201 : {
202 0 : amrex::Box domain(a_geom.Domain());
203 0 : domain.convert(amrex::IntVect::TheNodeVector());
204 0 : const amrex::Dim3 lo= amrex::lbound(domain), hi = amrex::ubound(domain);
205 0 : for (amrex::MFIter mfi(*a_rhs, amrex::TilingIfNotGPU()); mfi.isValid(); ++mfi)
206 : {
207 0 : amrex::Box bx = mfi.tilebox();
208 0 : bx.grow(2);
209 0 : bx = bx & domain;
210 0 : amrex::Array4<amrex::Real> const& rhs = a_rhs->array(mfi);
211 0 : amrex::ParallelFor (bx,[=] AMREX_GPU_DEVICE(int i, int j, int k) {
212 :
213 0 : for (int dir = 0; dir < AMREX_SPACEDIM; dir++)
214 : {
215 0 : Face face = Face::INT;
216 :
217 : #if AMREX_SPACEDIM == 2
218 :
219 0 : if (i==lo.x && j==lo.y) face = Face::XLO_YLO;
220 0 : else if (i==lo.x && j==hi.y) face = Face::XLO_YHI;
221 0 : else if (i==hi.x && j==lo.y) face = Face::XHI_YLO;
222 0 : else if (i==hi.x && j==hi.y) face = Face::XHI_YHI;
223 :
224 0 : else if (i==lo.x) face = Face::XLO;
225 0 : else if (i==hi.x) face = Face::XHI;
226 0 : else if (j==lo.y) face = Face::YLO;
227 0 : else if (j==hi.y) face = Face::YHI;
228 :
229 : #elif AMREX_SPACEDIM == 3
230 :
231 0 : if (i==lo.x && j==lo.y && k==lo.z) face = Face::XLO_YLO_ZLO;
232 0 : else if (i==lo.x && j==lo.y && k==hi.z) face = Face::XLO_YLO_ZHI;
233 0 : else if (i==lo.x && j==hi.y && k==lo.z) face = Face::XLO_YHI_ZLO;
234 0 : else if (i==lo.x && j==hi.y && k==hi.z) face = Face::XLO_YHI_ZHI;
235 0 : else if (i==hi.x && j==lo.y && k==lo.z) face = Face::XHI_YLO_ZLO;
236 0 : else if (i==hi.x && j==lo.y && k==hi.z) face = Face::XHI_YLO_ZHI;
237 0 : else if (i==hi.x && j==hi.y && k==lo.z) face = Face::XHI_YHI_ZLO;
238 0 : else if (i==hi.x && j==hi.y && k==hi.z) face = Face::XHI_YHI_ZHI;
239 :
240 0 : else if (j==lo.y && k==lo.z) face = Face::YLO_ZLO;
241 0 : else if (j==lo.y && k==hi.z) face = Face::YLO_ZHI;
242 0 : else if (j==hi.y && k==lo.z) face = Face::YHI_ZLO;
243 0 : else if (j==hi.y && k==hi.z) face = Face::YHI_ZHI;
244 0 : else if (k==lo.z && i==lo.x) face = Face::ZLO_XLO;
245 0 : else if (k==lo.z && i==hi.x) face = Face::ZLO_XHI;
246 0 : else if (k==hi.z && i==lo.x) face = Face::ZHI_XLO;
247 0 : else if (k==hi.z && i==hi.x) face = Face::ZHI_XHI;
248 0 : else if (i==lo.x && j==lo.y) face = Face::XLO_YLO;
249 0 : else if (i==lo.x && j==hi.y) face = Face::XLO_YHI;
250 0 : else if (i==hi.x && j==lo.y) face = Face::XHI_YLO;
251 0 : else if (i==hi.x && j==hi.y) face = Face::XHI_YHI;
252 :
253 0 : else if (i==lo.x) face = Face::XLO;
254 0 : else if (i==hi.x) face = Face::XHI;
255 0 : else if (j==lo.y) face = Face::YLO;
256 0 : else if (j==hi.y) face = Face::YHI;
257 0 : else if (k==lo.z) face = Face::ZLO;
258 0 : else if (k==hi.z) face = Face::ZHI;
259 :
260 : #endif
261 :
262 0 : if (a_geom.isPeriodic(0) && a_geom.isPeriodic(1))
263 : {
264 0 : if (face == Face::XLO || face == Face::XHI ||
265 0 : face == Face::YLO || face == Face::YHI) continue;
266 : }
267 0 : if (a_geom.isPeriodic(0))
268 : {
269 0 : if (face == Face::XLO || face==Face::XHI) continue;
270 : }
271 0 : if (a_geom.isPeriodic(1))
272 : {
273 0 : if (face == Face::YLO || face==Face::YHI) continue;
274 : }
275 :
276 0 : if (!(face == Face::INT))
277 : {
278 0 : if (a_homogeneous && m_bc_type[face][dir] == Type::Displacement)
279 0 : rhs(i,j,k,dir) = 0.0;
280 : else
281 0 : rhs(i,j,k,dir) = m_bc_val[face][dir](m_time);
282 : }
283 : }
284 0 : });
285 0 : }
286 0 : }
287 :
288 : virtual
289 0 : 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 0 : 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 0 : if (m[0] == lo.x && m[1] == lo.y) return m_bc_type[Face::XLO_YLO];
299 0 : if (m[0] == lo.x && m[1] == hi.y) return m_bc_type[Face::XLO_YHI];
300 0 : if (m[0] == hi.x && m[1] == lo.y) return m_bc_type[Face::XHI_YLO];
301 0 : if (m[0] == hi.x && m[1] == hi.y) return m_bc_type[Face::XHI_YHI];
302 0 : if (m[0] == lo.x) return m_bc_type[Face::XLO];
303 0 : if (m[0] == hi.x) return m_bc_type[Face::XHI];
304 0 : if (m[1] == lo.y) return m_bc_type[Face::YLO];
305 0 : if (m[1] == hi.y) return m_bc_type[Face::YHI];
306 :
307 : #elif AMREX_SPACEDIM == 3
308 :
309 0 : if (m[0] == lo.x && m[1] == lo.y && m[2] == lo.z) return m_bc_type[Face::XLO_YLO_ZLO];
310 0 : if (m[0] == lo.x && m[1] == lo.y && m[2] == hi.z) return m_bc_type[Face::XLO_YLO_ZHI];
311 0 : if (m[0] == lo.x && m[1] == hi.y && m[2] == lo.z) return m_bc_type[Face::XLO_YHI_ZLO];
312 0 : if (m[0] == lo.x && m[1] == hi.y && m[2] == hi.z) return m_bc_type[Face::XLO_YHI_ZHI];
313 0 : if (m[0] == hi.x && m[1] == lo.y && m[2] == lo.z) return m_bc_type[Face::XHI_YLO_ZLO];
314 0 : if (m[0] == hi.x && m[1] == lo.y && m[2] == hi.z) return m_bc_type[Face::XHI_YLO_ZHI];
315 0 : if (m[0] == hi.x && m[1] == hi.y && m[2] == lo.z) return m_bc_type[Face::XHI_YHI_ZLO];
316 0 : if (m[0] == hi.x && m[1] == hi.y && m[2] == hi.z) return m_bc_type[Face::XHI_YHI_ZHI];
317 :
318 0 : if (m[1] == lo.y && m[2] == lo.z) return m_bc_type[Face::YLO_ZLO];
319 0 : if (m[1] == lo.y && m[2] == hi.z) return m_bc_type[Face::YLO_ZHI];
320 0 : if (m[1] == hi.y && m[2] == lo.z) return m_bc_type[Face::YHI_ZLO];
321 0 : if (m[1] == hi.y && m[2] == hi.z) return m_bc_type[Face::YHI_ZHI];
322 0 : if (m[2] == lo.z && m[0] == lo.x) return m_bc_type[Face::ZLO_XLO];
323 0 : if (m[2] == lo.z && m[0] == hi.x) return m_bc_type[Face::ZLO_XHI];
324 0 : if (m[2] == hi.z && m[0] == lo.x) return m_bc_type[Face::ZHI_XLO];
325 0 : if (m[2] == hi.z && m[0] == hi.x) return m_bc_type[Face::ZHI_XHI];
326 0 : if (m[0] == lo.x && m[1] == lo.y) return m_bc_type[Face::XLO_YLO];
327 0 : if (m[0] == lo.x && m[1] == hi.y) return m_bc_type[Face::XLO_YHI];
328 0 : if (m[0] == hi.x && m[1] == lo.y) return m_bc_type[Face::XHI_YLO];
329 0 : if (m[0] == hi.x && m[1] == hi.y) return m_bc_type[Face::XHI_YHI];
330 :
331 0 : if (m[0] == lo.x) return m_bc_type[Face::XLO];
332 0 : if (m[0] == hi.x) return m_bc_type[Face::XHI];
333 0 : if (m[1] == lo.y) return m_bc_type[Face::YLO];
334 0 : if (m[1] == hi.y) return m_bc_type[Face::YHI];
335 0 : if (m[2] == lo.z) return m_bc_type[Face::ZLO];
336 0 : if (m[2] == hi.z) return m_bc_type[Face::ZHI];
337 :
338 : #endif
339 :
340 0 : return {AMREX_D_DECL(Type::None,Type::None,Type::None)};
341 : }
342 :
343 :
344 : AMREX_FORCE_INLINE
345 7896071 : 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 7896071 : amrex::IntVect m(AMREX_D_DECL(i,j,k));
357 :
358 : // Corners
359 : #if AMREX_SPACEDIM == 2
360 :
361 2217515 : 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 2016278 : 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 2000288 : 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 1799043 : 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 1513183 : if (m[0] == lo.x) return set(m_bc_type[Face::XLO], u, gradu, sigma, Set::Vector(-1, 0));
367 1163040 : if (m[0] == hi.x) return set(m_bc_type[Face::XHI], u, gradu, sigma, Set::Vector(+1, 0));
368 790149 : if (m[1] == lo.y) return set(m_bc_type[Face::YLO], u, gradu, sigma, Set::Vector( 0,-1));
369 526478 : 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 8799666 : 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 8381304 : 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 8241850 : 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 7823488 : 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 8241850 : 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 7823488 : 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 7684034 : 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 7265672 : 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 6707856 : 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 6279572 : 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 6279572 : 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 5851288 : 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 5423004 : 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 4994720 : 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 4994720 : 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 4566436 : 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 4566436 : 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 4138152 : 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 4138152 : 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 3709868 : 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 3067442 : if (m[0] == lo.x) return set(m_bc_type[Face::XLO], u, gradu, sigma, Set::Vector(-1, 0, 0));
396 2629236 : if (m[0] == hi.x) return set(m_bc_type[Face::XHI], u, gradu, sigma, Set::Vector(+1, 0, 0));
397 2191030 : if (m[1] == lo.y) return set(m_bc_type[Face::YLO], u, gradu, sigma, Set::Vector( 0,-1, 0));
398 1752824 : if (m[1] == hi.y) return set(m_bc_type[Face::YHI], u, gradu, sigma, Set::Vector( 0,+1, 0));
399 1314618 : if (m[2] == lo.z) return set(m_bc_type[Face::ZLO], u, gradu, sigma, Set::Vector( 0, 0,-1));
400 876412 : 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 0 : Util::Abort(INFO,"Boundary condition error");
405 0 : 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 7896071 : Set::Vector f = Set::Vector::Zero();
413 30002785 : for (int i = 0; i < AMREX_SPACEDIM; i++)
414 : {
415 22106714 : if (bc_type[i] == Type::Displacement)
416 9312018 : f(i) = u(i);
417 12794696 : else if (bc_type[i] == Type::Traction)
418 12793214 : f(i) = (sigma*n)(i);
419 1482 : else if (bc_type[i] == Type::Neumann)
420 0 : f(i) = (gradu*n)(i);
421 : else
422 1482 : continue;
423 : }
424 7896071 : return f;
425 : }
426 :
427 : protected:
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 :
438 : public:
439 17 : static void Parse(Constant & value, IO::ParmParse & pp)
440 : {
441 17 : std::map<std::string, Type> bcmap;
442 34 : bcmap["displacement"] = Type::Displacement;
443 34 : bcmap["disp"] = Type::Displacement;
444 34 : bcmap["neumann"] = Type::Neumann;
445 34 : bcmap["traction"] = Type::Traction;
446 34 : bcmap["trac"] = Type::Traction;
447 17 : bcmap["periodic"] = Type::Periodic;
448 :
449 17 : std::vector<std::string> str;
450 :
451 : // TYPES
452 :
453 : #if AMREX_SPACEDIM==3
454 0 : if (pp.contains("type.xloylozlo")) {
455 0 : pp_queryarr("type.xloylozlo",str); // 3D Corner
456 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YLO_ZLO][i] = bcmap[str[i]]; }
457 0 : if (pp.contains("type.xloylozhi")) {
458 0 : pp_queryarr("type.xloylozhi",str); // 3D Corner
459 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YLO_ZHI][i] = bcmap[str[i]]; }
460 0 : if (pp.contains("type.xloyhizlo")) {
461 0 : pp_queryarr("type.xloyhizlo",str); // 3D Corner
462 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YHI_ZLO][i] = bcmap[str[i]]; }
463 0 : if (pp.contains("type.xloyhizhi")) {
464 0 : pp_queryarr("type.xloyhizhi",str); // 3D Corner
465 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YHI_ZHI][i] = bcmap[str[i]]; }
466 0 : if (pp.contains("type.xhiylozlo")) {
467 0 : pp_queryarr("type.xhiylozlo",str); // 3D Corner
468 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YLO_ZLO][i] = bcmap[str[i]]; }
469 0 : if (pp.contains("type.xhiylozhi")) {
470 0 : pp_queryarr("type.xhiylozhi",str); // 3D Corner
471 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YLO_ZHI][i] = bcmap[str[i]]; }
472 0 : if (pp.contains("type.xhiyhizlo")) {
473 0 : pp_queryarr("type.xhiyhizlo",str); // 3D Corner
474 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YHI_ZLO][i] = bcmap[str[i]]; }
475 0 : if (pp.contains("type.xhiyhizhi")) {
476 0 : pp_queryarr("type.xhiyhizhi",str); // 3D Corner
477 0 : 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 136 : pp.ignore("type.xloylozlo"); pp.ignore("type.xloylozhi"); pp.ignore("type.xloyhizlo"); pp.ignore("type.xloyhizhi");
480 136 : 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 0 : if (pp.contains("type.ylozlo")) {
485 0 : pp_queryarr("type.ylozlo",str); // 3D Edge
486 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YLO_ZLO][i] = bcmap[str[i]]; }
487 0 : if (pp.contains("type.ylozhi")) {
488 0 : pp_queryarr("type.ylozhi",str); // 3D Edge
489 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YLO_ZHI][i] = bcmap[str[i]]; }
490 0 : if (pp.contains("type.yhizlo")) {
491 0 : pp_queryarr("type.yhizlo",str); // 3D Edge
492 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YHI_ZLO][i] = bcmap[str[i]]; }
493 0 : if (pp.contains("type.yhizhi")) {
494 0 : pp_queryarr("type.yhizhi",str); // 3D Edge
495 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YHI_ZHI][i] = bcmap[str[i]]; }
496 0 : if (pp.contains("type.zloxlo")) {
497 0 : pp_queryarr("type.zloxlo",str); // 3D Edge
498 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZLO_XLO][i] = bcmap[str[i]]; }
499 0 : if (pp.contains("type.zloxhi")) {
500 0 : pp_queryarr("type.zloxhi",str); // 3D Edge
501 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZLO_XHI][i] = bcmap[str[i]]; }
502 0 : if (pp.contains("type.zhixlo")) {
503 0 : pp_queryarr("type.zhixlo",str); // 3D Edge
504 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZHI_XLO][i] = bcmap[str[i]]; }
505 0 : if (pp.contains("type.zhixhi")) {
506 0 : pp_queryarr("type.zhixhi",str); // 3D Edge
507 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZHI_XHI][i] = bcmap[str[i]]; }
508 : #else
509 136 : pp.ignore("type.ylozlo"); pp.ignore("type.ylozhi"); pp.ignore("type.yhizlo"); pp.ignore("type.yhizhi");
510 136 : pp.ignore("type.zloxlo"); pp.ignore("type.zloxhi"); pp.ignore("type.zhixlo"); pp.ignore("type.zhixhi");
511 : #endif
512 :
513 34 : if (pp.contains("type.xloylo")) {
514 6 : pp_queryarr("type.xloylo",str); // 3D Edge / 2D Corner
515 18 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YLO][i] = bcmap[str[i]]; }
516 34 : if (pp.contains("type.xloyhi")) {
517 6 : pp_queryarr("type.xloyhi",str); // 3D Edge / 2D Corner
518 18 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO_YHI][i] = bcmap[str[i]]; }
519 34 : if (pp.contains("type.xhiylo")) {
520 6 : pp_queryarr("type.xhiylo",str); // 3D Edge / 2D Corner
521 18 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YLO][i] = bcmap[str[i]]; }
522 34 : if (pp.contains("type.xhiyhi")) {
523 6 : pp_queryarr("type.xhiyhi",str); // 3D Edge / 2D Corner
524 18 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI_YHI][i] = bcmap[str[i]]; }
525 :
526 34 : if (pp.contains("type.xlo")) {
527 6 : pp_queryarr("type.xlo",str); // 3D Face / 2D Edge
528 18 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XLO][i] = bcmap[str[i]]; }
529 34 : if (pp.contains("type.xhi")) {
530 6 : pp_queryarr("type.xhi",str); // 3D Face / 2D Edge
531 18 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::XHI][i] = bcmap[str[i]]; }
532 34 : if (pp.contains("type.ylo")) {
533 11 : pp_queryarr("type.ylo",str); // 3D Face / 2D Edge
534 33 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YLO][i] = bcmap[str[i]]; }
535 34 : if (pp.contains("type.yhi")) {
536 11 : pp_queryarr("type.yhi",str); // 3D Face / 2D Edge
537 33 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::YHI][i] = bcmap[str[i]]; }
538 : #if AMREX_SPACEDIM==3
539 0 : if (pp.contains("type.zlo")) {
540 0 : pp_queryarr("type.zlo",str); // 3D Face
541 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZLO][i] = bcmap[str[i]]; }
542 0 : if (pp.contains("type.zhi")) {
543 0 : pp_queryarr("type.zhi",str); // 3D Face
544 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_type[Face::ZHI][i] = bcmap[str[i]]; }
545 : #else
546 51 : pp.ignore("type.zlo"); pp.ignore("type.zhi");
547 : #endif
548 :
549 : // VALS
550 : //std::vector<Set::Scalar> val;
551 17 : std::vector<std::string> val;
552 :
553 : #if AMREX_SPACEDIM==3
554 0 : if (pp.contains("val.xloylozlo")) {
555 0 : pp_queryarr("val.xloylozlo",val); // 3D Corner
556 0 : 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 0 : if (pp.contains("val.xloylozhi")) {
558 0 : pp_queryarr("val.xloylozhi",val); // 3D Corner
559 0 : 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 0 : if (pp.contains("val.xloyhizlo")) {
561 0 : pp_queryarr("val.xloyhizlo",val); // 3D Corner
562 0 : 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 0 : if (pp.contains("val.xloyhizhi")) {
564 0 : pp_queryarr("val.xloyhizhi",val); // 3D Corner
565 0 : 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 0 : if (pp.contains("val.xhiylozlo")) {
567 0 : pp_queryarr("val.xhiylozlo",val); // 3D Corner
568 0 : 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 0 : if (pp.contains("val.xhiylozhi")) {
570 0 : pp_queryarr("val.xhiylozhi",val); // 3D Corner
571 0 : 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 0 : if (pp.contains("val.xhiyhizlo")) {
573 0 : pp_queryarr("val.xhiyhizlo",val); // 3D Corner
574 0 : 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 0 : if (pp.contains("val.xhiyhizhi")) {
576 0 : pp_queryarr("val.xhiyhizhi",val); // 3D Corner
577 0 : 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 136 : pp.ignore("val.xloylozlo"); pp.ignore("val.xloylozhi"); pp.ignore("val.xloyhizlo"); pp.ignore("val.xloyhizhi");
580 136 : 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 0 : if (pp.contains("val.ylozlo")) {
585 0 : pp_queryarr("val.ylozlo",val); // 3D Edge
586 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YLO_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
587 0 : if (pp.contains("val.ylozhi")) {
588 0 : pp_queryarr("val.ylozhi",val); // 3D Edge
589 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YLO_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
590 0 : if (pp.contains("val.yhizlo")) {
591 0 : pp_queryarr("val.yhizlo",val); // 3D Edge
592 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YHI_ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
593 0 : if (pp.contains("val.yhizhi")) {
594 0 : pp_queryarr("val.yhizhi",val); // 3D Edge
595 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YHI_ZHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
596 0 : if (pp.contains("val.zloxlo")) {
597 0 : pp_queryarr("val.zloxlo",val); // 3D Edge
598 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZLO_XLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
599 0 : if (pp.contains("val.zloxhi")) {
600 0 : pp_queryarr("val.zloxhi",val); // 3D Edge
601 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZLO_XHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
602 0 : if (pp.contains("val.zhixlo")) {
603 0 : pp_queryarr("val.zhixlo",val); // 3D Edge
604 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZHI_XLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
605 0 : if (pp.contains("val.zhixhi")) {
606 0 : pp_queryarr("val.zhixhi",val); // 3D Edge
607 0 : 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 136 : pp.ignore("val.ylozlo"); pp.ignore("val.ylozhi"); pp.ignore("val.yhizlo"); pp.ignore("val.yhizhi");
610 136 : pp.ignore("val.zloxlo"); pp.ignore("val.zloxhi"); pp.ignore("val.zhixlo"); pp.ignore("val.zhixhi");
611 : #endif
612 34 : if (pp.contains("val.xloylo")) {
613 0 : pp_queryarr("val.xloylo",val); // 3D Edge / 2D Corner
614 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
615 34 : if (pp.contains("val.xloyhi")) {
616 0 : pp_queryarr("val.xloyhi",val); // 3D Edge / 2D Corner
617 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO_YHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
618 34 : if (pp.contains("val.xhiylo")) {
619 6 : pp_queryarr("val.xhiylo",val); // 3D Edge / 2D Corner
620 18 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI_YLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
621 34 : if (pp.contains("val.xhiyhi")) {
622 6 : pp_queryarr("val.xhiyhi",val); // 3D Edge / 2D Corner
623 18 : 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 34 : if (pp.contains("val.xlo")) {
626 0 : pp_queryarr("val.xlo",val); // 3D Face / 2D Edge
627 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
628 34 : if (pp.contains("val.xhi")) {
629 6 : pp_queryarr("val.xhi",val); // 3D Face / 2D Edge
630 18 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::XHI][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
631 34 : if (pp.contains("val.ylo")) {
632 0 : pp_queryarr("val.ylo",val); // 3D Face / 2D Edge
633 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::YLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
634 34 : if (pp.contains("val.yhi")) {
635 0 : pp_queryarr("val.yhi",val); // 3D Face / 2D Edge
636 0 : 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 0 : if (pp.contains("val.zlo")) {
639 0 : pp_queryarr("val.zlo",val); // 3D Face
640 0 : for (int i = 0; i < AMREX_SPACEDIM; i++) value.m_bc_val[Face::ZLO][i] = Numeric::Interpolator::Linear<Set::Scalar>(val[i]); }
641 0 : if (pp.contains("val.zhi")) {
642 0 : pp_queryarr("val.zhi",val); // 3D Face
643 0 : 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 51 : pp.ignore("val.zlo"); pp.ignore("val.zhi");
646 : #endif
647 :
648 17 : }
649 : };
650 : }
651 : }
652 : }
653 : #endif
|
