From 45bcae87e36241e75d08fc0ce5420cc35259178f Mon Sep 17 00:00:00 2001 From: Justin Angus Date: Tue, 9 Jul 2024 17:43:40 -0700 Subject: [PATCH] fixed ComputeTemperature routine called from Coulomb collisions to work for particles of varying weight. --- .../BinaryCollision/Coulomb/ComputeTemperature.H | 15 +++++++++------ .../Coulomb/ElasticCollisionPerez.H | 4 ++-- 2 files changed, 11 insertions(+), 8 deletions(-) diff --git a/Source/Particles/Collision/BinaryCollision/Coulomb/ComputeTemperature.H b/Source/Particles/Collision/BinaryCollision/Coulomb/ComputeTemperature.H index 04440288ba3..6d0b48615a4 100644 --- a/Source/Particles/Collision/BinaryCollision/Coulomb/ComputeTemperature.H +++ b/Source/Particles/Collision/BinaryCollision/Coulomb/ComputeTemperature.H @@ -14,6 +14,7 @@ template 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 w, T_R const * AMREX_RESTRICT ux, T_R const * AMREX_RESTRICT uy, T_R const * AMREX_RESTRICT uz, T_R const m ) { @@ -26,6 +27,7 @@ T_R ComputeTemperature ( T_R vx = T_R(0.0); T_R vy = T_R(0.0); T_R vz = T_R(0.0); T_R vs = T_R(0.0); T_R gm = T_R(0.0); T_R us = T_R(0.0); + T_R wtot = T_R(0.0); for (int i = Is; i < static_cast(Ie); ++i) { @@ -33,14 +35,15 @@ T_R ComputeTemperature ( uy[ I[i] ] * uy[ I[i] ] + uz[ I[i] ] * uz[ I[i] ] ); gm = std::sqrt( T_R(1.0) + us*inv_c2 ); - vx += ux[ I[i] ] / gm; - vy += uy[ I[i] ] / gm; - vz += uz[ I[i] ] / gm; - vs += us / gm / gm; + wtot += w[ I[i] ]; + vx += w[ I[i] ] * ux[ I[i] ] / gm; + vy += w[ I[i] ] * uy[ I[i] ] / gm; + vz += w[ I[i] ] * uz[ I[i] ] / gm; + vs += w[ I[i] ] * us / gm / gm; } - vx = vx / N; vy = vy / N; - vz = vz / N; vs = vs / N; + vx = vx / wtot; vy = vy / wtot; + vz = vz / wtot; vs = vs / wtot; return m/T_R(3.0)*(vs-(vx*vx+vy*vy+vz*vz)); } diff --git a/Source/Particles/Collision/BinaryCollision/Coulomb/ElasticCollisionPerez.H b/Source/Particles/Collision/BinaryCollision/Coulomb/ElasticCollisionPerez.H index af4e6ba38c9..a782fac5f6e 100644 --- a/Source/Particles/Collision/BinaryCollision/Coulomb/ElasticCollisionPerez.H +++ b/Source/Particles/Collision/BinaryCollision/Coulomb/ElasticCollisionPerez.H @@ -74,12 +74,12 @@ void ElasticCollisionPerez ( T_PR T1t; T_PR T2t; if ( T1 <= T_PR(0.0) && L <= T_PR(0.0) ) { - T1t = ComputeTemperature(I1s,I1e,I1,u1x,u1y,u1z,m1); + T1t = ComputeTemperature(I1s,I1e,I1,w1,u1x,u1y,u1z,m1); } else { T1t = T1; } if ( T2 <= T_PR(0.0) && L <= T_PR(0.0) ) { - T2t = ComputeTemperature(I2s,I2e,I2,u2x,u2y,u2z,m2); + T2t = ComputeTemperature(I2s,I2e,I2,w2,u2x,u2y,u2z,m2); } else { T2t = T2; }