190 BL_PROFILE_REGION(
"realbdy_compute_interior_ghost_RHS()");
198 if (width > set_width+1) width -= 1;
201 Real F1 = 1./(10.*delta_t);
202 Real F2 = 1./(50.*delta_t);
205 Real dT = bdy_time_interval;
206 Real time_since_start = time - start_bdy_time;
207 int n_time =
static_cast<int>( time_since_start / dT);
208 Real alpha = (time_since_start - n_time * dT) / dT;
209 AMREX_ALWAYS_ASSERT( alpha >= 0. && alpha <= 1.0);
210 Real oma = 1.0 - alpha;
218 FArrayBox U_xlo, U_xhi, U_ylo, U_yhi;
219 FArrayBox V_xlo, V_xhi, V_ylo, V_yhi;
220 FArrayBox T_xlo, T_xhi, T_ylo, T_yhi;
241 for (
int ivar(ivarU); ivar < BdyEnd; ivar++) {
242 int ivar_idx = var_map[ivar];
243 Box domain = geom.Domain();
244 auto ixtype = S_cur_data[ivar_idx].boxArray().ixType();
245 domain.convert(ixtype);
250 IntVect ng_vect{2,2,0};
251 Box gdom(domain); gdom.grow(ng_vect);
252 Box bx_xlo, bx_xhi, bx_ylo, bx_yhi;
260 U_xlo.resize(bx_xlo,1,The_Async_Arena()); U_xhi.resize(bx_xhi,1,The_Async_Arena());
261 U_ylo.resize(bx_ylo,1,The_Async_Arena()); U_yhi.resize(bx_yhi,1,The_Async_Arena());
262 }
else if (ivar == ivarV) {
263 V_xlo.resize(bx_xlo,1,The_Async_Arena()); V_xhi.resize(bx_xhi,1,The_Async_Arena());
264 V_ylo.resize(bx_ylo,1,The_Async_Arena()); V_yhi.resize(bx_yhi,1,The_Async_Arena());
265 }
else if (ivar == ivarT){
266 T_xlo.resize(bx_xlo,1,The_Async_Arena()); T_xhi.resize(bx_xhi,1,The_Async_Arena());
267 T_ylo.resize(bx_ylo,1,The_Async_Arena()); T_yhi.resize(bx_yhi,1,The_Async_Arena());
282 for (
int ivar(ivarU); ivar < BdyEnd; ivar++) {
283 int ivar_idx = var_map[ivar];
284 Box domain = geom.Domain();
285 auto ixtype = S_cur_data[ivar_idx].boxArray().ixType();
286 domain.convert(ixtype);
287 const auto& dom_lo = lbound(domain);
288 const auto& dom_hi = ubound(domain);
291 #pragma omp parallel if (Gpu::notInLaunchRegion())
293 for (MFIter mfi(S_cur_data[ivar_idx],TilingIfNotGPU()); mfi.isValid(); ++mfi) {
296 IntVect ng_vect{2,2,0};
297 Box gtbx = grow(mfi.tilebox(ixtype.toIntVect()),ng_vect);
298 Box tbx_xlo, tbx_xhi, tbx_ylo, tbx_yhi;
304 Array4<Real> arr_xlo; Array4<Real> arr_xhi;
305 Array4<Real> arr_ylo; Array4<Real> arr_yhi;
307 arr_xlo = U_xlo.array(); arr_xhi = U_xhi.array();
308 arr_ylo = U_ylo.array(); arr_yhi = U_yhi.array();
309 }
else if (ivar == ivarV) {
310 arr_xlo = V_xlo.array(); arr_xhi = V_xhi.array();
311 arr_ylo = V_ylo.array(); arr_yhi = V_yhi.array();
312 }
else if (ivar == ivarT){
313 arr_xlo = T_xlo.array(); arr_xhi = T_xhi.array();
314 arr_ylo = T_ylo.array(); arr_yhi = T_yhi.array();
320 const auto& bdatxlo_n = bdy_data_xlo[n_time ][ivar].const_array();
321 const auto& bdatxlo_np1 = bdy_data_xlo[n_time+1][ivar].const_array();
322 const auto& bdatxhi_n = bdy_data_xhi[n_time ][ivar].const_array();
323 const auto& bdatxhi_np1 = bdy_data_xhi[n_time+1][ivar].const_array();
324 const auto& bdatylo_n = bdy_data_ylo[n_time ][ivar].const_array();
325 const auto& bdatylo_np1 = bdy_data_ylo[n_time+1][ivar].const_array();
326 const auto& bdatyhi_n = bdy_data_yhi[n_time ][ivar].const_array();
327 const auto& bdatyhi_np1 = bdy_data_yhi[n_time+1][ivar].const_array();
336 int offset = set_width - 1;
337 if (width > set_width)
offset = width;
340 ParallelFor(tbx_xlo, tbx_xhi,
341 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
343 int ii = std::max(i , dom_lo.x);
344 ii = std::min(ii, dom_lo.x+
offset);
345 int jj = std::max(j , dom_lo.y);
346 jj = std::min(jj, dom_hi.y);
350 rho_interp = 0.5 * ( r_arr(i-1,j ,k) + r_arr(i,j,k) );
351 }
else if (ivar==ivarV) {
352 rho_interp = 0.5 * ( r_arr(i ,j-1,k) + r_arr(i,j,k) );
354 rho_interp = r_arr(i,j,k);
356 arr_xlo(i,j,k) = rho_interp * ( oma * bdatxlo_n (ii,jj,k,0)
357 + alpha * bdatxlo_np1(ii,jj,k,0) );
359 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
361 int ii = std::max(i , dom_hi.x-
offset);
362 ii = std::min(ii, dom_hi.x);
363 int jj = std::max(j , dom_lo.y);
364 jj = std::min(jj, dom_hi.y);
368 rho_interp = 0.5 * ( r_arr(i-1,j ,k) + r_arr(i,j,k) );
369 }
else if (ivar==ivarV) {
370 rho_interp = 0.5 * ( r_arr(i ,j-1,k) + r_arr(i,j,k) );
372 rho_interp = r_arr(i,j,k);
374 arr_xhi(i,j,k) = rho_interp * ( oma * bdatxhi_n (ii,jj,k,0)
375 + alpha * bdatxhi_np1(ii,jj,k,0) );
378 ParallelFor(tbx_ylo, tbx_yhi,
379 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
381 int ii = std::max(i , dom_lo.x);
382 ii = std::min(ii, dom_hi.x);
383 int jj = std::max(j , dom_lo.y);
384 jj = std::min(jj, dom_lo.y+
offset);
388 rho_interp = 0.5 * ( r_arr(i-1,j ,k) + r_arr(i,j,k) );
389 }
else if (ivar==ivarV) {
390 rho_interp = 0.5 * ( r_arr(i ,j-1,k) + r_arr(i,j,k) );
392 rho_interp = r_arr(i,j,k);
394 arr_ylo(i,j,k) = rho_interp * ( oma * bdatylo_n (ii,jj,k,0)
395 + alpha * bdatylo_np1(ii,jj,k,0) );
397 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
399 int ii = std::max(i , dom_lo.x);
400 ii = std::min(ii, dom_hi.x);
401 int jj = std::max(j , dom_hi.y-
offset);
402 jj = std::min(jj, dom_hi.y);
406 rho_interp = 0.5 * ( r_arr(i-1,j ,k) + r_arr(i,j,k) );
407 }
else if (ivar==ivarV) {
408 rho_interp = 0.5 * ( r_arr(i ,j-1,k) + r_arr(i,j,k) );
410 rho_interp = r_arr(i,j,k);
412 arr_yhi(i,j,k) = rho_interp * ( oma * bdatyhi_n (ii,jj,k,0)
413 + alpha * bdatyhi_np1(ii,jj,k,0) );
425 for (
int ivar(ivarU); ivar < BdyEnd; ivar++) {
426 int ivar_idx = ivar_map[ivar];
427 int icomp = comp_map[ivar];
429 Box domain = geom.Domain();
430 auto ix_type = S_old_data[ivar_idx].boxArray().ixType();
431 auto iv_type = ix_type.toIntVect();
432 domain.convert(ix_type);
433 const auto& dom_hi = ubound(domain);
434 const auto& dom_lo = lbound(domain);
436 int set_width_x = (iv_type[0]) ? set_width : set_width-1;
437 int set_width_y = (iv_type[1]) ? set_width : set_width-1;
440 #pragma omp parallel if (Gpu::notInLaunchRegion())
442 for (MFIter mfi(S_old_data[ivar_idx],TilingIfNotGPU()); mfi.isValid(); ++mfi) {
443 Box tbx = mfi.tilebox();
444 Box tbx_xlo, tbx_xhi, tbx_ylo, tbx_yhi;
449 Array4<Real> rhs_arr; Array4<Real> data_arr;
450 Array4<Real> arr_xlo; Array4<Real> arr_xhi;
451 Array4<Real> arr_ylo; Array4<Real> arr_yhi;
453 arr_xlo = U_xlo.array(); arr_xhi = U_xhi.array();
454 arr_ylo = U_ylo.array(); arr_yhi = U_yhi.array();
457 }
else if (ivar == ivarV) {
458 arr_xlo = V_xlo.array(); arr_xhi = V_xhi.array();
459 arr_ylo = V_ylo.array(); arr_yhi = V_yhi.array();
462 }
else if (ivar == ivarT){
463 arr_xlo = T_xlo.array(); arr_xhi = T_xhi.array();
464 arr_ylo = T_ylo.array(); arr_yhi = T_yhi.array();
472 width, set_width_x, set_width_y,
474 tbx_xlo, tbx_xhi, tbx_ylo, tbx_yhi,
475 arr_xlo, arr_xhi, arr_ylo, arr_yhi,
487 if (width > set_width) {
488 for (
int ivar(ivarU); ivar < BdyEnd; ivar++) {
489 int ivar_idx = ivar_map[ivar];
490 int icomp = comp_map[ivar];
492 Box domain = geom.Domain();
493 domain.convert(S_cur_data[ivar_idx].boxArray().ixType());
494 const auto& dom_hi = ubound(domain);
495 const auto& dom_lo = lbound(domain);
496 IntVect ng_vect = S_cur_data[ivar_idx].nGrowVect();
499 #pragma omp parallel if (Gpu::notInLaunchRegion())
501 for (MFIter mfi(S_cur_data[ivar_idx],TilingIfNotGPU()); mfi.isValid(); ++mfi) {
502 Box tbx = mfi.tilebox();
503 Box tbx_xlo, tbx_xhi, tbx_ylo, tbx_yhi;
509 Array4<Real> rhs_arr; Array4<Real> data_arr;
510 Array4<Real> arr_xlo; Array4<Real> arr_xhi;
511 Array4<Real> arr_ylo; Array4<Real> arr_yhi;
513 arr_xlo = U_xlo.array(); arr_xhi = U_xhi.array();
514 arr_ylo = U_ylo.array(); arr_yhi = U_yhi.array();
517 }
else if (ivar == ivarV) {
518 arr_xlo = V_xlo.array(); arr_xhi = V_xhi.array();
519 arr_ylo = V_ylo.array(); arr_yhi = V_yhi.array();
522 }
else if (ivar == ivarT){
523 arr_xlo = T_xlo.array(); arr_xhi = T_xhi.array();
524 arr_ylo = T_ylo.array(); arr_yhi = T_yhi.array();
532 width, set_width, dom_lo, dom_hi, F1, F2,
533 tbx_xlo, tbx_xhi, tbx_ylo, tbx_yhi,
534 arr_xlo, arr_xhi, arr_ylo, arr_yhi,
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:37
void realbdy_interior_bxs_xy(const Box &bx, const Box &domain, const int &width, Box &bx_xlo, Box &bx_xhi, Box &bx_ylo, Box &bx_yhi, const int &set_width, const IntVect &ng_vect, const bool get_int_ng)
Definition: ERF_InteriorGhostCells.cpp:23
AMREX_FORCE_INLINE IntVect offset(const int face_dir, const int normal)
Definition: ERF_ReadBndryPlanes.cpp:28
AMREX_GPU_HOST AMREX_FORCE_INLINE void realbdy_set_rhs_in_spec_region(const amrex::Real &dt, const int &icomp, const int &num_var, const int &width, const int &set_width_x, const int &set_width_y, const amrex::Dim3 &dom_lo, const amrex::Dim3 &dom_hi, const amrex::Box &bx_xlo, const amrex::Box &bx_xhi, const amrex::Box &bx_ylo, const amrex::Box &bx_yhi, const amrex::Array4< const amrex::Real > &arr_xlo, const amrex::Array4< const amrex::Real > &arr_xhi, const amrex::Array4< const amrex::Real > &arr_ylo, const amrex::Array4< const amrex::Real > &arr_yhi, const amrex::Array4< const amrex::Real > &data_arr, const amrex::Array4< amrex::Real > &rhs_arr)
Definition: ERF_Utils.H:169
AMREX_GPU_HOST AMREX_FORCE_INLINE void realbdy_compute_laplacian_relaxation(const int &icomp, const int &num_var, const int &width, const int &set_width, const amrex::Dim3 &dom_lo, const amrex::Dim3 &dom_hi, const amrex::Real &F1, const amrex::Real &F2, const amrex::Box &bx_xlo, const amrex::Box &bx_xhi, const amrex::Box &bx_ylo, const amrex::Box &bx_yhi, const amrex::Array4< const amrex::Real > &arr_xlo, const amrex::Array4< const amrex::Real > &arr_xhi, const amrex::Array4< const amrex::Real > &arr_ylo, const amrex::Array4< const amrex::Real > &arr_yhi, const amrex::Array4< const amrex::Real > &data_arr, const amrex::Array4< amrex::Real > &rhs_arr)
Definition: ERF_Utils.H:263
@ U
Definition: ERF_IndexDefines.H:108
@ NumTypes
Definition: ERF_IndexDefines.H:112
@ T
Definition: ERF_IndexDefines.H:110
@ V
Definition: ERF_IndexDefines.H:109
@ xvel
Definition: ERF_IndexDefines.H:141
@ cons
Definition: ERF_IndexDefines.H:140
@ yvel
Definition: ERF_IndexDefines.H:142