REMORA
Regional Modeling of Oceans Refined Adaptively
Loading...
Searching...
No Matches
REMORA_init_bcs.cpp
Go to the documentation of this file.
1#include <AMReX_Vector.H>
2#include <AMReX_BC_TYPES.H>
3#include <AMReX_ParmParse.H>
4
5#include <REMORA.H>
6
7using namespace amrex;
8
9void REMORA::init_bcs ()
10{
11 const int xvel_bc_idx = xvel_bc();
12 const int yvel_bc_idx = yvel_bc();
13 const int zvel_bc_idx = zvel_bc();
14 const int ubar_bc_idx = ubar_bc();
15 const int vbar_bc_idx = vbar_bc();
16 const int zeta_bc_idx = zeta_bc();
17 const int tke_bc_idx = tke_bc();
18 const int foextrap_periodic_bc_idx = foextrap_periodic_bc();
19 const int foextrap_bc_idx = foextrap_bc();
20 const int u2d_simple_bc_idx = u2d_simple_bc();
21 const int v2d_simple_bc_idx = v2d_simple_bc();
22
23 phys_bc_type.assign(num_bc_vars(), {});
24 m_bc_extdir_vals.assign(num_bc_vars(), {});
25
26 auto uses_velocity_input = [=] (int bcvar_type) noexcept {
27 return bcvar_type == xvel_bc_idx || bcvar_type == yvel_bc_idx || bcvar_type == zvel_bc_idx;
28 };
29
30 auto uses_scalar_input = [=] (int bcvar_type) noexcept {
31 return bcvar_type >= Tracer_comp && bcvar_type < ncons;
32 };
33
34 auto f_set_var_bc = [this, uses_velocity_input, uses_scalar_input, xvel_bc_idx, zeta_bc_idx, ubar_bc_idx, vbar_bc_idx]
35 (ParmParse& pp, int bcvar_type, Orientation ori, std::string bc_type_string) {
36 // const bool requires_file_data =
37 // (bc_type_string == "clamped") || (bc_type_string == "chapman") ||
38 // (bc_type_string == "flather") || (bc_type_string == "orlanski_rad_nudg");
39
40 // if (bcvar_type > Tracer_comp && requires_file_data) {
41 // amrex::Abort("Additional passive scalars do not yet support NetCDF-driven boundary conditions");
42 // }
43
44 if (bc_type_string == "symmetry")
45 {
46 phys_bc_type[bcvar_type][ori] = REMORA_BC::symmetry;
47 phys_bc_need_data[bdy_index[bcvar_type]][ori] = false;
48 domain_bc_type[ori] = "Symmetry";
49 }
50 else if (bc_type_string == "outflow")
51 {
52 phys_bc_type[bcvar_type][ori] = REMORA_BC::outflow;
53 phys_bc_need_data[bdy_index[bcvar_type]][ori] = false;
54 domain_bc_type[ori] = "Outflow";
55 }
56 else if (bc_type_string == "inflow")
57 {
58 phys_bc_type[bcvar_type][ori] = REMORA_BC::inflow;
59 phys_bc_need_data[bdy_index[bcvar_type]][ori] = false;
60 domain_bc_type[ori] = "Inflow";
61
62 if (uses_velocity_input(bcvar_type)) {
63 std::vector<Real> v;
64 pp.getarr("velocity", v, 0, AMREX_SPACEDIM);
65 m_bc_extdir_vals[bcvar_type][ori] = v[bcvar_type - xvel_bc_idx];
66 } else if (uses_scalar_input(bcvar_type)) {
67 Real scalar_in = 0.;
68 if (pp.queryAdd("scalar", scalar_in)) {
69 m_bc_extdir_vals[bcvar_type][ori] = scalar_in;
70 }
71 }
72 }
73 else if (bc_type_string == "noslipwall")
74 {
75 phys_bc_type[bcvar_type][ori] = REMORA_BC::no_slip_wall;
76 phys_bc_need_data[bdy_index[bcvar_type]][ori] = false;
77 domain_bc_type[ori] = "NoSlipWall";
78
79 if (uses_velocity_input(bcvar_type)) {
80 std::vector<Real> v;
81
82 // The values of m_bc_extdir_vals default to 0.
83 // But if we find "velocity" in the inputs file, use those values instead.
84 if (pp.queryarr("velocity", v, 0, AMREX_SPACEDIM))
85 {
86 v[ori.coordDir()] = 0.0_rt;
87 m_bc_extdir_vals[bcvar_type][ori] = v[bcvar_type - xvel_bc_idx];
88 }
89 }
90 }
91 else if (bc_type_string == "slipwall")
92 {
93 phys_bc_type[bcvar_type][ori] = REMORA_BC::slip_wall;
94 phys_bc_need_data[bdy_index[bcvar_type]][ori] = false;
95 domain_bc_type[ori] = "SlipWall";
96 }
97 else if (bc_type_string == "clamped")
98 {
99 phys_bc_type[bcvar_type][ori] = REMORA_BC::clamped;
100 phys_bc_need_data[bdy_index[bcvar_type]][ori] = true;
101 domain_bc_type[ori] = "Clamped";
102 solverChoice.boundary_from_netcdf = true;
103 }
104 else if (bc_type_string == "chapman")
105 {
106 phys_bc_type[bcvar_type][ori] = REMORA_BC::chapman;
107 phys_bc_need_data[bdy_index[bcvar_type]][ori] = true;
108 domain_bc_type[ori] = "Chapman";
109 solverChoice.boundary_from_netcdf = true;
110
111 if (bcvar_type != zeta_bc_idx) {
112 amrex::Abort("Chapman BC can only be applied to zeta");
113 }
114 }
115 else if (bc_type_string == "flather")
116 {
117 phys_bc_type[bcvar_type][ori] = REMORA_BC::flather;
118 phys_bc_need_data[bdy_index[bcvar_type]][ori] = true;
119 domain_bc_type[ori] = "Flather";
120 solverChoice.boundary_from_netcdf = true;
121
122 if (!(bcvar_type == ubar_bc_idx || bcvar_type == vbar_bc_idx)) {
123 amrex::Abort("Flather BC can only be applied to ubar or vbar");
124 }
125 }
126 else if (bc_type_string == "orlanski_rad")
127 {
128 phys_bc_type[bcvar_type][ori] = REMORA_BC::orlanski_rad;
129 phys_bc_need_data[bdy_index[bcvar_type]][ori] = false;
130 domain_bc_type[ori] = "Orlanski Radiation";
131 }
132 else if (bc_type_string == "orlanski_rad_nudg")
133 {
134 phys_bc_type[bcvar_type][ori] = REMORA_BC::orlanski_rad_nudge;
135 phys_bc_need_data[bdy_index[bcvar_type]][ori] = true;
136 domain_bc_type[ori] = "Orlanski Radiation with nudging";
137 solverChoice.boundary_from_netcdf = true;
138 }
139 else if (bc_type_string == "periodic")
140 {
141 if (!geom[0].isPeriodic(ori.coordDir())) {
142 amrex::Abort("Periodic boundary specified in a non-periodic direction");
143 }
144 phys_bc_type[bcvar_type][ori] = REMORA_BC::periodic;
145 phys_bc_need_data[bdy_index[bcvar_type]][ori] = false;
146 domain_bc_type[ori] = "Periodic";
147 }
148 else
149 {
150 phys_bc_type[bcvar_type][ori] = REMORA_BC::undefined;
151 phys_bc_need_data[bdy_index[bcvar_type]][ori] = false;
152 }
153
154 if (geom[0].isPeriodic(ori.coordDir()))
155 {
156 domain_bc_type[ori] = "Periodic";
157 if (phys_bc_type[bcvar_type][ori] == REMORA_BC::undefined)
158 {
159 phys_bc_type[bcvar_type][ori] = REMORA_BC::periodic;
160 phys_bc_need_data[bdy_index[bcvar_type]][ori] = false;
161 } else if (phys_bc_type[bcvar_type][ori] != REMORA_BC::periodic) {
162 amrex::Abort("Wrong BC type for periodic boundary");
163 }
164 }
165
166 if (phys_bc_type[bcvar_type][ori] == REMORA_BC::undefined)
167 {
168 amrex::Print() << "BC Type specified for fac is " << bc_type_string << std::endl;
169 amrex::Abort("This BC type is unknown");
170 }
171 };
172
173 auto f_by_side = [this, &f_set_var_bc] (std::string const& bcid, Orientation ori)
174 {
175 ParmParse pp("remora.bc."+bcid);
176 std::string bc_type_in = "null";
177 // Default z directions to slipwall
178 if (bcid=="zlo" or bcid=="zhi") bc_type_in = "slipwall";
179 pp.queryAdd("type", bc_type_in);
180 std::string bc_type = amrex::toLower(bc_type_in);
181
182 for (int icomp = 0; icomp < num_bc_vars(); ++icomp) {
183 f_set_var_bc(pp, icomp, ori, bc_type);
184 }
185 };
186
187 auto f_by_var = [this, &f_set_var_bc, zvel_bc_idx] (std::string const& varname, int bcvar_type)
188 {
189 amrex::Vector<Orientation> orientations = {Orientation(Direction::x,Orientation::low), Orientation(Direction::y,Orientation::high),Orientation(Direction::x,Orientation::high),Orientation(Direction::y,Orientation::low)}; // west, south, east, north [matches ROMS]
190 std::vector<std::string> bc_types = {"null","null","null","null"};
191 ParmParse pp("remora.bc."+varname);
192 std::string bc_type_in = "null";
193 // default zvel to outflow
194 if (bcvar_type == zvel_bc_idx) {
195 bc_types = {"outflow","outflow","outflow","outflow"};
196 for (int i=0; i<4; i++) {
197 auto ori = orientations[i];
198 if (geom[0].isPeriodic(ori.coordDir())) {
199 bc_types[i] = "periodic";
200 }
201 }
202 }
203 pp.queryarr("type", bc_types);
204 AMREX_ASSERT(bc_types.size() == 4);
205 for (int i=0; i<4; i++) {
206 std::string bc_type = amrex::toLower(bc_types[i]);
207 auto ori = orientations[i];
208 f_set_var_bc(pp, bcvar_type, ori, bc_type);
209 }
210 };
211
212 bdy_index.assign(num_bc_vars(), BdyVars::NumTypes);
213 bdy_index[BCVars::Temp_bc_comp] = BdyVars::t;
214 bdy_index[BCVars::Salt_bc_comp] = BdyVars::s;
215 bdy_index[xvel_bc_idx] = BdyVars::u;
216 bdy_index[yvel_bc_idx] = BdyVars::v;
217 bdy_index[ubar_bc_idx] = BdyVars::ubar;
218 bdy_index[vbar_bc_idx] = BdyVars::vbar;
219 bdy_index[zeta_bc_idx] = BdyVars::zeta;
220
221 for (OrientationIter oit; oit; ++oit) {
222 Orientation ori = oit();
223 // These are simply defaults for Dirichlet faces -- they should be over-written below if needed
224 m_bc_extdir_vals[BCVars::Temp_bc_comp ][ori] = 1.e19_rt;
225 m_bc_extdir_vals[BCVars::Salt_bc_comp ][ori] = 1.e20_rt;
226 for (int icomp = Tracer_comp; icomp < ncons; ++icomp) {
227 m_bc_extdir_vals[icomp][ori] = 1.e21_rt + static_cast<Real>(icomp - Tracer_comp);
228 }
229
230 m_bc_extdir_vals[xvel_bc_idx][ori] = 0.0_rt; // default
231 m_bc_extdir_vals[yvel_bc_idx][ori] = 0.0_rt;
232 m_bc_extdir_vals[zvel_bc_idx][ori] = 0.0_rt;
233
234 m_bc_extdir_vals[ubar_bc_idx][ori] = 0.0_rt; // default
235 m_bc_extdir_vals[vbar_bc_idx][ori] = 0.0_rt;
236 m_bc_extdir_vals[u2d_simple_bc_idx][ori] = 0.0_rt;
237 m_bc_extdir_vals[v2d_simple_bc_idx][ori] = 0.0_rt;
238 }
239
240 // Whether to specify boundary conditions by variable (then side).
241 // Alternative is to do it by side by indicating keywords that indicate multiple variables
242 set_bcs_by_var = false;
243
244 ParmParse pp("remora");
245 pp.queryAdd("boundary_per_variable", set_bcs_by_var);
246 if (!set_bcs_by_var) {
247 f_by_side("xlo", Orientation(Direction::x,Orientation::low));
248 f_by_side("xhi", Orientation(Direction::x,Orientation::high));
249 f_by_side("ylo", Orientation(Direction::y,Orientation::low));
250 f_by_side("yhi", Orientation(Direction::y,Orientation::high));
251 } else {
252 f_by_var("temp", BCVars::Temp_bc_comp);
253 f_by_var("salt", BCVars::Salt_bc_comp);
254 for (int icomp = Tracer_comp; icomp < ncons; ++icomp) {
255 f_by_var("scalar", icomp);
256 }
257 f_by_var("u", xvel_bc_idx);
258 f_by_var("v", yvel_bc_idx);
259 f_by_var("w", zvel_bc_idx);
260 f_by_var("ubar", ubar_bc_idx);
261 f_by_var("vbar", vbar_bc_idx);
262 f_by_var("zeta", zeta_bc_idx);
263 f_by_var("tke", tke_bc_idx);
264 }
265
266 // Always specify z direction by side keyword
267 f_by_side("zlo", Orientation(Direction::z,Orientation::low));
268 f_by_side("zhi", Orientation(Direction::z,Orientation::high));
269
270 // *****************************************************************************
271 //
272 // Here we translate the physical boundary conditions -- one type per face --
273 // into logical boundary conditions for each velocity component
274 //
275 // *****************************************************************************
276 {
277 domain_bcs_type.resize(num_bc_vars());
278 domain_bcs_type_d.resize(num_bc_vars());
279
280 for (OrientationIter oit; oit; ++oit) {
281 Orientation ori = oit();
282 int dir = ori.coordDir();
283 Orientation::Side side = ori.faceDir();
284 // only do this for xvel and yvel
285 for (int i = 0; i < 2; i++) {
286 auto const bct = phys_bc_type[xvel_bc_idx+i][ori];
287 if ( bct == REMORA_BC::symmetry )
288 {
289 if (side == Orientation::low) {
290 domain_bcs_type[xvel_bc_idx+i].setLo(dir, REMORABCType::reflect_even);
291 if (i==1)
292 domain_bcs_type[xvel_bc_idx+dir].setLo(dir, REMORABCType::reflect_odd);
293 } else {
294 domain_bcs_type[xvel_bc_idx+i].setHi(dir, REMORABCType::reflect_even);
295 if (i==1)
296 domain_bcs_type[xvel_bc_idx+dir].setHi(dir, REMORABCType::reflect_odd);
297 }
298 }
299 else if (bct == REMORA_BC::outflow)
300 {
301 if (side == Orientation::low) {
302 domain_bcs_type[xvel_bc_idx+i].setLo(dir, REMORABCType::foextrap);
303 } else {
304 domain_bcs_type[xvel_bc_idx+i].setHi(dir, REMORABCType::foextrap);
305 }
306 }
307 else if (bct == REMORA_BC::inflow)
308 {
309 if (side == Orientation::low) {
310 domain_bcs_type[xvel_bc_idx+i].setLo(dir, REMORABCType::ext_dir);
311 } else {
312 domain_bcs_type[xvel_bc_idx+i].setHi(dir, REMORABCType::ext_dir);
313 }
314 }
315 else if (bct == REMORA_BC::no_slip_wall)
316 {
317 if (side == Orientation::low) {
318 domain_bcs_type[xvel_bc_idx+i].setLo(dir, REMORABCType::ext_dir);
319 } else {
320 domain_bcs_type[xvel_bc_idx+i].setHi(dir, REMORABCType::ext_dir);
321 }
322 }
323 else if (bct == REMORA_BC::slip_wall)
324 {
325 if (side == Orientation::low) {
326 domain_bcs_type[xvel_bc_idx+i].setLo(dir, REMORABCType::foextrap);
327 if (i==1) {
328 // Only normal direction has ext_dir
329 domain_bcs_type[xvel_bc_idx+dir].setLo(dir, REMORABCType::ext_dir);
330 }
331
332 } else {
333 domain_bcs_type[xvel_bc_idx+i].setHi(dir, REMORABCType::foextrap);
334 if (i==1) {
335 // Only normal direction has ext_dir
336 domain_bcs_type[xvel_bc_idx+dir].setHi(dir, REMORABCType::ext_dir);
337 }
338 }
339 }
340 else if (bct == REMORA_BC::periodic)
341 {
342 if (side == Orientation::low) {
343 domain_bcs_type[xvel_bc_idx+i].setLo(dir, REMORABCType::int_dir);
344 } else {
345 domain_bcs_type[xvel_bc_idx+i].setHi(dir, REMORABCType::int_dir);
346 }
347 }
348 else if (bct == REMORA_BC::clamped)
349 {
350 if (side == Orientation::low) {
351 domain_bcs_type[xvel_bc_idx+i].setLo(dir, REMORABCType::clamped);
352 } else {
353 domain_bcs_type[xvel_bc_idx+i].setHi(dir, REMORABCType::clamped);
354 }
355 }
356 else if (bct == REMORA_BC::orlanski_rad)
357 {
358 if (side == Orientation::low) {
359 domain_bcs_type[xvel_bc_idx+i].setLo(dir, REMORABCType::orlanski_rad);
360 } else {
361 domain_bcs_type[xvel_bc_idx+i].setHi(dir, REMORABCType::orlanski_rad);
362 }
363 }
364 else if (bct == REMORA_BC::orlanski_rad_nudge)
365 {
366 if (side == Orientation::low) {
367 domain_bcs_type[xvel_bc_idx+i].setLo(dir, REMORABCType::orlanski_rad_nudge);
368 } else {
369 domain_bcs_type[xvel_bc_idx+i].setHi(dir, REMORABCType::orlanski_rad_nudge);
370 }
371 }
372 else
373 {
374 amrex::Abort("Velocity boundary condition not validly specified");
375 }
376 }
377
378 // Always set zvel_bc to foextrap
379 if (side == Orientation::low) {
380 domain_bcs_type[zvel_bc_idx].setLo(dir, REMORABCType::foextrap);
381 } else {
382 domain_bcs_type[zvel_bc_idx].setHi(dir, REMORABCType::foextrap);
383 }
384 }
385 }
386
387 // *****************************************************************************
388 //
389 // Here we translate the physical boundary conditions -- one type per face --
390 // into logical boundary conditions for each cell-centered variable
391 //
392 // *****************************************************************************
393 {
394 for (OrientationIter oit; oit; ++oit) {
395 Orientation ori = oit();
396 int dir = ori.coordDir();
397 Orientation::Side side = ori.faceDir();
398 for (int i = 0; i < ncons; i++) {
399 auto const bct = phys_bc_type[BCVars::cons_bc+i][ori];
400 if ( bct == REMORA_BC::symmetry )
401 {
402 if (side == Orientation::low) {
403 domain_bcs_type[BCVars::cons_bc+i].setLo(dir, REMORABCType::reflect_even);
404 } else {
405 domain_bcs_type[BCVars::cons_bc+i].setHi(dir, REMORABCType::reflect_even);
406 }
407 }
408 else if ( bct == REMORA_BC::outflow)
409 {
410 if (side == Orientation::low) {
411 domain_bcs_type[BCVars::cons_bc+i].setLo(dir, REMORABCType::foextrap);
412 } else {
413 domain_bcs_type[BCVars::cons_bc+i].setHi(dir, REMORABCType::foextrap);
414 }
415 }
416 else if ( bct == REMORA_BC::no_slip_wall)
417 {
418 if (side == Orientation::low) {
419 domain_bcs_type[BCVars::cons_bc+i].setLo(dir, REMORABCType::foextrap);
420 } else {
421 domain_bcs_type[BCVars::cons_bc+i].setHi(dir, REMORABCType::foextrap);
422 }
423 }
424 else if (bct == REMORA_BC::slip_wall)
425 {
426 if (side == Orientation::low) {
427 domain_bcs_type[BCVars::cons_bc+i].setLo(dir, REMORABCType::foextrap);
428 } else {
429 domain_bcs_type[BCVars::cons_bc+i].setHi(dir, REMORABCType::foextrap);
430 }
431 }
432 else if (bct == REMORA_BC::inflow)
433 {
434 if (side == Orientation::low) {
435 domain_bcs_type[BCVars::cons_bc+i].setLo(dir, REMORABCType::ext_dir);
436 } else {
437 domain_bcs_type[BCVars::cons_bc+i].setHi(dir, REMORABCType::ext_dir);
438 }
439 }
440 else if (bct == REMORA_BC::periodic)
441 {
442 if (side == Orientation::low) {
443 domain_bcs_type[BCVars::cons_bc+i].setLo(dir, REMORABCType::int_dir);
444 } else {
445 domain_bcs_type[BCVars::cons_bc+i].setHi(dir, REMORABCType::int_dir);
446 }
447 }
448 else if ( bct == REMORA_BC::clamped)
449 {
450 if (side == Orientation::low) {
451 domain_bcs_type[BCVars::cons_bc+i].setLo(dir, REMORABCType::clamped);
452 } else {
453 domain_bcs_type[BCVars::cons_bc+i].setHi(dir, REMORABCType::clamped);
454 }
455 }
456 else if ( bct == REMORA_BC::orlanski_rad)
457 {
458 if (side == Orientation::low) {
459 domain_bcs_type[BCVars::cons_bc+i].setLo(dir, REMORABCType::orlanski_rad);
460 } else {
461 domain_bcs_type[BCVars::cons_bc+i].setHi(dir, REMORABCType::orlanski_rad);
462 }
463 }
464 else if ( bct == REMORA_BC::orlanski_rad_nudge)
465 {
466 if (side == Orientation::low) {
467 domain_bcs_type[BCVars::cons_bc+i].setLo(dir, REMORABCType::orlanski_rad_nudge);
468 } else {
469 domain_bcs_type[BCVars::cons_bc+i].setHi(dir, REMORABCType::orlanski_rad_nudge);
470 }
471 }
472 else
473 {
474 amrex::Abort("Scalar/tracer boundary condition not validly specified");
475 }
476 }
477 }
478 }
479
480 // *****************************************************************************
481 //
482 // Here we translate the physical boundary conditions -- one type per face --
483 // into logical boundary conditions for ubar and vbar. Also add simplified
484 // 2d boundary condition (corresponding to BCs in bc_2d.F
485 //
486 // *****************************************************************************
487 {
488 for (OrientationIter oit; oit; ++oit) {
489 Orientation ori = oit();
490 int dir = ori.coordDir();
491 Orientation::Side side = ori.faceDir();
492 for (int i = 0; i < 2; i++) {
493 auto const bct = phys_bc_type[ubar_bc_idx+i][ori];
494 if ( bct == REMORA_BC::symmetry )
495 {
496 if (side == Orientation::low) {
497 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::reflect_even);
498 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::reflect_even);
499 if (i==1 and dir!=2) {
500 domain_bcs_type[ubar_bc_idx+dir].setLo(dir, REMORABCType::reflect_odd);
501 domain_bcs_type[u2d_simple_bc_idx+dir].setLo(dir, REMORABCType::reflect_odd);
502 }
503 } else {
504 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::reflect_even);
505 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::reflect_even);
506 if (i==1 and dir!=2) {
507 domain_bcs_type[ubar_bc_idx+dir].setHi(dir, REMORABCType::reflect_odd);
508 domain_bcs_type[u2d_simple_bc_idx+dir].setHi(dir, REMORABCType::reflect_odd);
509 }
510 }
511 }
512 else if (bct == REMORA_BC::outflow)
513 {
514 if (side == Orientation::low) {
515 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::foextrap);
516 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::foextrap);
517 } else {
518 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::foextrap);
519 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::foextrap);
520 }
521 }
522 else if (bct == REMORA_BC::inflow)
523 {
524 if (side == Orientation::low) {
525 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::ext_dir);
526 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::ext_dir);
527 } else {
528 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::ext_dir);
529 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::ext_dir);
530 }
531 }
532 else if (bct == REMORA_BC::no_slip_wall)
533 {
534 if (side == Orientation::low) {
535 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::ext_dir);
536 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::ext_dir);
537 } else {
538 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::ext_dir);
539 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::ext_dir);
540 }
541 }
542 else if (bct == REMORA_BC::slip_wall)
543 {
544 if (side == Orientation::low) {
545 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::foextrap);
546 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::foextrap);
547 if (i==1 and dir!=2) {
548 // Only normal direction has ext_dir
549 domain_bcs_type[ubar_bc_idx+dir].setLo(dir, REMORABCType::ext_dir);
550 domain_bcs_type[u2d_simple_bc_idx+dir].setLo(dir, REMORABCType::ext_dir);
551 }
552
553 } else {
554 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::foextrap);
555 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::foextrap);
556 if (i==1 and dir!=2) {
557 // Only normal direction has ext_dir
558 domain_bcs_type[ubar_bc_idx+dir].setHi(dir, REMORABCType::ext_dir);
559 domain_bcs_type[u2d_simple_bc_idx+dir].setHi(dir, REMORABCType::ext_dir);
560 }
561 }
562 }
563 else if (bct == REMORA_BC::periodic)
564 {
565 if (side == Orientation::low) {
566 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::int_dir);
567 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::int_dir);
568 } else {
569 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::int_dir);
570 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::int_dir);
571 }
572 }
573 else if (bct == REMORA_BC::clamped)
574 {
575 if (side == Orientation::low) {
576 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::clamped);
577 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::foextrap);
578 } else {
579 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::clamped);
580 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::foextrap);
581 }
582 }
583 else if (bct == REMORA_BC::flather)
584 {
585 if (side == Orientation::low) {
586 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::chapman);
587 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::foextrap);
588 if (i==1 and dir!=2) {
589 // Only normal direction has Flather
590 domain_bcs_type[ubar_bc_idx+dir].setLo(dir, REMORABCType::flather);
591 domain_bcs_type[u2d_simple_bc_idx+dir].setLo(dir, REMORABCType::foextrap);
592 }
593
594 } else {
595 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::chapman);
596 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::foextrap);
597 if (i==1 and dir!=2) {
598 // Only normal direction has Flather
599 domain_bcs_type[ubar_bc_idx+dir].setHi(dir, REMORABCType::flather);
600 domain_bcs_type[u2d_simple_bc_idx+dir].setHi(dir, REMORABCType::foextrap);
601 }
602 }
603 }
604 else if (bct == REMORA_BC::orlanski_rad)
605 {
606 if (side == Orientation::low) {
607 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::orlanski_rad);
608 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::foextrap);
609 } else {
610 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::orlanski_rad);
611 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::foextrap);
612 }
613 }
614 else if (bct == REMORA_BC::orlanski_rad_nudge)
615 {
616 if (side == Orientation::low) {
617 domain_bcs_type[ubar_bc_idx+i].setLo(dir, REMORABCType::orlanski_rad_nudge);
618 domain_bcs_type[u2d_simple_bc_idx+i].setLo(dir, REMORABCType::foextrap);
619 } else {
620 domain_bcs_type[ubar_bc_idx+i].setHi(dir, REMORABCType::orlanski_rad_nudge);
621 domain_bcs_type[u2d_simple_bc_idx+i].setHi(dir, REMORABCType::foextrap);
622 }
623 }
624 else
625 {
626 amrex::Abort("ubar or vbar boundary condition not validly specified");
627 }
628 }
629 }
630 }
631
632 // *****************************************************************************
633 //
634 // Here we translate the physical boundary conditions -- one type per face --
635 // into logical boundary conditions for zeta and tke
636 //
637 // *****************************************************************************
638 {
639 for (OrientationIter oit; oit; ++oit) {
640 Orientation ori = oit();
641 int dir = ori.coordDir();
642 Orientation::Side side = ori.faceDir();
643 for (int i = 0; i < 2; i++) {
644 auto const bct = phys_bc_type[zeta_bc_idx+i][ori];
645 if ( bct == REMORA_BC::symmetry )
646 {
647 if (side == Orientation::low) {
648 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::reflect_even);
649 } else {
650 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::reflect_even);
651 }
652 }
653 else if ( bct == REMORA_BC::outflow)
654 {
655 if (side == Orientation::low) {
656 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::foextrap);
657 } else {
658 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::foextrap);
659 }
660 }
661 else if ( bct == REMORA_BC::no_slip_wall)
662 {
663 if (side == Orientation::low) {
664 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::foextrap);
665 } else {
666 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::foextrap);
667 }
668 }
669 else if (bct == REMORA_BC::slip_wall)
670 {
671 if (side == Orientation::low) {
672 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::foextrap);
673 } else {
674 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::foextrap);
675 }
676 }
677 else if (bct == REMORA_BC::inflow)
678 {
679 if (side == Orientation::low) {
680 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::ext_dir);
681 } else {
682 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::ext_dir);
683 }
684 }
685 else if (bct == REMORA_BC::periodic)
686 {
687 if (side == Orientation::low) {
688 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::int_dir);
689 } else {
690 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::int_dir);
691 }
692 }
693 else if (bct == REMORA_BC::chapman)
694 {
695 if (side == Orientation::low) {
696 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::chapman);
697 } else {
698 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::chapman);
699 }
700 }
701 else if ( bct == REMORA_BC::clamped)
702 {
703 if (side == Orientation::low) {
704 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::clamped);
705 } else {
706 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::clamped);
707 }
708 }
709 else if ( bct == REMORA_BC::orlanski_rad)
710 {
711 if (side == Orientation::low) {
712 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::orlanski_rad);
713 } else {
714 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::orlanski_rad);
715 }
716 }
717 else if ( bct == REMORA_BC::orlanski_rad_nudge)
718 {
719 if (side == Orientation::low) {
720 domain_bcs_type[zeta_bc_idx+i].setLo(dir, REMORABCType::orlanski_rad_nudge);
721 } else {
722 domain_bcs_type[zeta_bc_idx+i].setHi(dir, REMORABCType::orlanski_rad_nudge);
723 }
724 }
725 else
726 {
727 amrex::Abort("Free surface (zeta) boundary condition not validly specified");
728 }
729 }
730 }
731 }
732
733 // *****************************************************************************
734 //
735 // Here we define a boundary condition that will foextrap while respecting periodicity
736 // This is used as a "null BC"
737 //
738 // *****************************************************************************
739 {
740 for (OrientationIter oit; oit; ++oit) {
741 Orientation ori = oit();
742 int dir = ori.coordDir();
743 Orientation::Side side = ori.faceDir();
744 if (side == Orientation::low) {
745 domain_bcs_type[foextrap_periodic_bc_idx].setLo(dir, REMORABCType::foextrap);
746 } else {
747 domain_bcs_type[foextrap_periodic_bc_idx].setHi(dir, REMORABCType::foextrap);
748 }
749 }
750 }
751
752 // *****************************************************************************
753 //
754 // Here we define a boundary condition that will unconditionally foextrap
755 //
756 // *****************************************************************************
757 {
758 for (OrientationIter oit; oit; ++oit) {
759 Orientation ori = oit();
760 int dir = ori.coordDir();
761 Orientation::Side side = ori.faceDir();
762 if (side == Orientation::low) {
763 domain_bcs_type[foextrap_bc_idx].setLo(dir, REMORABCType::foextrap);
764 } else {
765 domain_bcs_type[foextrap_bc_idx].setHi(dir, REMORABCType::foextrap);
766 }
767 }
768 }
769
770
771#ifdef AMREX_USE_GPU
772 Gpu::htod_memcpy
773 (domain_bcs_type_d.data(), domain_bcs_type.data(),
774 sizeof(amrex::BCRec)*num_bc_vars());
775#else
776 std::memcpy
777 (domain_bcs_type_d.data(), domain_bcs_type.data(),
778 sizeof(amrex::BCRec)*num_bc_vars());
779#endif
780}
REMORA_BC::slip_wall
@ slip_wall
REMORA_BC::orlanski_rad
@ orlanski_rad
REMORA_BC::flather
@ flather
REMORA_BC::symmetry
@ symmetry
REMORA_BC::clamped
@ clamped
REMORA_BC::undefined
@ undefined
REMORA_BC::chapman
@ chapman
REMORA_BC::outflow
@ outflow
REMORA_BC::no_slip_wall
@ no_slip_wall
REMORA_BC::inflow
@ inflow
REMORA_BC::periodic
@ periodic
REMORA_BC::orlanski_rad_nudge
@ orlanski_rad_nudge
Tracer_comp
#define Tracer_comp
Definition REMORA_IndexDefines.H:10
REMORA::foextrap_periodic_bc
int foextrap_periodic_bc() const noexcept
Definition REMORA.H:1059
REMORA::ncons
int ncons
Number of conserved scalars in the state (temperature + salt + passive scalars)
Definition REMORA.H:1346
REMORA::zvel_bc
int zvel_bc() const noexcept
Definition REMORA.H:1054
REMORA::init_bcs
void init_bcs()
Read in boundary parameters from input file and set up data structures.
Definition REMORA_init_bcs.cpp:9
REMORA::xvel_bc
int xvel_bc() const noexcept
Definition REMORA.H:1052
REMORA::domain_bcs_type
amrex::Vector< amrex::BCRec > domain_bcs_type
vector (over BCVars) of BCRecs
Definition REMORA.H:1295
REMORA::set_bcs_by_var
bool set_bcs_by_var
whether to set boundary conditions by variable rather than just by side
Definition REMORA.H:1280
REMORA::tke_bc
int tke_bc() const noexcept
Definition REMORA.H:1058
REMORA::num_bc_vars
int num_bc_vars() const noexcept
Definition REMORA.H:1063
REMORA::yvel_bc
int yvel_bc() const noexcept
Definition REMORA.H:1053
REMORA::foextrap_bc
int foextrap_bc() const noexcept
Definition REMORA.H:1060
REMORA::v2d_simple_bc
int v2d_simple_bc() const noexcept
Definition REMORA.H:1062
REMORA::zeta_bc
int zeta_bc() const noexcept
Definition REMORA.H:1057
REMORA::phys_bc_type
amrex::Vector< amrex::GpuArray< REMORA_BC, AMREX_SPACEDIM *2 > > phys_bc_type
Array holding the "physical" boundary condition types (e.g. "inflow")
Definition REMORA.H:1306
REMORA::bdy_index
amrex::Vector< int > bdy_index
Container to connect boundary data being read in boundary condition containers.
Definition REMORA.H:1315
REMORA::u2d_simple_bc
int u2d_simple_bc() const noexcept
Definition REMORA.H:1061
REMORA::solverChoice
static SolverChoice solverChoice
Container for algorithmic choices.
Definition REMORA.H:1403
REMORA::ubar_bc
int ubar_bc() const noexcept
Definition REMORA.H:1055
REMORA::vbar_bc
int vbar_bc() const noexcept
Definition REMORA.H:1056
REMORA::domain_bc_type
amrex::Array< std::string, 2 *AMREX_SPACEDIM > domain_bc_type
Array of strings describing domain boundary conditions.
Definition REMORA.H:1300
REMORA::phys_bc_need_data
amrex::GpuArray< amrex::GpuArray< bool, AMREX_SPACEDIM *2 >, BdyVars::NumTypes+1 > phys_bc_need_data
These are flags that indicate whether we need to read in boundary data from file.
Definition REMORA.H:1309
REMORA::m_bc_extdir_vals
amrex::Vector< amrex::GpuArray< amrex::Real, AMREX_SPACEDIM *2 > > m_bc_extdir_vals
Array holding the Dirichlet values at walls which need them.
Definition REMORA.H:1303
REMORA::domain_bcs_type_d
amrex::Gpu::DeviceVector< amrex::BCRec > domain_bcs_type_d
GPU vector (over BCVars) of BCRecs.
Definition REMORA.H:1297
BCVars::cons_bc
static constexpr int cons_bc
Definition REMORA_IndexDefines.H:24
BCVars::Salt_bc_comp
static constexpr int Salt_bc_comp
Definition REMORA_IndexDefines.H:26
BCVars::Temp_bc_comp
static constexpr int Temp_bc_comp
Definition REMORA_IndexDefines.H:25