WarpX
WarpX_PML_kernels.H
Go to the documentation of this file.
1 /* Copyright 2019 Remi Lehe, Revathi Jambunathan, Revathi Jambunathan
2  *
3  *
4  * This file is part of WarpX.
5  *
6  * License: BSD-3-Clause-LBNL
7  */
8 #ifndef WARPX_PML_KERNELS_H_
9 #define WARPX_PML_KERNELS_H_
10 
12 #include "Utils/TextMsg.H"
13 
14 #include <AMReX.H>
15 #include <AMReX_FArrayBox.H>
16 
18 void warpx_damp_pml_ex (int i, int j, int k, amrex::Array4<amrex::Real> const& Ex,
19  const amrex::IntVect& Ex_stag,
20  const amrex::Real* const sigma_fac_x,
21  const amrex::Real* const sigma_fac_y,
22  const amrex::Real* const sigma_fac_z,
23  const amrex::Real* const sigma_star_fac_x,
24  const amrex::Real* const sigma_star_fac_y,
25  const amrex::Real* const sigma_star_fac_z,
26  int xlo, int ylo, int zlo,
27  const bool dive_cleaning)
28 {
29 #if defined(WARPX_DIM_1D_Z)
30  amrex::ignore_unused(i, j, k, Ex, Ex_stag, sigma_fac_x, sigma_fac_y, sigma_fac_z,
31  sigma_star_fac_x, sigma_star_fac_y, sigma_star_fac_z, xlo, ylo, zlo, dive_cleaning);
32  amrex::Abort("PML not implemented in Cartesian 1D geometry");
33 #endif
34 
35 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
36 
37  amrex::ignore_unused(sigma_fac_y, sigma_star_fac_y, ylo);
38 
39  // sx = 0 means that Ex is staggered in x, while sx = 1 means that Ex is nodal in x (same for z)
40  const int sx = Ex_stag[0];
41  const int sz = Ex_stag[1];
42 
43  if (dive_cleaning)
44  {
45  // Exx
46  if (sx == 0) {
47  Ex(i,j,k,PMLComp::xx) *= sigma_star_fac_x[i-xlo];
48  } else {
49  Ex(i,j,k,PMLComp::xx) *= sigma_fac_x[i-xlo];
50  }
51  }
52 
53  // Exz
54  if (sz == 0) {
55  Ex(i,j,k,PMLComp::xz) *= sigma_star_fac_z[j-zlo];
56  } else {
57  Ex(i,j,k,PMLComp::xz) *= sigma_fac_z[j-zlo];
58  }
59 
60 #elif defined(WARPX_DIM_3D)
61 
62  // sx = 0 means that Ex is staggered in x, while sx = 1 means that Ex is nodal in x (same for y, z)
63  const int sx = Ex_stag[0];
64  const int sy = Ex_stag[1];
65  const int sz = Ex_stag[2];
66 
67  if (dive_cleaning)
68  {
69  // Exx
70  if (sx == 0) {
71  Ex(i,j,k,PMLComp::xx) *= sigma_star_fac_x[i-xlo];
72  } else {
73  Ex(i,j,k,PMLComp::xx) *= sigma_fac_x[i-xlo];
74  }
75  }
76 
77  // Exy
78  if (sy == 0) {
79  Ex(i,j,k,PMLComp::xy) *= sigma_star_fac_y[j-ylo];
80  } else {
81  Ex(i,j,k,PMLComp::xy) *= sigma_fac_y[j-ylo];
82  }
83 
84  // Exz
85  if (sz == 0) {
86  Ex(i,j,k,PMLComp::xz) *= sigma_star_fac_z[k-zlo];
87  } else {
88  Ex(i,j,k,PMLComp::xz) *= sigma_fac_z[k-zlo];
89  }
90 
91 #endif
92 }
93 
95 void warpx_damp_pml_ey (int i, int j, int k, amrex::Array4<amrex::Real> const& Ey,
96  const amrex::IntVect& Ey_stag,
97  const amrex::Real* const sigma_fac_x,
98  const amrex::Real* const sigma_fac_y,
99  const amrex::Real* const sigma_fac_z,
100  const amrex::Real* const sigma_star_fac_x,
101  const amrex::Real* const sigma_star_fac_y,
102  const amrex::Real* const sigma_star_fac_z,
103  int xlo, int ylo, int zlo,
104  const bool dive_cleaning)
105 {
106 #if defined(WARPX_DIM_1D_Z)
107  amrex::ignore_unused(i, j, k, Ey, Ey_stag, sigma_fac_x, sigma_fac_y, sigma_fac_z,
108  sigma_star_fac_x, sigma_star_fac_y, sigma_star_fac_z, xlo, ylo, zlo, dive_cleaning);
109  amrex::Abort("PML not implemented in Cartesian 1D geometry");
110 #endif
111 
112 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
113 
114  amrex::ignore_unused(sigma_fac_y, sigma_star_fac_y, ylo, dive_cleaning);
115 
116  // sx = 0 means that Ey is staggered in x, while sx = 1 means that Ey is nodal in x (same for z)
117  const int sx = Ey_stag[0];
118  const int sz = Ey_stag[1];
119 
120  // Eyx
121  if (sx == 0) {
122  Ey(i,j,k,PMLComp::yx) *= sigma_star_fac_x[i-xlo];
123  } else {
124  Ey(i,j,k,PMLComp::yx) *= sigma_fac_x[i-xlo];
125  }
126 
127  // Eyz
128  if (sz == 0) {
129  Ey(i,j,k,PMLComp::yz) *= sigma_star_fac_z[j-zlo];
130  } else {
131  Ey(i,j,k,PMLComp::yz) *= sigma_fac_z[j-zlo];
132  }
133 
134 #elif defined(WARPX_DIM_3D)
135 
136  // sx = 0 means that Ey is staggered in x, while sx = 1 means that Ey is nodal in x (same for y, z)
137  const int sx = Ey_stag[0];
138  const int sy = Ey_stag[1];
139  const int sz = Ey_stag[2];
140 
141  // Eyx
142  if (sx == 0) {
143  Ey(i,j,k,PMLComp::yx) *= sigma_star_fac_x[i-xlo];
144  } else {
145  Ey(i,j,k,PMLComp::yx) *= sigma_fac_x[i-xlo];
146  }
147 
148  if (dive_cleaning)
149  {
150  // Eyy
151  if (sy == 0) {
152  Ey(i,j,k,PMLComp::yy) *= sigma_star_fac_y[j-ylo];
153  } else {
154  Ey(i,j,k,PMLComp::yy) *= sigma_fac_y[j-ylo];
155  }
156  }
157 
158  // Eyz
159  if (sz == 0) {
160  Ey(i,j,k,PMLComp::yz) *= sigma_star_fac_z[k-zlo];
161  } else {
162  Ey(i,j,k,PMLComp::yz) *= sigma_fac_z[k-zlo];
163  }
164 
165 #endif
166 }
167 
169 void warpx_damp_pml_ez (int i, int j, int k, amrex::Array4<amrex::Real> const& Ez,
170  const amrex::IntVect& Ez_stag,
171  const amrex::Real* const sigma_fac_x,
172  const amrex::Real* const sigma_fac_y,
173  const amrex::Real* const sigma_fac_z,
174  const amrex::Real* const sigma_star_fac_x,
175  const amrex::Real* const sigma_star_fac_y,
176  const amrex::Real* const sigma_star_fac_z,
177  int xlo, int ylo, int zlo,
178  const bool dive_cleaning)
179 {
180 #if defined(WARPX_DIM_1D_Z)
181  amrex::ignore_unused(i, j, k, Ez, Ez_stag, sigma_fac_x, sigma_fac_y, sigma_fac_z,
182  sigma_star_fac_x, sigma_star_fac_y, sigma_star_fac_z, xlo, ylo, zlo, dive_cleaning);
183  amrex::Abort("PML not implemented in Cartesian 1D geometry");
184 #endif
185 
186 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
187 
188  amrex::ignore_unused(sigma_fac_y, sigma_star_fac_y, ylo);
189 
190  // sx = 0 means that Ez is staggered in x, while sx = 1 means that Ez is nodal in x (same for z)
191  const int sx = Ez_stag[0];
192  const int sz = Ez_stag[1];
193 
194  // Ezx
195  if (sx == 0) {
196  Ez(i,j,k,PMLComp::zx) *= sigma_star_fac_x[i-xlo];
197  } else {
198  Ez(i,j,k,PMLComp::zx) *= sigma_fac_x[i-xlo];
199  }
200 
201  if (dive_cleaning)
202  {
203  // Ezz
204  if (sz == 0) {
205  Ez(i,j,k,PMLComp::zz) *= sigma_star_fac_z[j-zlo];
206  } else {
207  Ez(i,j,k,PMLComp::zz) *= sigma_fac_z[j-zlo];
208  }
209  }
210 
211 #elif defined(WARPX_DIM_3D)
212 
213  // sx = 0 means that Ez is staggered in x, while sx = 1 means that Ez is nodal in x (same for y, z)
214  const int sx = Ez_stag[0];
215  const int sy = Ez_stag[1];
216  const int sz = Ez_stag[2];
217 
218  // Ezx
219  if (sx == 0) {
220  Ez(i,j,k,PMLComp::zx) *= sigma_star_fac_x[i-xlo];
221  } else {
222  Ez(i,j,k,PMLComp::zx) *= sigma_fac_x[i-xlo];
223  }
224 
225  // Ezy
226  if (sy == 0) {
227  Ez(i,j,k,PMLComp::zy) *= sigma_star_fac_y[j-ylo];
228  } else {
229  Ez(i,j,k,PMLComp::zy) *= sigma_fac_y[j-ylo];
230  }
231 
232  if (dive_cleaning)
233  {
234  // Ezz
235  if (sz == 0) {
236  Ez(i,j,k,PMLComp::zz) *= sigma_star_fac_z[k-zlo];
237  } else {
238  Ez(i,j,k,PMLComp::zz) *= sigma_fac_z[k-zlo];
239  }
240  }
241 
242 #endif
243 }
244 
246 void warpx_damp_pml_bx (int i, int j, int k, amrex::Array4<amrex::Real> const& Bx,
247  const amrex::IntVect& Bx_stag,
248  const amrex::Real* const sigma_fac_x,
249  const amrex::Real* const sigma_fac_y,
250  const amrex::Real* const sigma_fac_z,
251  const amrex::Real* const sigma_star_fac_x,
252  const amrex::Real* const sigma_star_fac_y,
253  const amrex::Real* const sigma_star_fac_z,
254  int xlo, int ylo, int zlo,
255  const bool divb_cleaning)
256 {
257 #if defined(WARPX_DIM_1D_Z)
258  amrex::ignore_unused(i, j, k, Bx, Bx_stag, sigma_fac_x, sigma_fac_y, sigma_fac_z,
259  sigma_star_fac_x, sigma_star_fac_y, sigma_star_fac_z, xlo, ylo, zlo, divb_cleaning);
260  amrex::Abort("PML not implemented in Cartesian 1D geometry");
261 #endif
262 
263 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
264 
265  amrex::ignore_unused(sigma_fac_y, sigma_star_fac_y, ylo);
266 
267  // sx = 0 means that Bx is staggered in x, while sx = 1 means that Bx is nodal in x (same for z)
268  const int sx = Bx_stag[0];
269  const int sz = Bx_stag[1];
270 
271  if (divb_cleaning)
272  {
273  // Bxx
274  if (sx == 0) {
275  Bx(i,j,k,PMLComp::xx) *= sigma_star_fac_x[i-xlo];
276  } else {
277  Bx(i,j,k,PMLComp::xx) *= sigma_fac_x[i-xlo];
278  }
279  }
280 
281  // Bxz
282  if (sz == 0) {
283  Bx(i,j,k,PMLComp::xz) *= sigma_star_fac_z[j-zlo];
284  } else {
285  Bx(i,j,k,PMLComp::xz) *= sigma_fac_z[j-zlo];
286  }
287 
288 #elif defined(WARPX_DIM_3D)
289 
290  // sx = 0 means that Bx is staggered in x, while sx = 1 means that Bx is nodal in x (same for y, z)
291  const int sx = Bx_stag[0];
292  const int sy = Bx_stag[1];
293  const int sz = Bx_stag[2];
294 
295  if (divb_cleaning)
296  {
297  // Bxx
298  if (sx == 0) {
299  Bx(i,j,k,PMLComp::xx) *= sigma_star_fac_x[i-xlo];
300  } else {
301  Bx(i,j,k,PMLComp::xx) *= sigma_fac_x[i-xlo];
302  }
303  }
304 
305  // Bxy
306  if (sy == 0) {
307  Bx(i,j,k,PMLComp::xy) *= sigma_star_fac_y[j-ylo];
308  } else {
309  Bx(i,j,k,PMLComp::xy) *= sigma_fac_y[j-ylo];
310  }
311 
312  // Bxz
313  if (sz == 0) {
314  Bx(i,j,k,PMLComp::xz) *= sigma_star_fac_z[k-zlo];
315  } else {
316  Bx(i,j,k,PMLComp::xz) *= sigma_fac_z[k-zlo];
317  }
318 
319 #endif
320 }
321 
323 void warpx_damp_pml_by (int i, int j, int k, amrex::Array4<amrex::Real> const& By,
324  const amrex::IntVect& By_stag,
325  const amrex::Real* const sigma_fac_x,
326  const amrex::Real* const sigma_fac_y,
327  const amrex::Real* const sigma_fac_z,
328  const amrex::Real* const sigma_star_fac_x,
329  const amrex::Real* const sigma_star_fac_y,
330  const amrex::Real* const sigma_star_fac_z,
331  int xlo, int ylo, int zlo,
332  const bool divb_cleaning)
333 {
334 #if defined(WARPX_DIM_1D_Z)
335  amrex::ignore_unused(i, j, k, By, By_stag, sigma_fac_x, sigma_fac_y, sigma_fac_z,
336  sigma_star_fac_x, sigma_star_fac_y, sigma_star_fac_z, xlo, ylo, zlo, divb_cleaning);
337  amrex::Abort("PML not implemented in Cartesian 1D geometry");
338 #endif
339 
340 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
341 
342  amrex::ignore_unused(sigma_fac_y, sigma_star_fac_y, ylo, divb_cleaning);
343 
344  // sx = 0 means that By is staggered in x, while sx = 1 means that By is nodal in x (same for z)
345  const int sx = By_stag[0];
346  const int sz = By_stag[1];
347 
348  // Byx
349  if (sx == 0) {
350  By(i,j,k,PMLComp::yx) *= sigma_star_fac_x[i-xlo];
351  } else {
352  By(i,j,k,PMLComp::yx) *= sigma_fac_x[i-xlo];
353  }
354 
355  // Byz
356  if (sz == 0) {
357  By(i,j,k,PMLComp::yz) *= sigma_star_fac_z[j-zlo];
358  } else {
359  By(i,j,k,PMLComp::yz) *= sigma_fac_z[j-zlo];
360  }
361 
362 #elif defined(WARPX_DIM_3D)
363 
364  // sx = 0 means that By is staggered in x, while sx = 1 means that By is nodal in x (same for y, z)
365  const int sx = By_stag[0];
366  const int sy = By_stag[1];
367  const int sz = By_stag[2];
368 
369  // Byx
370  if (sx == 0) {
371  By(i,j,k,PMLComp::yx) *= sigma_star_fac_x[i-xlo];
372  } else {
373  By(i,j,k,PMLComp::yx) *= sigma_fac_x[i-xlo];
374  }
375 
376  if (divb_cleaning)
377  {
378  // Byy
379  if (sy == 0) {
380  By(i,j,k,PMLComp::yy) *= sigma_star_fac_y[j-ylo];
381  } else {
382  By(i,j,k,PMLComp::yy) *= sigma_fac_y[j-ylo];
383  }
384  }
385 
386  // Byz
387  if (sz == 0) {
388  By(i,j,k,PMLComp::yz) *= sigma_star_fac_z[k-zlo];
389  } else {
390  By(i,j,k,PMLComp::yz) *= sigma_fac_z[k-zlo];
391  }
392 
393 #endif
394 }
395 
397 void warpx_damp_pml_bz (int i, int j, int k, amrex::Array4<amrex::Real> const& Bz,
398  const amrex::IntVect& Bz_stag,
399  const amrex::Real* const sigma_fac_x,
400  const amrex::Real* const sigma_fac_y,
401  const amrex::Real* const sigma_fac_z,
402  const amrex::Real* const sigma_star_fac_x,
403  const amrex::Real* const sigma_star_fac_y,
404  const amrex::Real* const sigma_star_fac_z,
405  int xlo, int ylo, int zlo,
406  const bool divb_cleaning)
407 {
408 #if defined(WARPX_DIM_1D_Z)
409  amrex::ignore_unused(i, j, k, Bz, Bz_stag, sigma_fac_x, sigma_fac_y, sigma_fac_z,
410  sigma_star_fac_x, sigma_star_fac_y, sigma_star_fac_z, xlo, ylo, zlo, divb_cleaning);
411  amrex::Abort("PML not implemented in Cartesian 1D geometry");
412 #endif
413 
414 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
415 
416  amrex::ignore_unused(sigma_fac_y, sigma_star_fac_y, ylo);
417 
418  // sx = 0 means that Bz is staggered in x, while sx = 1 means that Bz is nodal in x (same for z)
419  const int sx = Bz_stag[0];
420  const int sz = Bz_stag[1];
421 
422  // Bzx
423  if (sx == 0) {
424  Bz(i,j,k,PMLComp::zx) *= sigma_star_fac_x[i-xlo];
425  } else {
426  Bz(i,j,k,PMLComp::zx) *= sigma_fac_x[i-xlo];
427  }
428 
429  if (divb_cleaning)
430  {
431  // Bzz
432  if (sz == 0) {
433  Bz(i,j,k,PMLComp::zz) *= sigma_star_fac_z[j-zlo];
434  } else {
435  Bz(i,j,k,PMLComp::zz) *= sigma_fac_z[j-zlo];
436  }
437  }
438 
439 #elif defined(WARPX_DIM_3D)
440 
441  // sx = 0 means that Bz is staggered in x, while sx = 1 means that Bz is nodal in x (same for y, z)
442  const int sx = Bz_stag[0];
443  const int sy = Bz_stag[1];
444  const int sz = Bz_stag[2];
445 
446  // Bzx
447  if (sx == 0) {
448  Bz(i,j,k,PMLComp::zx) *= sigma_star_fac_x[i-xlo];
449  } else {
450  Bz(i,j,k,PMLComp::zx) *= sigma_fac_x[i-xlo];
451  }
452 
453  // Bzy
454  if (sy == 0) {
455  Bz(i,j,k,PMLComp::zy) *= sigma_star_fac_y[j-ylo];
456  } else {
457  Bz(i,j,k,PMLComp::zy) *= sigma_fac_y[j-ylo];
458  }
459 
460  if (divb_cleaning)
461  {
462  // Bzz
463  if (sz == 0) {
464  Bz(i,j,k,PMLComp::zz) *= sigma_star_fac_z[k-zlo];
465  } else {
466  Bz(i,j,k,PMLComp::zz) *= sigma_fac_z[k-zlo];
467  }
468  }
469 
470 #endif
471 }
472 
474 void warpx_damp_pml_scalar (int i, int j, int k, amrex::Array4<amrex::Real> const& arr,
475  const amrex::IntVect& arr_stag,
476  const amrex::Real* const sigma_fac_x,
477  const amrex::Real* const sigma_fac_y,
478  const amrex::Real* const sigma_fac_z,
479  const amrex::Real* const sigma_star_fac_x,
480  const amrex::Real* const sigma_star_fac_y,
481  const amrex::Real* const sigma_star_fac_z,
482  int xlo, int ylo, int zlo)
483 {
484 #if defined(WARPX_DIM_1D_Z)
485  amrex::ignore_unused(i, j, k, arr, arr_stag, sigma_fac_x, sigma_fac_y, sigma_fac_z,
486  sigma_star_fac_x, sigma_star_fac_y, sigma_star_fac_z, xlo, ylo, zlo);
487  amrex::Abort("PML not implemented in Cartesian 1D geometry");
488 #endif
489 
490 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
491 
492  amrex::ignore_unused(sigma_fac_y, sigma_star_fac_y, ylo);
493 
494  // sx = 0 means that arr is staggered in x, while sx = 1 means that arr is nodal in x (same for z)
495  const int sx = arr_stag[0];
496  const int sz = arr_stag[1];
497 
498  // Component along x
499  if (sx == 0) {
500  arr(i,j,k,PMLComp::x) *= sigma_star_fac_x[i-xlo];
501  } else {
502  arr(i,j,k,PMLComp::x) *= sigma_fac_x[i-xlo];
503  }
504 
505  // Component along z
506  if (sz == 0) {
507  arr(i,j,k,PMLComp::z) *= sigma_star_fac_z[j-zlo];
508  } else {
509  arr(i,j,k,PMLComp::z) *= sigma_fac_z[j-zlo];
510  }
511 
512 #elif defined(WARPX_DIM_3D)
513 
514  // sx = 0 means that arr is staggered in x, while sx = 1 means that arr is nodal in x (same for y, z)
515  const int sx = arr_stag[0];
516  const int sy = arr_stag[1];
517  const int sz = arr_stag[2];
518 
519  // Component along x
520  if (sx == 0) {
521  arr(i,j,k,PMLComp::x) *= sigma_star_fac_x[i-xlo];
522  } else {
523  arr(i,j,k,PMLComp::x) *= sigma_fac_x[i-xlo];
524  }
525 
526  // Component along y
527  if (sy == 0) {
528  arr(i,j,k,PMLComp::y) *= sigma_star_fac_y[j-ylo];
529  } else {
530  arr(i,j,k,PMLComp::y) *= sigma_fac_y[j-ylo];
531  }
532 
533  // Component along z
534  if (sz == 0) {
535  arr(i,j,k,PMLComp::z) *= sigma_star_fac_z[k-zlo];
536  } else {
537  arr(i,j,k,PMLComp::z) *= sigma_fac_z[k-zlo];
538  }
539 
540 #endif
541 }
542 
543 #endif
#define AMREX_FORCE_INLINE
#define AMREX_GPU_HOST_DEVICE
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void warpx_damp_pml_ex(int i, int j, int k, amrex::Array4< amrex::Real > const &Ex, const amrex::IntVect &Ex_stag, const amrex::Real *const sigma_fac_x, const amrex::Real *const sigma_fac_y, const amrex::Real *const sigma_fac_z, const amrex::Real *const sigma_star_fac_x, const amrex::Real *const sigma_star_fac_y, const amrex::Real *const sigma_star_fac_z, int xlo, int ylo, int zlo, const bool dive_cleaning)
Definition: WarpX_PML_kernels.H:18
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void warpx_damp_pml_scalar(int i, int j, int k, amrex::Array4< amrex::Real > const &arr, const amrex::IntVect &arr_stag, const amrex::Real *const sigma_fac_x, const amrex::Real *const sigma_fac_y, const amrex::Real *const sigma_fac_z, const amrex::Real *const sigma_star_fac_x, const amrex::Real *const sigma_star_fac_y, const amrex::Real *const sigma_star_fac_z, int xlo, int ylo, int zlo)
Definition: WarpX_PML_kernels.H:474
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void warpx_damp_pml_ey(int i, int j, int k, amrex::Array4< amrex::Real > const &Ey, const amrex::IntVect &Ey_stag, const amrex::Real *const sigma_fac_x, const amrex::Real *const sigma_fac_y, const amrex::Real *const sigma_fac_z, const amrex::Real *const sigma_star_fac_x, const amrex::Real *const sigma_star_fac_y, const amrex::Real *const sigma_star_fac_z, int xlo, int ylo, int zlo, const bool dive_cleaning)
Definition: WarpX_PML_kernels.H:95
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void warpx_damp_pml_bz(int i, int j, int k, amrex::Array4< amrex::Real > const &Bz, const amrex::IntVect &Bz_stag, const amrex::Real *const sigma_fac_x, const amrex::Real *const sigma_fac_y, const amrex::Real *const sigma_fac_z, const amrex::Real *const sigma_star_fac_x, const amrex::Real *const sigma_star_fac_y, const amrex::Real *const sigma_star_fac_z, int xlo, int ylo, int zlo, const bool divb_cleaning)
Definition: WarpX_PML_kernels.H:397
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void warpx_damp_pml_by(int i, int j, int k, amrex::Array4< amrex::Real > const &By, const amrex::IntVect &By_stag, const amrex::Real *const sigma_fac_x, const amrex::Real *const sigma_fac_y, const amrex::Real *const sigma_fac_z, const amrex::Real *const sigma_star_fac_x, const amrex::Real *const sigma_star_fac_y, const amrex::Real *const sigma_star_fac_z, int xlo, int ylo, int zlo, const bool divb_cleaning)
Definition: WarpX_PML_kernels.H:323
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void warpx_damp_pml_bx(int i, int j, int k, amrex::Array4< amrex::Real > const &Bx, const amrex::IntVect &Bx_stag, const amrex::Real *const sigma_fac_x, const amrex::Real *const sigma_fac_y, const amrex::Real *const sigma_fac_z, const amrex::Real *const sigma_star_fac_x, const amrex::Real *const sigma_star_fac_y, const amrex::Real *const sigma_star_fac_z, int xlo, int ylo, int zlo, const bool divb_cleaning)
Definition: WarpX_PML_kernels.H:246
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void warpx_damp_pml_ez(int i, int j, int k, amrex::Array4< amrex::Real > const &Ez, const amrex::IntVect &Ez_stag, const amrex::Real *const sigma_fac_x, const amrex::Real *const sigma_fac_y, const amrex::Real *const sigma_fac_z, const amrex::Real *const sigma_star_fac_x, const amrex::Real *const sigma_star_fac_y, const amrex::Real *const sigma_star_fac_z, int xlo, int ylo, int zlo, const bool dive_cleaning)
Definition: WarpX_PML_kernels.H:169
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
void Abort(const std::string &msg)
i
Definition: check_interp_points_and_weights.py:174
@ zz
Definition: PMLComponent.H:18
@ xz
Definition: PMLComponent.H:16
@ z
Definition: PMLComponent.H:19
@ yy
Definition: PMLComponent.H:17
@ zx
Definition: PMLComponent.H:18
@ zy
Definition: PMLComponent.H:18
@ yz
Definition: PMLComponent.H:17
@ x
Definition: PMLComponent.H:19
@ y
Definition: PMLComponent.H:19
@ xx
Definition: PMLComponent.H:16
@ yx
Definition: PMLComponent.H:17
@ xy
Definition: PMLComponent.H:16