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