WarpX
WarpXPushFieldsEM_K.H
Go to the documentation of this file.
1 /* Copyright 2019-2020
2  *
3  * This file is part of WarpX.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 
8 #ifndef WARPX_PushFieldsEM_K_h
9 #define WARPX_PushFieldsEM_K_h
10 
11 #include "Utils/WarpXConst.H"
12 
13 #include <AMReX.H>
14 
15 /*
16  * \brief Return a tilebox that only covers the outer half of the guard cells.
17  * For tileboxes that don't include cells beyond the whole domain,
18  * an empty box is returned.
19  *
20  * \param[in,out] input_tb tilebox to be modified
21  * \param[in] dir direction where the tilebox smallEnd/bigEnd is modified
22  * \param[in] n_domain number of valid cells in the whole simulation domain
23  * \param[in] tb_smallEnd the lowest index of the tilebox, including guard cells
24  * \param[in] tb_bigEnd the highest index of the tilebox, including guard cells
25  */
28  const amrex::Box& input_tb,
29  const int dir,
30  const int iside,
31  const int n_domain,
32  const int tb_smallEnd,
33  const int tb_bigEnd)
34 {
35  using namespace amrex;
36 
37  amrex::Box constrained_tb = input_tb;
38 
39  // If the input_tb does not overlap either the lower or upper guard,
40  // an empty box is returned.
41 
42  if (iside == 0)
43  {
44  // Lower guard
45  const int n_guard = -tb_smallEnd;
46  const int upper_bound = (tb_smallEnd < 0 ? -n_guard/2 : tb_smallEnd);
47  constrained_tb.setBig(dir, upper_bound - 1);
48  }
49  else if (iside == 1)
50  {
51  // Upper guard
52  const int n_guard = tb_bigEnd - n_domain;
53  const int lower_bound = (tb_bigEnd > n_domain ? n_domain + n_guard/2 : tb_bigEnd);
54  constrained_tb.setSmall(dir, lower_bound + 1);
55  }
56 
57  return constrained_tb;
58 }
59 
60 /*
61  * \brief Damp a given field in the guard cells along a given direction
62  *
63  * \param[in,out] mf_arr array that contains the field values to be damped
64  * \oaram[in] i index along x
65  * \oaram[in] j index along y (in 3D) or z (in 2D/RZ)
66  * \oaram[in] k index along z (in 3D, \c k = 0 in 2D/RZ)
67  * \param[in] icomp index along the fourth component of the array
68  * \param]in] dir direction where the field will be damped
69  * \param[in] n_domain number of valid cells in the whole simulation domain
70  * \param[in] tb_smallEnd the lowest index of the tilebox, including guard cells
71  * \param[in] tb_bigEnd the highest index of the tilebox, including guard cells
72  */
75  amrex::Array4<amrex::Real> const& mf_arr,
76  const int i,
77  const int j,
78  const int k,
79  const int icomp,
80  const int dir,
81  const int n_domain,
82  const int tb_smallEnd,
83  const int tb_bigEnd)
84 {
85  using namespace amrex;
86 
87  // dir = 0: idx = i (x)
88  // dir = 1: idx = j (y in 3D, z in 2D/RZ)
89  // dir = 2: idx = k (z in 3D)
90  const int idx = ((dir == 0) ? i : ((dir == 1) ? j : k));
91 
92  if (idx < 0)
93  {
94  // Apply damping factor in guards cells below the lower end of the domain
95  const int n_guard = -tb_smallEnd;
96 
97  const auto cell = static_cast<amrex::Real>(idx + n_guard);
98 
99  const amrex::Real phase = MathConst::pi * cell / n_guard;
100  const amrex::Real sin_phase = std::sin(phase);
101  const amrex::Real damp_factor = sin_phase * sin_phase;
102 
103  mf_arr(i,j,k,icomp) *= damp_factor;
104  }
105  else if (idx > n_domain)
106  {
107  // Apply damping factor in guards cells above the upper end of the domain
108  const int n_guard = tb_bigEnd - n_domain;
109 
110  const auto cell = static_cast<amrex::Real>(tb_bigEnd - idx);
111 
112  const amrex::Real phase = MathConst::pi * cell / n_guard;
113  const amrex::Real sin_phase = std::sin(phase);
114  const amrex::Real damp_factor = sin_phase * sin_phase;
115 
116  mf_arr(i,j,k,icomp) *= damp_factor;
117  }
118 }
119 
120 #endif //WARPX_PushFieldsEM_K_h
#define AMREX_INLINE
#define AMREX_GPU_HOST_DEVICE
AMREX_GPU_HOST_DEVICE AMREX_INLINE void damp_field_in_guards(amrex::Array4< amrex::Real > const &mf_arr, const int i, const int j, const int k, const int icomp, const int dir, const int n_domain, const int tb_smallEnd, const int tb_bigEnd)
Definition: WarpXPushFieldsEM_K.H:74
AMREX_GPU_HOST_DEVICE AMREX_INLINE amrex::Box constrain_tilebox_to_guards(const amrex::Box &input_tb, const int dir, const int iside, const int n_domain, const int tb_smallEnd, const int tb_bigEnd)
Definition: WarpXPushFieldsEM_K.H:27
AMREX_GPU_HOST_DEVICE Box & setSmall(const IntVect &sm) noexcept
AMREX_GPU_HOST_DEVICE Box & setBig(const IntVect &bg) noexcept
AMREX_GPU_HOST_DEVICE ItType upper_bound(ItType first, ItType last, const ValType &val)
AMREX_GPU_HOST_DEVICE ItType lower_bound(ItType first, ItType last, const ValType &val)
i
Definition: check_interp_points_and_weights.py:174