WarpX
WarpXComm_K.H
Go to the documentation of this file.
1 /* Copyright 2019 Axel Huebl, Weiqun Zhang
2  *
3  * This file is part of WarpX.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 #ifndef WARPX_COMM_K_H_
8 #define WARPX_COMM_K_H_
9 
10 #include <AMReX.H>
11 #include <AMReX_FArrayBox.H>
12 
13 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
14 void warpx_interp (int j, int k, int l,
15  amrex::Array4<amrex::Real > const& arr_aux,
16  amrex::Array4<amrex::Real const> const& arr_fine,
17  amrex::Array4<amrex::Real const> const& arr_coarse,
18  const amrex::IntVect& arr_stag,
19  const amrex::IntVect& rr)
20 {
21  using namespace amrex;
22 
23  // NOTE Indices (j,k,l) in the following refer to (z,-,-) in 1D, (x,z,-) in 2D, and (x,y,z) in 3D
24 
25  // Refinement ratio
26  const int rj = rr[0];
27  const int rk = (AMREX_SPACEDIM == 1) ? 1 : rr[1];
28  const int rl = (AMREX_SPACEDIM <= 2) ? 1 : rr[2];
29 
30  // Staggering (0: cell-centered; 1: nodal)
31  const int sj = arr_stag[0];
32  const int sk = (AMREX_SPACEDIM == 1) ? 0 : arr_stag[1];
33  const int sl = (AMREX_SPACEDIM <= 2) ? 0 : arr_stag[2];
34 
35  // Number of points used for interpolation from coarse grid to fine grid
36  const int nj = (sj == 0) ? 1 : 2;
37  const int nk = (sk == 0) ? 1 : 2;
38  const int nl = (sl == 0) ? 1 : 2;
39 
40  const int jc = amrex::coarsen(j, rj);
41  const int kc = amrex::coarsen(k, rk);
42  const int lc = amrex::coarsen(l, rl);
43 
44  amrex::Real wj;
45  amrex::Real wk;
46  amrex::Real wl;
47 
48  // Interpolate from coarse grid to fine grid using either 1 point with weight 1, if both grids
49  // are cell-centered, or 2 points with weights depending on the distance, if both grids are nodal
50  amrex::Real res = 0.0_rt;
51  for (int jj = 0; jj < nj; jj++) {
52  for (int kk = 0; kk < nk; kk++) {
53  for (int ll = 0; ll < nl; ll++) {
54  wj = (sj == 0) ? 1.0_rt : (rj - amrex::Math::abs(j - (jc + jj) * rj))
55  / static_cast<amrex::Real>(rj);
56  wk = (sk == 0) ? 1.0_rt : (rk - amrex::Math::abs(k - (kc + kk) * rk))
57  / static_cast<amrex::Real>(rk);
58  wl = (sl == 0) ? 1.0_rt : (rl - amrex::Math::abs(l - (lc + ll) * rl))
59  / static_cast<amrex::Real>(rl);
60  res += wj * wk * wl * arr_coarse(jc+jj,kc+kk,lc+ll);
61  }
62  }
63  }
64  arr_aux(j,k,l) = arr_fine(j,k,l) + res;
65 }
66 
67 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
68 void warpx_interp_nd_bfield_x (int j, int k, int l,
69  amrex::Array4<amrex::Real> const& Bxa,
70  amrex::Array4<amrex::Real const> const& Bxf,
71  amrex::Array4<amrex::Real const> const& Bxc,
72  amrex::Array4<amrex::Real const> const& Bxg)
73 {
74  using namespace amrex;
75 
76  // Pad Bxf with zeros beyond ghost cells for out-of-bound accesses
77  const auto Bxf_zeropad = [Bxf] (const int jj, const int kk, const int ll) noexcept
78  {
79  return Bxf.contains(jj,kk,ll) ? Bxf(jj,kk,ll) : 0.0_rt;
80  };
81 
82  int jg = amrex::coarsen(j,2);
83  Real wx = (j == jg*2) ? 0.0_rt : 0.5_rt;
84  Real owx = 1.0_rt-wx;
85 
86  int kg = amrex::coarsen(k,2);
87  Real wy = (k == kg*2) ? 0.0_rt : 0.5_rt;
88  Real owy = 1.0_rt-wy;
89 
90 #if defined(WARPX_DIM_1D_Z)
91 
92  // interp from coarse nodal to fine nodal
93  Real bg = owx * Bxg(jg ,0,0)
94  + wx * Bxg(jg+1,0,0);
95 
96  // interp from coarse staggered to fine nodal
97  Real bc = owx * Bxc(jg ,0,0)
98  + wx * Bxc(jg+1,0,0);
99 
100  // interp from fine staggered to fine nodal
101  Real bf = 0.5_rt*(Bxf_zeropad(j-1,0,0) + Bxf_zeropad(j,0,0));
102  amrex::ignore_unused(owy);
103 
104 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
105 
106  // interp from coarse nodal to fine nodal
107  Real bg = owx * owy * Bxg(jg ,kg ,0)
108  + owx * wy * Bxg(jg ,kg+1,0)
109  + wx * owy * Bxg(jg+1,kg ,0)
110  + wx * wy * Bxg(jg+1,kg+1,0);
111 
112  // interp from coarse staggered to fine nodal
113  wy = 0.5_rt-wy; owy = 1.0_rt-wy;
114  Real bc = owx * owy * Bxc(jg ,kg ,0)
115  + owx * wy * Bxc(jg ,kg-1,0)
116  + wx * owy * Bxc(jg+1,kg ,0)
117  + wx * wy * Bxc(jg+1,kg-1,0);
118 
119  // interp from fine staggered to fine nodal
120  Real bf = 0.5_rt*(Bxf_zeropad(j,k-1,0) + Bxf_zeropad(j,k,0));
121 
122 #else
123 
124  int lg = amrex::coarsen(l,2);
125  Real wz = (l == lg*2) ? 0.0_rt : 0.5_rt;
126  Real owz = 1.0_rt-wz;
127 
128  // interp from coarse nodal to fine nodal
129  Real bg = owx * owy * owz * Bxg(jg ,kg ,lg )
130  + wx * owy * owz * Bxg(jg+1,kg ,lg )
131  + owx * wy * owz * Bxg(jg ,kg+1,lg )
132  + wx * wy * owz * Bxg(jg+1,kg+1,lg )
133  + owx * owy * wz * Bxg(jg ,kg ,lg+1)
134  + wx * owy * wz * Bxg(jg+1,kg ,lg+1)
135  + owx * wy * wz * Bxg(jg ,kg+1,lg+1)
136  + wx * wy * wz * Bxg(jg+1,kg+1,lg+1);
137 
138  // interp from coarse staggered to fine nodal
139  wy = 0.5_rt-wy; owy = 1.0_rt-wy;
140  wz = 0.5_rt-wz; owz = 1.0_rt-wz;
141  Real bc = owx * owy * owz * Bxc(jg ,kg ,lg )
142  + wx * owy * owz * Bxc(jg+1,kg ,lg )
143  + owx * wy * owz * Bxc(jg ,kg-1,lg )
144  + wx * wy * owz * Bxc(jg+1,kg-1,lg )
145  + owx * owy * wz * Bxc(jg ,kg ,lg-1)
146  + wx * owy * wz * Bxc(jg+1,kg ,lg-1)
147  + owx * wy * wz * Bxc(jg ,kg-1,lg-1)
148  + wx * wy * wz * Bxc(jg+1,kg-1,lg-1);
149 
150  // interp from fine stagged to fine nodal
151  Real bf = 0.25_rt*(Bxf_zeropad(j,k-1,l-1) + Bxf_zeropad(j,k,l-1) + Bxf_zeropad(j,k-1,l) + Bxf_zeropad(j,k,l));
152 #endif
153 
154  Bxa(j,k,l) = bg + (bf-bc);
155 }
156 
157 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
158 void warpx_interp_nd_bfield_y (int j, int k, int l,
159  amrex::Array4<amrex::Real> const& Bya,
160  amrex::Array4<amrex::Real const> const& Byf,
161  amrex::Array4<amrex::Real const> const& Byc,
162  amrex::Array4<amrex::Real const> const& Byg)
163 {
164  using namespace amrex;
165 
166  // Pad Byf with zeros beyond ghost cells for out-of-bound accesses
167  const auto Byf_zeropad = [Byf] (const int jj, const int kk, const int ll) noexcept
168  {
169  return Byf.contains(jj,kk,ll) ? Byf(jj,kk,ll) : 0.0_rt;
170  };
171 
172  int jg = amrex::coarsen(j,2);
173  Real wx = (j == jg*2) ? 0.0_rt : 0.5_rt;
174  Real owx = 1.0_rt-wx;
175 
176  int kg = amrex::coarsen(k,2);
177  Real wy = (k == kg*2) ? 0.0_rt : 0.5_rt;
178  Real owy = 1.0_rt-wy;
179 
180 #if defined(WARPX_DIM_1D_Z)
181 
182  // interp from coarse nodal to fine nodal
183  Real bg = owx * Byg(jg ,0,0)
184  + wx * Byg(jg+1,0,0);
185 
186  // interp from coarse staggered to fine nodal
187  Real bc = owx * Byc(jg ,0,0)
188  + wx * Byc(jg+1,0,0);
189 
190  // interp from fine staggered to fine nodal
191  Real bf = 0.5_rt*(Byf_zeropad(j-1,0,0) + Byf_zeropad(j,0,0));
192  amrex::ignore_unused(owy);
193 
194 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
195 
196  // interp from coarse nodal to fine nodal
197  Real bg = owx * owy * Byg(jg ,kg ,0)
198  + owx * wy * Byg(jg ,kg+1,0)
199  + wx * owy * Byg(jg+1,kg ,0)
200  + wx * wy * Byg(jg+1,kg+1,0);
201 
202  // interp from coarse stagged (cell-centered for By) to fine nodal
203  wx = 0.5_rt-wx; owx = 1.0_rt-wx;
204  wy = 0.5_rt-wy; owy = 1.0_rt-wy;
205  Real bc = owx * owy * Byc(jg ,kg ,0)
206  + owx * wy * Byc(jg ,kg-1,0)
207  + wx * owy * Byc(jg-1,kg ,0)
208  + wx * wy * Byc(jg-1,kg-1,0);
209 
210  // interp form fine stagger (cell-centered for By) to fine nodal
211  Real bf = 0.25_rt*(Byf_zeropad(j,k,0) + Byf_zeropad(j-1,k,0) + Byf_zeropad(j,k-1,0) + Byf_zeropad(j-1,k-1,0));
212 
213 #else
214 
215  int lg = amrex::coarsen(l,2);
216  Real wz = (l == lg*2) ? 0.0_rt : 0.5_rt;
217  Real owz = 1.0_rt-wz;
218 
219  // interp from coarse nodal to fine nodal
220  Real bg = owx * owy * owz * Byg(jg ,kg ,lg )
221  + wx * owy * owz * Byg(jg+1,kg ,lg )
222  + owx * wy * owz * Byg(jg ,kg+1,lg )
223  + wx * wy * owz * Byg(jg+1,kg+1,lg )
224  + owx * owy * wz * Byg(jg ,kg ,lg+1)
225  + wx * owy * wz * Byg(jg+1,kg ,lg+1)
226  + owx * wy * wz * Byg(jg ,kg+1,lg+1)
227  + wx * wy * wz * Byg(jg+1,kg+1,lg+1);
228 
229  // interp from coarse staggered to fine nodal
230  wx = 0.5_rt-wx; owx = 1.0_rt-wx;
231  wz = 0.5_rt-wz; owz = 1.0_rt-wz;
232  Real bc = owx * owy * owz * Byc(jg ,kg ,lg )
233  + wx * owy * owz * Byc(jg-1,kg ,lg )
234  + owx * wy * owz * Byc(jg ,kg+1,lg )
235  + wx * wy * owz * Byc(jg-1,kg+1,lg )
236  + owx * owy * wz * Byc(jg ,kg ,lg-1)
237  + wx * owy * wz * Byc(jg-1,kg ,lg-1)
238  + owx * wy * wz * Byc(jg ,kg+1,lg-1)
239  + wx * wy * wz * Byc(jg-1,kg+1,lg-1);
240 
241  // interp from fine stagged to fine nodal
242  Real bf = 0.25_rt*(Byf_zeropad(j-1,k,l-1) + Byf_zeropad(j,k,l-1) + Byf_zeropad(j-1,k,l) + Byf_zeropad(j,k,l));
243 
244 #endif
245 
246  Bya(j,k,l) = bg + (bf-bc);
247 }
248 
249 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
250 void warpx_interp_nd_bfield_z (int j, int k, int l,
251  amrex::Array4<amrex::Real> const& Bza,
252  amrex::Array4<amrex::Real const> const& Bzf,
253  amrex::Array4<amrex::Real const> const& Bzc,
254  amrex::Array4<amrex::Real const> const& Bzg)
255 {
256  using namespace amrex;
257 
258  // Pad Bzf with zeros beyond ghost cells for out-of-bound accesses
259  const auto Bzf_zeropad = [Bzf] (const int jj, const int kk, const int ll) noexcept
260  {
261  return Bzf.contains(jj,kk,ll) ? Bzf(jj,kk,ll) : 0.0_rt;
262  };
263 
264  int jg = amrex::coarsen(j,2);
265  Real wx = (j == jg*2) ? 0.0_rt : 0.5_rt;
266  Real owx = 1.0_rt-wx;
267 
268  int kg = amrex::coarsen(k,2);
269  Real wy = (k == kg*2) ? 0.0_rt : 0.5_rt;
270  Real owy = 1.0_rt-wy;
271 
272 #if defined(WARPX_DIM_1D_Z)
273 
274  // interp from coarse nodal to fine nodal
275  Real bg = owx * Bzg(jg ,0,0)
276  + wx * Bzg(jg+1,0,0);
277 
278  // interp from coarse staggered to fine nodal
279  Real bc = owx * Bzc(jg ,0,0)
280  + wx * Bzc(jg+1,0,0);
281 
282  // interp from fine staggered to fine nodal
283  Real bf = Bzf_zeropad(j,0,0);
284  amrex::ignore_unused(owy);
285 
286 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
287 
288  // interp from coarse nodal to fine nodal
289  Real bg = owx * owy * Bzg(jg ,kg ,0)
290  + owx * wy * Bzg(jg ,kg+1,0)
291  + wx * owy * Bzg(jg+1,kg ,0)
292  + wx * wy * Bzg(jg+1,kg+1,0);
293 
294  // interp from coarse staggered to fine nodal
295  wx = 0.5_rt-wx; owx = 1.0_rt-wx;
296  Real bc = owx * owy * Bzc(jg ,kg ,0)
297  + owx * wy * Bzc(jg ,kg+1,0)
298  + wx * owy * Bzc(jg-1,kg ,0)
299  + wx * wy * Bzc(jg-1,kg+1,0);
300 
301  // interp from fine staggered to fine nodal
302  Real bf = 0.5_rt*(Bzf_zeropad(j-1,k,0) + Bzf_zeropad(j,k,0));
303 
304 #else
305 
306  int lg = amrex::coarsen(l,2);
307  Real wz = (l == lg*2) ? 0.0_rt : 0.5_rt;
308  Real owz = 1.0_rt-wz;
309 
310  // interp from coarse nodal to fine nodal
311  Real bg = owx * owy * owz * Bzg(jg ,kg ,lg )
312  + wx * owy * owz * Bzg(jg+1,kg ,lg )
313  + owx * wy * owz * Bzg(jg ,kg+1,lg )
314  + wx * wy * owz * Bzg(jg+1,kg+1,lg )
315  + owx * owy * wz * Bzg(jg ,kg ,lg+1)
316  + wx * owy * wz * Bzg(jg+1,kg ,lg+1)
317  + owx * wy * wz * Bzg(jg ,kg+1,lg+1)
318  + wx * wy * wz * Bzg(jg+1,kg+1,lg+1);
319 
320  // interp from coarse staggered to fine nodal
321  wx = 0.5_rt-wx; owx = 1.0_rt-wx;
322  wy = 0.5_rt-wy; owy = 1.0_rt-wy;
323  Real bc = owx * owy * owz * Bzc(jg ,kg ,lg )
324  + wx * owy * owz * Bzc(jg-1,kg ,lg )
325  + owx * wy * owz * Bzc(jg ,kg-1,lg )
326  + wx * wy * owz * Bzc(jg-1,kg-1,lg )
327  + owx * owy * wz * Bzc(jg ,kg ,lg+1)
328  + wx * owy * wz * Bzc(jg-1,kg ,lg+1)
329  + owx * wy * wz * Bzc(jg ,kg-1,lg+1)
330  + wx * wy * wz * Bzc(jg-1,kg-1,lg+1);
331 
332  // interp from fine stagged to fine nodal
333  Real bf = 0.25_rt*(Bzf_zeropad(j-1,k-1,l) + Bzf_zeropad(j,k-1,l) + Bzf_zeropad(j-1,k,l) + Bzf_zeropad(j,k,l));
334 
335 #endif
336 
337  Bza(j,k,l) = bg + (bf-bc);
338 }
339 
340 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
341 void warpx_interp_nd_efield_x (int j, int k, int l,
342  amrex::Array4<amrex::Real> const& Exa,
343  amrex::Array4<amrex::Real const> const& Exf,
344  amrex::Array4<amrex::Real const> const& Exc,
345  amrex::Array4<amrex::Real const> const& Exg)
346 {
347  using namespace amrex;
348 
349  // Pad Exf with zeros beyond ghost cells for out-of-bound accesses
350  const auto Exf_zeropad = [Exf] (const int jj, const int kk, const int ll) noexcept
351  {
352  return Exf.contains(jj,kk,ll) ? Exf(jj,kk,ll) : 0.0_rt;
353  };
354 
355  int jg = amrex::coarsen(j,2);
356  Real wx = (j == jg*2) ? 0.0_rt : 0.5_rt;
357  Real owx = 1.0_rt-wx;
358 
359  int kg = amrex::coarsen(k,2);
360  Real wy = (k == kg*2) ? 0.0_rt : 0.5_rt;
361  Real owy = 1.0_rt-wy;
362 
363 #if defined(WARPX_DIM_1D_Z)
364 
365  // interp from coarse nodal to fine nodal
366  Real eg = owx * Exg(jg ,0,0)
367  + wx * Exg(jg+1,0,0);
368 
369  // interp from coarse staggered to fine nodal
370  Real ec = owx * Exc(jg ,0,0)
371  + wx * Exc(jg+1,0,0);
372 
373  // interp from fine staggered to fine nodal
374  Real ef = Exf_zeropad(j,0,0);
375  amrex::ignore_unused(owy);
376 
377 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
378 
379  // interp from coarse nodal to fine nodal
380  Real eg = owx * owy * Exg(jg ,kg ,0)
381  + owx * wy * Exg(jg ,kg+1,0)
382  + wx * owy * Exg(jg+1,kg ,0)
383  + wx * wy * Exg(jg+1,kg+1,0);
384 
385  // interp from coarse staggered to fine nodal
386  wx = 0.5_rt-wx; owx = 1.0_rt-wx;
387  Real ec = owx * owy * Exc(jg ,kg ,0)
388  + owx * wy * Exc(jg ,kg+1,0)
389  + wx * owy * Exc(jg-1,kg ,0)
390  + wx * wy * Exc(jg-1,kg+1,0);
391 
392  // interp from fine staggered to fine nodal
393  Real ef = 0.5_rt*(Exf_zeropad(j-1,k,0) + Exf_zeropad(j,k,0));
394 
395 #else
396 
397  int lg = amrex::coarsen(l,2);
398  Real wz = (l == lg*2) ? 0.0 : 0.5;
399  Real owz = 1.0_rt-wz;
400 
401  // interp from coarse nodal to fine nodal
402  Real eg = owx * owy * owz * Exg(jg ,kg ,lg )
403  + wx * owy * owz * Exg(jg+1,kg ,lg )
404  + owx * wy * owz * Exg(jg ,kg+1,lg )
405  + wx * wy * owz * Exg(jg+1,kg+1,lg )
406  + owx * owy * wz * Exg(jg ,kg ,lg+1)
407  + wx * owy * wz * Exg(jg+1,kg ,lg+1)
408  + owx * wy * wz * Exg(jg ,kg+1,lg+1)
409  + wx * wy * wz * Exg(jg+1,kg+1,lg+1);
410 
411  // interp from coarse staggered to fine nodal
412  wx = 0.5_rt-wx; owx = 1.0_rt-wx;
413  Real ec = owx * owy * owz * Exc(jg ,kg ,lg )
414  + wx * owy * owz * Exc(jg-1,kg ,lg )
415  + owx * wy * owz * Exc(jg ,kg+1,lg )
416  + wx * wy * owz * Exc(jg-1,kg+1,lg )
417  + owx * owy * wz * Exc(jg ,kg ,lg+1)
418  + wx * owy * wz * Exc(jg-1,kg ,lg+1)
419  + owx * wy * wz * Exc(jg ,kg+1,lg+1)
420  + wx * wy * wz * Exc(jg-1,kg+1,lg+1);
421 
422  // interp from fine staggered to fine nodal
423  Real ef = 0.5_rt*(Exf_zeropad(j-1,k,l) + Exf_zeropad(j,k,l));
424 
425 #endif
426 
427  Exa(j,k,l) = eg + (ef-ec);
428 }
429 
430 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
431 void warpx_interp_nd_efield_y (int j, int k, int l,
432  amrex::Array4<amrex::Real> const& Eya,
433  amrex::Array4<amrex::Real const> const& Eyf,
434  amrex::Array4<amrex::Real const> const& Eyc,
435  amrex::Array4<amrex::Real const> const& Eyg)
436 {
437  using namespace amrex;
438 
439  // Pad Eyf with zeros beyond ghost cells for out-of-bound accesses
440  const auto Eyf_zeropad = [Eyf] (const int jj, const int kk, const int ll) noexcept
441  {
442  return Eyf.contains(jj,kk,ll) ? Eyf(jj,kk,ll) : 0.0_rt;
443  };
444 
445  int jg = amrex::coarsen(j,2);
446  Real wx = (j == jg*2) ? 0.0_rt : 0.5_rt;
447  Real owx = 1.0_rt-wx;
448 
449  int kg = amrex::coarsen(k,2);
450  Real wy = (k == kg*2) ? 0.0_rt : 0.5_rt;
451  Real owy = 1.0_rt-wy;
452 
453 
454 
455 #if defined(WARPX_DIM_1D_Z)
456 
457  // interp from coarse nodal to fine nodal
458  Real eg = owx * Eyg(jg ,0,0)
459  + wx * Eyg(jg+1,0,0);
460 
461  // interp from coarse staggered to fine nodal
462  Real ec = owx * Eyc(jg ,0,0)
463  + wx * Eyc(jg+1,0,0);
464 
465  // interp from fine staggered to fine nodal
466  Real ef = Eyf_zeropad(j,0,0);
467  amrex::ignore_unused(owy);
468 
469 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
470 
471  // interp from coarse nodal to fine nodal
472  Real eg = owx * owy * Eyg(jg ,kg ,0)
473  + owx * wy * Eyg(jg ,kg+1,0)
474  + wx * owy * Eyg(jg+1,kg ,0)
475  + wx * wy * Eyg(jg+1,kg+1,0);
476 
477  // interp from coarse staggered to fine nodal (Eyc is actually nodal)
478  Real ec = owx * owy * Eyc(jg ,kg ,0)
479  + owx * wy * Eyc(jg ,kg+1,0)
480  + wx * owy * Eyc(jg+1,kg ,0)
481  + wx * wy * Eyc(jg+1,kg+1,0);
482 
483  // interp from fine staggered to fine nodal
484  Real ef = Eyf_zeropad(j,k,0);
485 
486 #else
487 
488  int lg = amrex::coarsen(l,2);
489  Real wz = (l == lg*2) ? 0.0 : 0.5;
490  Real owz = 1.0_rt-wz;
491 
492  // interp from coarse nodal to fine nodal
493  Real eg = owx * owy * owz * Eyg(jg ,kg ,lg )
494  + wx * owy * owz * Eyg(jg+1,kg ,lg )
495  + owx * wy * owz * Eyg(jg ,kg+1,lg )
496  + wx * wy * owz * Eyg(jg+1,kg+1,lg )
497  + owx * owy * wz * Eyg(jg ,kg ,lg+1)
498  + wx * owy * wz * Eyg(jg+1,kg ,lg+1)
499  + owx * wy * wz * Eyg(jg ,kg+1,lg+1)
500  + wx * wy * wz * Eyg(jg+1,kg+1,lg+1);
501 
502  // interp from coarse staggered to fine nodal
503  wy = 0.5_rt-wy; owy = 1.0_rt-wy;
504  Real ec = owx * owy * owz * Eyc(jg ,kg ,lg )
505  + wx * owy * owz * Eyc(jg+1,kg ,lg )
506  + owx * wy * owz * Eyc(jg ,kg-1,lg )
507  + wx * wy * owz * Eyc(jg+1,kg-1,lg )
508  + owx * owy * wz * Eyc(jg ,kg ,lg+1)
509  + wx * owy * wz * Eyc(jg+1,kg ,lg+1)
510  + owx * wy * wz * Eyc(jg ,kg-1,lg+1)
511  + wx * wy * wz * Eyc(jg+1,kg-1,lg+1);
512 
513  // interp from fine staggered to fine nodal
514  Real ef = 0.5_rt*(Eyf_zeropad(j,k-1,l) + Eyf_zeropad(j,k,l));
515 
516 #endif
517 
518  Eya(j,k,l) = eg + (ef-ec);
519 }
520 
521 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
522 void warpx_interp_nd_efield_z (int j, int k, int l,
523  amrex::Array4<amrex::Real> const& Eza,
524  amrex::Array4<amrex::Real const> const& Ezf,
525  amrex::Array4<amrex::Real const> const& Ezc,
526  amrex::Array4<amrex::Real const> const& Ezg)
527 {
528  using namespace amrex;
529 
530  // Pad Ezf with zeros beyond ghost cells for out-of-bound accesses
531  const auto Ezf_zeropad = [Ezf] (const int jj, const int kk, const int ll) noexcept
532  {
533  return Ezf.contains(jj,kk,ll) ? Ezf(jj,kk,ll) : 0.0_rt;
534  };
535 
536  int jg = amrex::coarsen(j,2);
537  Real wx = (j == jg*2) ? 0.0_rt : 0.5_rt;
538  Real owx = 1.0_rt-wx;
539 
540  int kg = amrex::coarsen(k,2);
541  Real wy = (k == kg*2) ? 0.0_rt : 0.5_rt;
542  Real owy = 1.0_rt-wy;
543 
544 #if defined(WARPX_DIM_1D_Z)
545 
546  // interp from coarse nodal to fine nodal
547  Real eg = owx * Ezg(jg ,0,0)
548  + wx * Ezg(jg+1,0,0);
549 
550  // interp from coarse staggered to fine nodal
551  Real ec = owx * Ezc(jg ,0,0)
552  + wx * Ezc(jg+1,0,0);
553 
554  // interp from fine staggered to fine nodal
555  Real ef = 0.5_rt*(Ezf_zeropad(j-1,0,0) + Ezf_zeropad(j,0,0));
556  amrex::ignore_unused(owy);
557 
558 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
559 
560  // interp from coarse nodal to fine nodal
561  Real eg = owx * owy * Ezg(jg ,kg ,0)
562  + owx * wy * Ezg(jg ,kg+1,0)
563  + wx * owy * Ezg(jg+1,kg ,0)
564  + wx * wy * Ezg(jg+1,kg+1,0);
565 
566  // interp from coarse stagged to fine nodal
567  wy = 0.5_rt-wy; owy = 1.0_rt-wy;
568  Real ec = owx * owy * Ezc(jg ,kg ,0)
569  + owx * wy * Ezc(jg ,kg-1,0)
570  + wx * owy * Ezc(jg+1,kg ,0)
571  + wx * wy * Ezc(jg+1,kg-1,0);
572 
573  // interp from fine staggered to fine nodal
574  Real ef = 0.5_rt*(Ezf_zeropad(j,k-1,0) + Ezf_zeropad(j,k,0));
575 
576 #else
577 
578  int lg = amrex::coarsen(l,2);
579  Real wz = (l == lg*2) ? 0.0_rt : 0.5_rt;
580  Real owz = 1.0_rt-wz;
581 
582  // interp from coarse nodal to fine nodal
583  Real eg = owx * owy * owz * Ezg(jg ,kg ,lg )
584  + wx * owy * owz * Ezg(jg+1,kg ,lg )
585  + owx * wy * owz * Ezg(jg ,kg+1,lg )
586  + wx * wy * owz * Ezg(jg+1,kg+1,lg )
587  + owx * owy * wz * Ezg(jg ,kg ,lg+1)
588  + wx * owy * wz * Ezg(jg+1,kg ,lg+1)
589  + owx * wy * wz * Ezg(jg ,kg+1,lg+1)
590  + wx * wy * wz * Ezg(jg+1,kg+1,lg+1);
591 
592  // interp from coarse staggered to fine nodal
593  wz = 0.5_rt-wz; owz = 1.0_rt-wz;
594  Real ec = owx * owy * owz * Ezc(jg ,kg ,lg )
595  + wx * owy * owz * Ezc(jg+1,kg ,lg )
596  + owx * wy * owz * Ezc(jg ,kg+1,lg )
597  + wx * wy * owz * Ezc(jg+1,kg+1,lg )
598  + owx * owy * wz * Ezc(jg ,kg ,lg-1)
599  + wx * owy * wz * Ezc(jg+1,kg ,lg-1)
600  + owx * wy * wz * Ezc(jg ,kg+1,lg-1)
601  + wx * wy * wz * Ezc(jg+1,kg+1,lg-1);
602 
603  // interp from fine staggered to fine nodal
604  Real ef = 0.5_rt*(Ezf_zeropad(j,k,l-1) + Ezf_zeropad(j,k,l));
605 
606 #endif
607 
608  Eza(j,k,l) = eg + (ef-ec);
609 }
610 
630 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
631 void warpx_interp (const int j,
632  const int k,
633  const int l,
634  amrex::Array4<amrex::Real > const& dst_arr,
635  amrex::Array4<amrex::Real const> const& src_arr,
636  const amrex::IntVect& dst_stag,
637  const amrex::IntVect& src_stag,
638  const int nox = 2,
639  const int noy = 2,
640  const int noz = 2,
641  amrex::Real const* stencil_coeffs_x = nullptr,
642  amrex::Real const* stencil_coeffs_y = nullptr,
643  amrex::Real const* stencil_coeffs_z = nullptr)
644 {
645  using namespace amrex;
646 
647  // Pad input array with zeros beyond ghost cells
648  // for out-of-bound accesses due to large-stencil operations
649  const auto src_arr_zeropad = [src_arr] (const int jj, const int kk, const int ll) noexcept
650  {
651  return src_arr.contains(jj,kk,ll) ? src_arr(jj,kk,ll) : 0.0_rt;
652  };
653 
654  // Avoid compiler warnings
655 #if defined(WARPX_DIM_1D_Z)
656  amrex::ignore_unused(nox, noy, stencil_coeffs_x, stencil_coeffs_y);
657 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
658  amrex::ignore_unused(noy, stencil_coeffs_y);
659 #endif
660 
661  // If dst_nodal = true , we are centering from a staggered grid to a nodal grid
662  // If dst_nodal = false, we are centering from a nodal grid to a staggered grid
663  const bool dst_nodal = (dst_stag == amrex::IntVect::TheNodeVector());
664 
665  // See 1D examples below to understand the meaning of this integer shift
666  const int shift = (dst_nodal) ? 0 : 1;
667 
668  // Staggering (s = 0 if cell-centered, s = 1 if nodal)
669  const int sj = (dst_nodal) ? src_stag[0] : dst_stag[0];
670 #if (AMREX_SPACEDIM >= 2)
671  const int sk = (dst_nodal) ? src_stag[1] : dst_stag[1];
672 #endif
673 #if defined(WARPX_DIM_3D)
674  const int sl = (dst_nodal) ? src_stag[2] : dst_stag[2];
675 #endif
676 
677  // Interpolate along j,k,l only if source MultiFab is staggered along j,k,l
678  const bool interp_j = (sj == 0);
679 #if (AMREX_SPACEDIM >= 2)
680  const bool interp_k = (sk == 0);
681 #endif
682 #if defined(WARPX_DIM_3D)
683  const bool interp_l = (sl == 0);
684 #endif
685 
686 #if defined(WARPX_DIM_1D_Z)
687  const int noj = noz;
688 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
689  const int noj = nox;
690  const int nok = noz;
691 #elif defined(WARPX_DIM_3D)
692  const int noj = nox;
693  const int nok = noy;
694  const int nol = noz;
695 #endif
696 
697  // Additional normalization factor
698  const amrex::Real wj = (interp_j) ? 0.5_rt : 1.0_rt;
699 #if defined(WARPX_DIM_1D_Z)
700  constexpr amrex::Real wk = 1.0_rt;
701  constexpr amrex::Real wl = 1.0_rt;
702 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
703  const amrex::Real wk = (interp_k) ? 0.5_rt : 1.0_rt;
704  constexpr amrex::Real wl = 1.0_rt;
705 #elif defined(WARPX_DIM_3D)
706  const amrex::Real wk = (interp_k) ? 0.5_rt : 1.0_rt;
707  const amrex::Real wl = (interp_l) ? 0.5_rt : 1.0_rt;
708 #endif
709 
710  // Min and max for interpolation loop along j
711  const int jmin = (interp_j) ? j - noj/2 + shift : j;
712  const int jmax = (interp_j) ? j + noj/2 + shift - 1 : j;
713 
714  // Min and max for interpolation loop along k
715 #if defined(WARPX_DIM_1D_Z)
716  // k = 0 always
717  const int kmin = k;
718  const int kmax = k;
719 #else
720  const int kmin = (interp_k) ? k - nok/2 + shift : k;
721  const int kmax = (interp_k) ? k + nok/2 + shift - 1 : k;
722 #endif
723 
724  // Min and max for interpolation loop along l
725 #if (AMREX_SPACEDIM <= 2)
726  // l = 0 always
727  const int lmin = l;
728  const int lmax = l;
729 #elif defined(WARPX_DIM_3D)
730  const int lmin = (interp_l) ? l - nol/2 + shift : l;
731  const int lmax = (interp_l) ? l + nol/2 + shift - 1 : l;
732 #endif
733 
734  // Number of interpolation points
735  const int nj = jmax - jmin;
736  const int nk = kmax - kmin;
737  const int nl = lmax - lmin;
738 
739  // Example of 1D centering from nodal grid to nodal grid (simple copy):
740  //
741  // j
742  // --o-----o-----o-- result(j) = f(j)
743  // --o-----o-----o--
744  // j-1 j j+1
745  //
746  // Example of 1D linear centering from staggered grid to nodal grid:
747  //
748  // j
749  // --o-----o-----o-- result(j) = (f(j-1) + f(j)) / 2
750  // -----x-----x-----
751  // j-1 j
752  //
753  // Example of 1D linear centering from nodal grid to staggered grid:
754  // (note the shift of +1 in the indices with respect to the case above, see variable "shift")
755  //
756  // j
757  // --x-----x-----x-- result(j) = (f(j) + f(j+1)) / 2
758  // -----o-----o-----
759  // j j+1
760  //
761  // Example of 1D finite-order centering from staggered grid to nodal grid:
762  //
763  // j
764  // --o-----o-----o-----o-----o-----o-----o-- result(j) = c_0 * (f(j-1) + f(j) ) / 2
765  // -----x-----x-----x-----x-----x-----x----- + c_1 * (f(j-2) + f(j+1)) / 2
766  // j-3 j-2 j-1 j j+1 j+2 + c_2 * (f(j-3) + f(j+2)) / 2
767  // c_2 c_1 c_0 c_0 c_1 c_2 + ...
768  //
769  // Example of 1D finite-order centering from nodal grid to staggered grid:
770  // (note the shift of +1 in the indices with respect to the case above, see variable "shift")
771  //
772  // j
773  // --x-----x-----x-----x-----x-----x-----x-- result(j) = c_0 * (f(j) + f(j+1)) / 2
774  // -----o-----o-----o-----o-----o-----o----- + c_1 * (f(j-1) + f(j+2)) / 2
775  // j-2 j-1 j j+1 j+2 j+3 + c_2 * (f(j-2) + f(j+3)) / 2
776  // c_2 c_1 c_0 c_0 c_1 c_2 + ...
777 
778  amrex::Real res = 0.0_rt;
779 
780  amrex::Real cj = 1.0_rt;
781  amrex::Real ck = 1.0_rt;
782  amrex::Real cl = 1.0_rt;
783 
784 #if defined(WARPX_DIM_1D_Z)
785  amrex::Real const* scj = stencil_coeffs_z;
786 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
787  amrex::Real const* scj = stencil_coeffs_x;
788  amrex::Real const* sck = stencil_coeffs_z;
789 #elif defined(WARPX_DIM_3D)
790  amrex::Real const* scj = stencil_coeffs_x;
791  amrex::Real const* sck = stencil_coeffs_y;
792  amrex::Real const* scl = stencil_coeffs_z;
793 #endif
794 
795  for (int ll = 0; ll <= nl; ll++)
796  {
797 #if defined(WARPX_DIM_3D)
798  if (interp_l) cl = scl[ll];
799 #endif
800  for (int kk = 0; kk <= nk; kk++)
801  {
802 #if (AMREX_SPACEDIM >= 2)
803  if (interp_k) ck = sck[kk];
804 #endif
805  for (int jj = 0; jj <= nj; jj++)
806  {
807  if (interp_j) cj = scj[jj];
808 
809  res += cj * ck * cl * src_arr_zeropad(jmin+jj,kmin+kk,lmin+ll);
810  }
811  }
812  }
813 
814  dst_arr(j,k,l) = wj * wk * wl * res;
815 }
816 
817 #endif
amrex::shift
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Box shift(const Box &b, int dir, int nzones) noexcept
Stencil.noz
int noz
Definition: Stencil.py:472
Stencil.noy
int noy
Definition: Stencil.py:471
Stencil.nox
int nox
Definition: Stencil.py:470
warpx_interp_nd_efield_y
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void warpx_interp_nd_efield_y(int j, int k, int l, amrex::Array4< amrex::Real > const &Eya, amrex::Array4< amrex::Real const > const &Eyf, amrex::Array4< amrex::Real const > const &Eyc, amrex::Array4< amrex::Real const > const &Eyg)
Definition: WarpXComm_K.H:431
amrex::ignore_unused
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
AMREX_FORCE_INLINE
#define AMREX_FORCE_INLINE
amrex::coarsen
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Box coarsen(const Box &b, int ref_ratio) noexcept
AMREX_GPU_DEVICE
#define AMREX_GPU_DEVICE
amrex::Array4::contains
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool contains(int i, int j, int k) const noexcept
warpx_interp_nd_bfield_y
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void warpx_interp_nd_bfield_y(int j, int k, int l, amrex::Array4< amrex::Real > const &Bya, amrex::Array4< amrex::Real const > const &Byf, amrex::Array4< amrex::Real const > const &Byc, amrex::Array4< amrex::Real const > const &Byg)
Definition: WarpXComm_K.H:158
warpx_interp_nd_bfield_x
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void warpx_interp_nd_bfield_x(int j, int k, int l, amrex::Array4< amrex::Real > const &Bxa, amrex::Array4< amrex::Real const > const &Bxf, amrex::Array4< amrex::Real const > const &Bxc, amrex::Array4< amrex::Real const > const &Bxg)
Definition: WarpXComm_K.H:68
warpx_interp_nd_efield_x
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void warpx_interp_nd_efield_x(int j, int k, int l, amrex::Array4< amrex::Real > const &Exa, amrex::Array4< amrex::Real const > const &Exf, amrex::Array4< amrex::Real const > const &Exc, amrex::Array4< amrex::Real const > const &Exg)
Definition: WarpXComm_K.H:341
AMReX_FArrayBox.H
check_interp_points_and_weights.jj
jj
Definition: check_interp_points_and_weights.py:159
warpx_interp_nd_efield_z
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void warpx_interp_nd_efield_z(int j, int k, int l, amrex::Array4< amrex::Real > const &Eza, amrex::Array4< amrex::Real const > const &Ezf, amrex::Array4< amrex::Real const > const &Ezc, amrex::Array4< amrex::Real const > const &Ezg)
Definition: WarpXComm_K.H:522
warpx_interp_nd_bfield_z
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void warpx_interp_nd_bfield_z(int j, int k, int l, amrex::Array4< amrex::Real > const &Bza, amrex::Array4< amrex::Real const > const &Bzf, amrex::Array4< amrex::Real const > const &Bzc, amrex::Array4< amrex::Real const > const &Bzg)
Definition: WarpXComm_K.H:250
amrex::IntVect::TheNodeVector
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE IntVect TheNodeVector() noexcept
warpx_interp
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void warpx_interp(int j, int k, int l, amrex::Array4< amrex::Real > const &arr_aux, amrex::Array4< amrex::Real const > const &arr_fine, amrex::Array4< amrex::Real const > const &arr_coarse, const amrex::IntVect &arr_stag, const amrex::IntVect &rr)
Definition: WarpXComm_K.H:14
Generated by   doxygen 1.8.13