1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
|
/* Copyright 2019-2021 Maxence Thevenet, Michael Rowan, Luca Fedeli, Axel Huebl
*
* This file is part of WarpX.
*
* License: BSD-3-Clause-LBNL
*/
#ifndef SHAPEFACTORS_H_
#define SHAPEFACTORS_H_
#include <AMReX.H>
#include <AMReX_GpuQualifiers.H>
/**
* Compute shape factor and return index of leftmost cell where
* particle writes.
* Specializations are defined for orders 0 to 3 (using "if constexpr").
* Shape factor functors may be evaluated with double arguments
* in current deposition to ensure that current deposited by
* particles that move only a small distance is still resolved.
* Without this safeguard, single and double precision versions
* can give disagreeing results in the time evolution for some
* problem setups.
*/
template <int depos_order>
struct Compute_shape_factor
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(
T* const sx,
T xmid) const
{
if constexpr (depos_order == 0){
const auto j = static_cast<int>(xmid + T(0.5));
sx[0] = T(1.0);
return j;
}
else if constexpr (depos_order == 1){
const auto j = static_cast<int>(xmid);
const T xint = xmid - T(j);
sx[0] = T(1.0) - xint;
sx[1] = xint;
return j;
}
else if constexpr (depos_order == 2){
const auto j = static_cast<int>(xmid + T(0.5));
const T xint = xmid - T(j);
sx[0] = T(0.5)*(T(0.5) - xint)*(T(0.5) - xint);
sx[1] = T(0.75) - xint*xint;
sx[2] = T(0.5)*(T(0.5) + xint)*(T(0.5) + xint);
// index of the leftmost cell where particle deposits
return j-1;
}
else if constexpr (depos_order == 3){
const auto j = static_cast<int>(xmid);
const T xint = xmid - T(j);
sx[0] = (T(1.0))/(T(6.0))*(T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - xint);
sx[1] = (T(2.0))/(T(3.0)) - xint*xint*(T(1.0) - xint/(T(2.0)));
sx[2] = (T(2.0))/(T(3.0)) - (T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - T(0.5)*(T(1.0) - xint));
sx[3] = (T(1.0))/(T(6.0))*xint*xint*xint;
// index of the leftmost cell where particle deposits
return j-1;
}
else{
amrex::Abort("Unknown particle shape selected in Compute_shape_factor");
amrex::ignore_unused(sx, xmid);
}
return 0;
}
};
/**
* Compute shifted shape factor and return index of leftmost cell where
* particle writes, for Esirkepov algorithm.
* Specializations are defined below for orders 1, 2 and 3 (using "if constexpr").
*/
template <int depos_order>
struct Compute_shifted_shape_factor
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(
T* const sx,
const T x_old,
const int i_new) const
{
if constexpr (depos_order == 1){
const auto i = static_cast<int>(x_old);
const int i_shift = i - i_new;
const T xint = x_old - T(i);
sx[1+i_shift] = T(1.0) - xint;
sx[2+i_shift] = xint;
return i;
}
else if constexpr (depos_order == 2){
const auto i = static_cast<int>(x_old + T(0.5));
const int i_shift = i - (i_new + 1);
const T xint = x_old - T(i);
sx[1+i_shift] = T(0.5)*(T(0.5) - xint)*(T(0.5) - xint);
sx[2+i_shift] = T(0.75) - xint*xint;
sx[3+i_shift] = T(0.5)*(T(0.5) + xint)*(T(0.5) + xint);
// index of the leftmost cell where particle deposits
return i - 1;
}
else if constexpr (depos_order == 3){
const auto i = static_cast<int>(x_old);
const int i_shift = i - (i_new + 1);
const T xint = x_old - i;
sx[1+i_shift] = (T(1.0))/(T(6.0))*(T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - xint);
sx[2+i_shift] = (T(2.0))/(T(3.0)) - xint*xint*(T(1.0) - xint/(T(2.0)));
sx[3+i_shift] = (T(2.0))/(T(3.0)) - (T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - T(0.5)*(T(1.0) - xint));
sx[4+i_shift] = (T(1.0))/(T(6.0))*xint*xint*xint;
// index of the leftmost cell where particle deposits
return i - 1;
}
else{
amrex::Abort("Unknown particle shape selected in Compute_shifted_shape_factor");
amrex::ignore_unused(sx, x_old, i_new);
}
return 0;
}
};
#endif // SHAPEFACTORS_H_
|
