7 #ifndef WARPX_PARTICLES_COLLISION_ELASTIC_COLLISION_PEREZ_H_
8 #define WARPX_PARTICLES_COLLISION_ELASTIC_COLLISION_PEREZ_H_
43 template <
typename T_index,
typename T_PR,
typename T_R,
typename SoaData_type>
46 T_index
const I1s, T_index
const I1e,
47 T_index
const I2s, T_index
const I2e,
48 T_index
const* AMREX_RESTRICT I1,
49 T_index
const* AMREX_RESTRICT I2,
50 SoaData_type soa_1, SoaData_type soa_2,
51 T_PR
const q1, T_PR
const q2,
52 T_PR
const m1, T_PR
const m2,
53 T_PR
const T1, T_PR
const T2,
54 T_R
const dt, T_PR
const L, T_R
const dV,
56 bool const isSameSpecies, T_index coll_idx)
58 const T_index NI1 = I1e - I1s;
59 const T_index NI2 = I2e - I2s;
75 if ( T1 <= T_PR(0.0) && L <= T_PR(0.0) )
80 if ( T2 <= T_PR(0.0) && L <= T_PR(0.0) )
90 for (T_index i1=I1s; i1<I1e; ++i1) { n1 += w1[ I1[i1] ]; }
91 for (T_index i2=I2s; i2<I2e; ++i2) { n2 += w2[ I2[i2] ]; }
98 n1 = n1 / dV; n2 = n2 / dV;
100 T_index i1 = I1s; T_index i2 = I2s;
101 for (T_index k = 0; k < max_N; ++k)
103 n12 +=
amrex::min( w1[ I1[i1] ], w2[ I2[i2] ] );
104 ++i1;
if ( i1 == I1e ) { i1 = I1s; }
105 ++i2;
if ( i2 == I2e ) { i2 = I2s; }
110 if (isSameSpecies) { n12 *= T_PR(2.0); }
114 if ( T1t < T_PR(0.0) || T2t < T_PR(0.0) ) {
121 T_PR rmin = std::pow( T_PR(4.0) * MathConst::pi / T_PR(3.0) *
125 #if (defined WARPX_DIM_RZ)
132 T_index i1 = I1s + coll_idx;
133 T_index i2 = I2s + coll_idx;
137 for (T_index k = coll_idx; k < max_N; k += min_N)
139 #if (defined WARPX_DIM_RZ)
148 T_PR
const theta = theta2[I2[i2]]-theta1[I1[i1]];
149 T_PR
const u1xbuf = u1x[I1[i1]];
150 u1x[I1[i1]] = u1xbuf*std::cos(theta) - u1y[I1[i1]]*std::sin(theta);
151 u1y[I1[i1]] = u1xbuf*std::sin(theta) + u1y[I1[i1]]*std::cos(theta);
155 u1x[ I1[i1] ], u1y[ I1[i1] ], u1z[ I1[i1] ],
156 u2x[ I2[i2] ], u2y[ I2[i2] ], u2z[ I2[i2] ],
158 q1, m1, w1[ I1[i1] ], q2, m2, w2[ I2[i2] ],
162 #if (defined WARPX_DIM_RZ)
163 T_PR
const u1xbuf_new = u1x[I1[i1]];
164 u1x[I1[i1]] = u1xbuf_new*std::cos(-theta) - u1y[I1[i1]]*std::sin(-theta);
165 u1y[I1[i1]] = u1xbuf_new*std::sin(-theta) + u1y[I1[i1]]*std::cos(-theta);
#define AMREX_GPU_HOST_DEVICE
AMREX_GPU_HOST_DEVICE T_R ComputeTemperature(T_index const Is, T_index const Ie, T_index const *AMREX_RESTRICT I, T_R const *AMREX_RESTRICT ux, T_R const *AMREX_RESTRICT uy, T_R const *AMREX_RESTRICT uz, T_R const m)
Definition: ComputeTemperature.H:15
AMREX_GPU_HOST_DEVICE AMREX_INLINE void ElasticCollisionPerez(T_index const I1s, T_index const I1e, T_index const I2s, T_index const I2e, T_index const *AMREX_RESTRICT I1, T_index const *AMREX_RESTRICT I2, SoaData_type soa_1, SoaData_type soa_2, T_PR const q1, T_PR const q2, T_PR const m1, T_PR const m2, T_PR const T1, T_PR const T2, T_R const dt, T_PR const L, T_R const dV, amrex::RandomEngine const &engine, bool const isSameSpecies, T_index coll_idx)
Definition: ElasticCollisionPerez.H:45
AMREX_GPU_HOST_DEVICE AMREX_INLINE void UpdateMomentumPerezElastic(T_PR &u1x, T_PR &u1y, T_PR &u1z, T_PR &u2x, T_PR &u2y, T_PR &u2z, T_PR const n1, T_PR const n2, T_PR const n12, T_PR const q1, T_PR const m1, T_PR const w1, T_PR const q2, T_PR const m2, T_PR const w2, T_R const dt, T_PR const L, T_PR const lmdD, amrex::RandomEngine const &engine)
Definition: UpdateMomentumPerezElastic.H:34
static constexpr auto ep0
vacuum permittivity: dielectric permittivity of vacuum [F/m]
Definition: constant.H:46
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE const T & max(const T &a, const T &b) noexcept
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE const T & min(const T &a, const T &b) noexcept
float dt
Definition: stencil.py:442
@ theta
RZ needs all three position components.
Definition: NamedComponentParticleContainer.H:36
@ uz
Definition: NamedComponentParticleContainer.H:34
@ w
weight
Definition: NamedComponentParticleContainer.H:33
@ uy
Definition: NamedComponentParticleContainer.H:34
@ ux
Definition: NamedComponentParticleContainer.H:34