Block-Structured AMR Software Framework
AMReX_Particle_mod_K.H
Go to the documentation of this file.
1 #ifndef AMREX_PARTICLE_MOD_K_H_
2 #define AMREX_PARTICLE_MOD_K_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_FArrayBox.H>
6 #include <AMReX_REAL.H>
7 #include <cmath>
8 
9 namespace amrex {
10 
11 template <typename P>
12 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
13 void amrex_deposit_cic (P const& p, int nc, amrex::Array4<amrex::Real> const& rho,
14  amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& plo,
15  amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dxi)
16 {
17 #if (AMREX_SPACEDIM == 1)
18  amrex::Real lx = (p.pos(0) - plo[0]) * dxi[0] + Real(0.5);
19 
20  int i = static_cast<int>(amrex::Math::floor(lx));
21 
22  amrex::Real xint = lx - static_cast<Real>(i);
23 
24  amrex::Real sx[2] = {Real(1.0) - xint, xint};
25 
26  for (int ii = 0; ii <= 1; ++ii) {
27  amrex::Gpu::Atomic::AddNoRet(&rho(i+ii-1, 0, 0, 0), static_cast<Real>(sx[ii]*p.rdata(0)));
28  }
29 
30  for (int comp=1; comp < nc; ++comp) {
31  for (int ii = 0; ii <= 1; ++ii) {
32  amrex::Gpu::Atomic::AddNoRet(&rho(i+ii-1, 0, 0, comp),
33  static_cast<Real>(sx[ii]*p.rdata(0)*p.rdata(comp)));
34  }
35  }
36 #elif (AMREX_SPACEDIM == 2)
37  amrex::Real lx = (p.pos(0) - plo[0]) * dxi[0] + Real(0.5);
38  amrex::Real ly = (p.pos(1) - plo[1]) * dxi[1] + Real(0.5);
39 
40  int i = static_cast<int>(amrex::Math::floor(lx));
41  int j = static_cast<int>(amrex::Math::floor(ly));
42 
43  amrex::Real xint = lx - static_cast<Real>(i);
44  amrex::Real yint = ly - static_cast<Real>(j);
45 
46  amrex::Real sx[2] = {Real(1.0) - xint, xint};
47  amrex::Real sy[2] = {Real(1.0) - yint, yint};
48 
49  for (int jj = 0; jj <= 1; ++jj) {
50  for (int ii = 0; ii <= 1; ++ii) {
51  amrex::Gpu::Atomic::AddNoRet(&rho(i+ii-1, j+jj-1, 0, 0),
52  static_cast<Real>(sx[ii]*sy[jj]*p.rdata(0)));
53  }
54  }
55 
56  for (int comp=1; comp < nc; ++comp) {
57  for (int jj = 0; jj <= 1; ++jj) {
58  for (int ii = 0; ii <= 1; ++ii) {
59  amrex::Gpu::Atomic::AddNoRet(&rho(i+ii-1, j+jj-1, 0, comp),
60  static_cast<Real>(sx[ii]*sy[jj]*p.rdata(0)*p.rdata(comp)));
61  }
62  }
63  }
64 
65 #elif (AMREX_SPACEDIM == 3)
66  amrex::Real lx = (p.pos(0) - plo[0]) * dxi[0] + Real(0.5);
67  amrex::Real ly = (p.pos(1) - plo[1]) * dxi[1] + Real(0.5);
68  amrex::Real lz = (p.pos(2) - plo[2]) * dxi[2] + Real(0.5);
69 
70  int i = static_cast<int>(amrex::Math::floor(lx));
71  int j = static_cast<int>(amrex::Math::floor(ly));
72  int k = static_cast<int>(amrex::Math::floor(lz));
73 
74  amrex::Real xint = lx - static_cast<Real>(i);
75  amrex::Real yint = ly - static_cast<Real>(j);
76  amrex::Real zint = lz - static_cast<Real>(k);
77 
78  amrex::Real sx[] = {Real(1.0) - xint, xint};
79  amrex::Real sy[] = {Real(1.0) - yint, yint};
80  amrex::Real sz[] = {Real(1.0) - zint, zint};
81 
82  for (int kk = 0; kk <= 1; ++kk) {
83  for (int jj = 0; jj <= 1; ++jj) {
84  for (int ii = 0; ii <= 1; ++ii) {
85  amrex::Gpu::Atomic::AddNoRet(&rho(i+ii-1, j+jj-1, k+kk-1, 0),
86  static_cast<Real>(sx[ii]*sy[jj]*sz[kk]*p.rdata(0)));
87  }
88  }
89  }
90 
91  for (int comp=1; comp < nc; ++comp) {
92  for (int kk = 0; kk <= 1; ++kk) {
93  for (int jj = 0; jj <= 1; ++jj) {
94  for (int ii = 0; ii <= 1; ++ii) {
95  amrex::Gpu::Atomic::AddNoRet(&rho(i+ii-1, j+jj-1, k+kk-1, comp),
96  static_cast<Real>(sx[ii]*sy[jj]*sz[kk]*p.rdata(0)*p.rdata(comp)));
97  }
98  }
99  }
100  }
101 #else
102  amrex::Abort("Not implemented.");
103 #endif
104 }
105 
106 template <typename P>
107 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
108 void amrex_deposit_particle_dx_cic (P const& p, int nc, amrex::Array4<amrex::Real> const& rho,
109  amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& plo,
110  amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dxi,
111  amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& pdxi)
112 {
113 #if (AMREX_SPACEDIM == 1)
114  amrex::Real factor = (pdxi[0]/dxi[0]);
115 
116  amrex::Real lx = (p.pos(0) - plo[0] - Real(0.5)/pdxi[0]) * dxi[0];
117 
118  amrex::Real hx = (p.pos(0) - plo[0] + Real(0.5)/pdxi[0]) * dxi[0];
119 
120  int lo_x = static_cast<int>(amrex::Math::floor(lx));
121 
122  int hi_x = static_cast<int>(amrex::Math::floor(hx));
123 
124  for (int i = lo_x; i <= hi_x; ++i) {
125  if (i < rho.begin.x || i >= rho.end.x) { continue; }
126  amrex::Real wx = amrex::min(hx - static_cast<Real>(i), amrex::Real(1.0)) - amrex::max(lx - static_cast<Real>(i), amrex::Real(0.0));
127  amrex::Real weight = wx*factor;
128  amrex::Gpu::Atomic::AddNoRet(&rho(i, 0, 0, 0), static_cast<Real>(weight*p.rdata(0)));
129  }
130 
131  for (int comp = 1; comp < nc; ++comp)
132  {
133  for (int i = lo_x; i <= hi_x; ++i) {
134  if (i < rho.begin.x || i >= rho.end.x) { continue; }
135  amrex::Real wx = amrex::min(hx - static_cast<Real>(i), amrex::Real(1.0)) - amrex::max(lx - static_cast<Real>(i), amrex::Real(0.0));
136  amrex::Real weight = wx*factor;
137  amrex::Gpu::Atomic::AddNoRet(&rho(i, 0, 0, comp), static_cast<Real>(weight*p.rdata(0)*p.rdata(comp)));
138  }
139  }
140 
141 #elif (AMREX_SPACEDIM == 2)
142  amrex::Real factor = (pdxi[0]/dxi[0])*(pdxi[1]/dxi[1]);
143 
144  amrex::Real lx = (p.pos(0) - plo[0] - Real(0.5)/pdxi[0]) * dxi[0];
145  amrex::Real ly = (p.pos(1) - plo[1] - Real(0.5)/pdxi[1]) * dxi[1];
146 
147  amrex::Real hx = (p.pos(0) - plo[0] + Real(0.5)/pdxi[0]) * dxi[0];
148  amrex::Real hy = (p.pos(1) - plo[1] + Real(0.5)/pdxi[1]) * dxi[1];
149 
150  int lo_x = static_cast<int>(amrex::Math::floor(lx));
151  int lo_y = static_cast<int>(amrex::Math::floor(ly));
152 
153  int hi_x = static_cast<int>(amrex::Math::floor(hx));
154  int hi_y = static_cast<int>(amrex::Math::floor(hy));
155 
156  for (int j = lo_y; j <= hi_y; ++j) {
157  if (j < rho.begin.y || j >= rho.end.y) { continue; }
158  amrex::Real wy = amrex::min(hy - static_cast<Real>(j), amrex::Real(1.0)) - amrex::max(ly - static_cast<Real>(j), amrex::Real(0.0));
159  for (int i = lo_x; i <= hi_x; ++i) {
160  if (i < rho.begin.x || i >= rho.end.x) { continue; }
161  amrex::Real wx = amrex::min(hx - static_cast<Real>(i), amrex::Real(1.0)) - amrex::max(lx - static_cast<Real>(i), amrex::Real(0.0));
162  amrex::Real weight = wx*wy*factor;
163  amrex::Gpu::Atomic::AddNoRet(&rho(i, j, 0, 0), static_cast<Real>(weight*p.rdata(0)));
164  }
165  }
166 
167  for (int comp = 1; comp < nc; ++comp) {
168  for (int j = lo_y; j <= hi_y; ++j) {
169  if (j < rho.begin.y || j >= rho.end.y) { continue; }
170  amrex::Real wy = amrex::min(hy - static_cast<Real>(j), amrex::Real(1.0)) - amrex::max(ly - static_cast<Real>(j), amrex::Real(0.0));
171  for (int i = lo_x; i <= hi_x; ++i) {
172  if (i < rho.begin.x || i >= rho.end.x) { continue; }
173  amrex::Real wx = amrex::min(hx - static_cast<Real>(i), amrex::Real(1.0)) - amrex::max(lx - static_cast<Real>(i), amrex::Real(0.0));
174  amrex::Real weight = wx*wy*factor;
175  amrex::Gpu::Atomic::AddNoRet(&rho(i, j, 0, comp), static_cast<Real>(weight*p.rdata(0)*p.rdata(comp)));
176  }
177  }
178  }
179 
180 #elif (AMREX_SPACEDIM == 3)
181  amrex::Real factor = (pdxi[0]/dxi[0])*(pdxi[1]/dxi[1])*(pdxi[2]/dxi[2]);
182 
183  amrex::Real lx = (p.pos(0) - plo[0] - Real(0.5)/pdxi[0]) * dxi[0];
184  amrex::Real ly = (p.pos(1) - plo[1] - Real(0.5)/pdxi[1]) * dxi[1];
185  amrex::Real lz = (p.pos(2) - plo[2] - Real(0.5)/pdxi[2]) * dxi[2];
186 
187  amrex::Real hx = (p.pos(0) - plo[0] + Real(0.5)/pdxi[0]) * dxi[0];
188  amrex::Real hy = (p.pos(1) - plo[1] + Real(0.5)/pdxi[1]) * dxi[1];
189  amrex::Real hz = (p.pos(2) - plo[2] + Real(0.5)/pdxi[2]) * dxi[2];
190 
191  int lo_x = static_cast<int>(amrex::Math::floor(lx));
192  int lo_y = static_cast<int>(amrex::Math::floor(ly));
193  int lo_z = static_cast<int>(amrex::Math::floor(lz));
194 
195  int hi_x = static_cast<int>(amrex::Math::floor(hx));
196  int hi_y = static_cast<int>(amrex::Math::floor(hy));
197  int hi_z = static_cast<int>(amrex::Math::floor(hz));
198 
199  for (int k = lo_z; k <= hi_z; ++k) {
200  if (k < rho.begin.z || k >= rho.end.z) { continue; }
201  amrex::Real wz = amrex::min(hz - static_cast<Real>(k), amrex::Real(1.0)) - amrex::max(lz - static_cast<Real>(k), amrex::Real(0.0));
202  for (int j = lo_y; j <= hi_y; ++j) {
203  if (j < rho.begin.y || j >= rho.end.y) { continue; }
204  amrex::Real wy = amrex::min(hy - static_cast<Real>(j), amrex::Real(1.0)) - amrex::max(ly - static_cast<Real>(j), amrex::Real(0.0));
205  for (int i = lo_x; i <= hi_x; ++i) {
206  if (i < rho.begin.x || i >= rho.end.x) { continue; }
207  amrex::Real wx = amrex::min(hx - static_cast<Real>(i), amrex::Real(1.0)) - amrex::max(lx - static_cast<Real>(i), amrex::Real(0.0));
208  amrex::Real weight = wx*wy*wz*factor;
209  amrex::Gpu::Atomic::AddNoRet(&rho(i, j, k, 0), static_cast<Real>(weight*p.rdata(0)));
210  }
211  }
212  }
213 
214  for (int comp = 1; comp < nc; ++comp) {
215  for (int k = lo_z; k <= hi_z; ++k) {
216  if (k < rho.begin.z || k >= rho.end.z) { continue; }
217  amrex::Real wz = amrex::min(hz - static_cast<Real>(k), amrex::Real(1.0)) - amrex::max(lz - static_cast<Real>(k), amrex::Real(0.0));
218  for (int j = lo_y; j <= hi_y; ++j) {
219  if (j < rho.begin.y || j >= rho.end.y) { continue; }
220  amrex::Real wy = amrex::min(hy - static_cast<Real>(j), amrex::Real(1.0)) - amrex::max(ly - static_cast<Real>(j), amrex::Real(0.0));
221  for (int i = lo_x; i <= hi_x; ++i) {
222  if (i < rho.begin.x || i >= rho.end.x) { continue; }
223  amrex::Real wx = amrex::min(hx - static_cast<Real>(i), amrex::Real(1.0)) - amrex::max(lx - static_cast<Real>(i), amrex::Real(0.0));
224  amrex::Real weight = wx*wy*wz*factor;
225  amrex::Gpu::Atomic::AddNoRet(&rho(i, j, k, comp), static_cast<Real>(weight*p.rdata(0)*p.rdata(comp)));
226  }
227  }
228  }
229  }
230 #else
231  amrex::Abort("Not implemented.")
232 #endif
233 }
234 
235 }
236 
237 #endif
AMREX_FORCE_INLINE
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
AMReX_FArrayBox.H
AMREX_GPU_HOST_DEVICE
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
AMReX_REAL.H
amrex::Gpu::Atomic::AddNoRet
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void AddNoRet(T *sum, T value) noexcept
Definition: AMReX_GpuAtomic.H:281
amrex::openbc::P
static constexpr int P
Definition: AMReX_OpenBC.H:14
amrex
Definition: AMReX_Amr.cpp:49
amrex::amrex_deposit_cic
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void amrex_deposit_cic(P const &p, int nc, amrex::Array4< amrex::Real > const &rho, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi)
Definition: AMReX_Particle_mod_K.H:13
amrex::max
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE const T & max(const T &a, const T &b) noexcept
Definition: AMReX_Algorithm.H:35
amrex::min
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE const T & min(const T &a, const T &b) noexcept
Definition: AMReX_Algorithm.H:21
amrex::amrex_deposit_particle_dx_cic
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void amrex_deposit_particle_dx_cic(P const &p, int nc, amrex::Array4< amrex::Real > const &rho, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &pdxi)
Definition: AMReX_Particle_mod_K.H:108
amrex::Abort
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:221
amrex::Array4
Definition: AMReX_Array4.H:61
amrex::Array4::end
Dim3 end
Definition: AMReX_Array4.H:67
amrex::Dim3::x
int x
Definition: AMReX_Dim3.H:12
amrex::Dim3::z
int z
Definition: AMReX_Dim3.H:12
amrex::Dim3::y
int y
Definition: AMReX_Dim3.H:12
amrex::GpuArray
Definition: AMReX_Array.H:33