WarpX
Loading...
Searching...
No Matches
ShapeFactors.H
Go to the documentation of this file.
1/* Copyright 2019-2021 Maxence Thevenet, Michael Rowan, Luca Fedeli, Axel Huebl
2 *
3 * This file is part of WarpX.
4 *
5 * License: BSD-3-Clause-LBNL
6 */
7#ifndef WARPX_SHAPEFACTORS_H_
8#define WARPX_SHAPEFACTORS_H_
9
10#include "Utils/TextMsg.H"
11
12#include <AMReX.H>
13#include <AMReX_GpuQualifiers.H>
14
15
27template <int depos_order>
29{
30 template< typename T >
33 T* const sx,
34 T xmid) const
35 {
36 if constexpr (depos_order == 0){
37 const auto j = static_cast<int>(xmid + T(0.5));
38 sx[0] = T(1.0);
39 return j;
40 }
41 else if constexpr (depos_order == 1){
42 const auto j = static_cast<int>(xmid);
43 const T xint = xmid - T(j);
44 sx[0] = T(1.0) - xint;
45 sx[1] = xint;
46 return j;
47 }
48 else if constexpr (depos_order == 2){
49 const auto j = static_cast<int>(xmid + T(0.5));
50 const T xint = xmid - T(j);
51 sx[0] = T(0.5)*(T(0.5) - xint)*(T(0.5) - xint);
52 sx[1] = T(0.75) - xint*xint;
53 sx[2] = T(0.5)*(T(0.5) + xint)*(T(0.5) + xint);
54 // index of the leftmost cell where particle deposits
55 return j-1;
56 }
57 else if constexpr (depos_order == 3){
58 const auto j = static_cast<int>(xmid);
59 const T xint = xmid - T(j);
60 sx[0] = (T(1.0))/(T(6.0))*(T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - xint);
61 sx[1] = (T(2.0))/(T(3.0)) - xint*xint*(T(1.0) - xint/(T(2.0)));
62 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));
63 sx[3] = (T(1.0))/(T(6.0))*xint*xint*xint;
64 // index of the leftmost cell where particle deposits
65 return j-1;
66 }
67 else if constexpr (depos_order == 4){
68 const auto j = static_cast<int>(xmid + T(0.5));
69 const T xint = xmid - T(j);
70 sx[0] = (T(1.0))/(T(24.0))*(T(0.5) - xint)*(T(0.5) - xint)*(T(0.5) - xint)*(T(0.5) - xint);
71 sx[1] = (T(1.0))/(T(24.0))*(T(4.75) - T(11.0)*xint + T(4.0)*xint*xint*(T(1.5) + xint - xint*xint));
72 sx[2] = (T(1.0))/(T(24.0))*(T(14.375) + T(6.0)*xint*xint*(xint*xint - T(2.5)));
73 sx[3] = (T(1.0))/(T(24.0))*(T(4.75) + T(11.0)*xint + T(4.0)*xint*xint*(T(1.5) - xint - xint*xint));
74 sx[4] = (T(1.0))/(T(24.0))*(T(0.5) + xint)*(T(0.5) + xint)*(T(0.5) + xint)*(T(0.5)+xint);
75 // index of the leftmost cell where particle deposits
76 return j-2;
77 }
78 else{
79 WARPX_ABORT_WITH_MESSAGE("Unknown particle shape selected in Compute_shape_factor");
80 amrex::ignore_unused(sx, xmid);
81 }
82 return 0;
83 }
84};
85
86
87
93template <int depos_order>
95{
96 template< typename T >
99 T* const sx,
100 const T x_old,
101 const int i_new) const
102 {
103 if constexpr (depos_order == 0){
104 const auto i = static_cast<int>(std::floor(x_old + T(0.5)));
105 const int i_shift = i - i_new;
106 sx[1+i_shift] = T(1.0);
107 return i;
108 }
109 else if constexpr (depos_order == 1){
110 const auto i = static_cast<int>(std::floor(x_old));
111 const int i_shift = i - i_new;
112 const T xint = x_old - T(i);
113 sx[1+i_shift] = T(1.0) - xint;
114 sx[2+i_shift] = xint;
115 return i;
116 }
117 else if constexpr (depos_order == 2){
118 const auto i = static_cast<int>(x_old + T(0.5));
119 const int i_shift = i - (i_new + 1);
120 const T xint = x_old - T(i);
121 sx[1+i_shift] = T(0.5)*(T(0.5) - xint)*(T(0.5) - xint);
122 sx[2+i_shift] = T(0.75) - xint*xint;
123 sx[3+i_shift] = T(0.5)*(T(0.5) + xint)*(T(0.5) + xint);
124 // index of the leftmost cell where particle deposits
125 return i - 1;
126 }
127 else if constexpr (depos_order == 3){
128 const auto i = static_cast<int>(x_old);
129 const int i_shift = i - (i_new + 1);
130 const T xint = x_old - T(i);
131 sx[1+i_shift] = (T(1.0))/(T(6.0))*(T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - xint);
132 sx[2+i_shift] = (T(2.0))/(T(3.0)) - xint*xint*(T(1.0) - xint/(T(2.0)));
133 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));
134 sx[4+i_shift] = (T(1.0))/(T(6.0))*xint*xint*xint;
135 // index of the leftmost cell where particle deposits
136 return i - 1;
137 }
138 else if constexpr (depos_order == 4){
139 const auto i = static_cast<int>(x_old + T(0.5));
140 const int i_shift = i - (i_new + 2);
141 const T xint = x_old - T(i);
142 sx[1+i_shift] = (T(1.0))/(T(24.0))*(T(0.5) - xint)*(T(0.5) - xint)*(T(0.5) - xint)*(T(0.5) - xint);
143 sx[2+i_shift] = (T(1.0))/(T(24.0))*(T(4.75) - T(11.0)*xint + T(4.0)*xint*xint*(T(1.5) + xint - xint*xint));
144 sx[3+i_shift] = (T(1.0))/(T(24.0))*(T(14.375) + T(6.0)*xint*xint*(xint*xint - T(2.5)));
145 sx[4+i_shift] = (T(1.0))/(T(24.0))*(T(4.75) + T(11.0)*xint + T(4.0)*xint*xint*(T(1.5) - xint - xint*xint));
146 sx[5+i_shift] = (T(1.0))/(T(24.0))*(T(0.5) + xint)*(T(0.5) + xint)*(T(0.5) + xint)*(T(0.5)+xint);
147 // index of the leftmost cell where particle deposits
148 return i - 2;
149 }
150 else{
151 WARPX_ABORT_WITH_MESSAGE("Unknown particle shape selected in Compute_shifted_shape_factor");
152 amrex::ignore_unused(sx, x_old, i_new);
153 }
154 return 0;
155 }
156};
157
166template <int depos_order>
168{
169 template< typename T >
172 T* const sx_old,
173 T* const sx_new,
174 T xold,
175 T xnew) const
176 {
177 const T xmid = T(0.5)*(xnew + xold);
178 if constexpr (depos_order == 1){
179 const auto j = static_cast<int>(xmid);
180 const T xint_old = xold - T(j);
181 sx_old[0] = T(1.0) - xint_old;
182 sx_old[1] = xint_old;
183 //
184 const T xint_new = xnew - T(j);
185 sx_new[0] = T(1.0) - xint_new;
186 sx_new[1] = xint_new;
187 return j;
188 }
189 else if constexpr (depos_order == 2){
190 const auto j = static_cast<int>(xmid + T(0.5));
191 const T xint_old = xold - T(j);
192 sx_old[0] = T(0.5)*(T(0.5) - xint_old)*(T(0.5) - xint_old);
193 sx_old[1] = T(0.75) - xint_old*xint_old;
194 sx_old[2] = T(0.5)*(T(0.5) + xint_old)*(T(0.5) + xint_old);
195 //
196 const T xint_new = xnew - T(j);
197 sx_new[0] = T(0.5)*(T(0.5) - xint_new)*(T(0.5) - xint_new);
198 sx_new[1] = T(0.75) - xint_new*xint_new;
199 sx_new[2] = T(0.5)*(T(0.5) + xint_new)*(T(0.5) + xint_new);
200 // index of the leftmost cell where particle deposits
201 return j-1;
202 }
203 else if constexpr (depos_order == 3){
204 const auto j = static_cast<int>(xmid);
205 const T xint_old = xold - T(j);
206 sx_old[0] = T(1.0)/T(6.0)*(T(1.0) - xint_old)*(T(1.0) - xint_old)*(T(1.0) - xint_old);
207 sx_old[1] = T(2.0)/T(3.0) - xint_old*xint_old*(T(1.0) - xint_old/(T(2.0)));
208 sx_old[2] = T(2.0)/T(3.0) - (T(1.0) - xint_old)*(T(1.0) - xint_old)*(T(1.0) - T(0.5)*(T(1.0) - xint_old));
209 sx_old[3] = T(1.0)/T(6.0)*xint_old*xint_old*xint_old;
210 //
211 const T xint_new = xnew - T(j);
212 sx_new[0] = T(1.0)/T(6.0)*(T(1.0) - xint_new)*(T(1.0) - xint_new)*(T(1.0) - xint_new);
213 sx_new[1] = T(2.0)/T(3.0) - xint_new*xint_new*(T(1.0) - xint_new/(T(2.0)));
214 sx_new[2] = T(2.0)/T(3.0) - (T(1.0) - xint_new)*(T(1.0) - xint_new)*(T(1.0) - T(0.5)*(T(1.0) - xint_new));
215 sx_new[3] = T(1.0)/T(6.0)*xint_new*xint_new*xint_new;
216 // index of the leftmost cell where particle deposits
217 return j-1;
218 }
219 else if constexpr (depos_order == 4){
220 const auto j = static_cast<int>(xmid + T(0.5));
221 const T xint_old = xold - T(j);
222 sx_old[0] = (T(1.0))/(T(24.0))*(T(0.5) - xint_old)*(T(0.5) - xint_old)*(T(0.5) - xint_old)*(T(0.5) - xint_old);
223 sx_old[1] = (T(1.0))/(T(24.0))*(T(4.75) - T(11.0)*xint_old + T(4.0)*xint_old*xint_old*(T(1.5) + xint_old - xint_old*xint_old));
224 sx_old[2] = (T(1.0))/(T(24.0))*(T(14.375) + T(6.0)*xint_old*xint_old*(xint_old*xint_old - T(2.5)));
225 sx_old[3] = (T(1.0))/(T(24.0))*(T(4.75) + T(11.0)*xint_old + T(4.0)*xint_old*xint_old*(T(1.5) - xint_old - xint_old*xint_old));
226 sx_old[4] = (T(1.0))/(T(24.0))*(T(0.5) + xint_old)*(T(0.5) + xint_old)*(T(0.5) + xint_old)*(T(0.5)+xint_old);
227 //
228 const T xint_new = xnew - T(j);
229 sx_new[0] = (T(1.0))/(T(24.0))*(T(0.5) - xint_new)*(T(0.5) - xint_new)*(T(0.5) - xint_new)*(T(0.5) - xint_new);
230 sx_new[1] = (T(1.0))/(T(24.0))*(T(4.75) - T(11.0)*xint_new + T(4.0)*xint_new*xint_new*(T(1.5) + xint_new - xint_new*xint_new));
231 sx_new[2] = (T(1.0))/(T(24.0))*(T(14.375) + T(6.0)*xint_new*xint_new*(xint_new*xint_new - T(2.5)));
232 sx_new[3] = (T(1.0))/(T(24.0))*(T(4.75) + T(11.0)*xint_new + T(4.0)*xint_new*xint_new*(T(1.5) - xint_new - xint_new*xint_new));
233 sx_new[4] = (T(1.0))/(T(24.0))*(T(0.5) + xint_new)*(T(0.5) + xint_new)*(T(0.5) + xint_new)*(T(0.5)+xint_new);
234 //
235 // index of the leftmost cell where particle deposits
236 return j-2;
237 }
238 else{
239 WARPX_ABORT_WITH_MESSAGE("Unknown particle shape selected in Compute_shape_factor_pair");
240 amrex::ignore_unused(sx_old, sx_new, xold, xnew);
241 }
242 return 0;
243 }
244};
245
246#endif // WARPX_SHAPEFACTORS_H_
#define AMREX_FORCE_INLINE
#define AMREX_GPU_HOST_DEVICE
#define WARPX_ABORT_WITH_MESSAGE(MSG)
Definition TextMsg.H:15
__host__ __device__ void ignore_unused(const Ts &...)
Definition ShapeFactors.H:168
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int operator()(T *const sx_old, T *const sx_new, T xold, T xnew) const
Definition ShapeFactors.H:171
Definition ShapeFactors.H:29
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int operator()(T *const sx, T xmid) const
Definition ShapeFactors.H:32
Definition ShapeFactors.H:95
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int operator()(T *const sx, const T x_old, const int i_new) const
Definition ShapeFactors.H:98