From 277dafddc5541252d7b90139ea06de2bc5b8332e Mon Sep 17 00:00:00 2001 From: Marco Acciarri Date: Mon, 25 Nov 2024 14:58:56 -0800 Subject: [PATCH] fixed double call to qdsmc_hybrid_electron_pc->InitParticles and added missing particle weight to HybridInitializeKe and HybridUpdateTe methods. --- Source/Fluids/QdsmcParticleContainer.H | 4 ++-- Source/Fluids/QdsmcParticleContainer.cpp | 10 +++------- Source/Fluids/Qdsmc_K.H | 13 ++++++++----- Source/Fluids/WarpXFluidContainer.cpp | 14 +++++++++++--- 4 files changed, 24 insertions(+), 17 deletions(-) diff --git a/Source/Fluids/QdsmcParticleContainer.H b/Source/Fluids/QdsmcParticleContainer.H index 76005bd790b..495f3cbe35a 100644 --- a/Source/Fluids/QdsmcParticleContainer.H +++ b/Source/Fluids/QdsmcParticleContainer.H @@ -21,8 +21,8 @@ #include "QdsmcParticleContainer_fwd.H" /** - * This enumerated struct is used to index the field probe particle - * values that are being stored as SoA data. Nattribs + * This enumerated struct is used to index the qdsmc fictitious + * particle values that are being stored as SoA data. Nattribs * is enumerated to give the number of attributes stored. */ struct QdsmcPIdx diff --git a/Source/Fluids/QdsmcParticleContainer.cpp b/Source/Fluids/QdsmcParticleContainer.cpp index d16fb39da01..9bb78546678 100644 --- a/Source/Fluids/QdsmcParticleContainer.cpp +++ b/Source/Fluids/QdsmcParticleContainer.cpp @@ -63,7 +63,7 @@ QdsmcParticleContainer::QdsmcParticleContainer (AmrCore* amr_core) : ParticleContainerPureSoA(amr_core->GetParGDB()) { SetParticleSize(); - InitParticles(0); // only level 0 is used in HybridSolver + //InitParticles(0); // only level 0 is used in HybridSolver, commented since this is already called in WarpX.cpp } @@ -214,10 +214,6 @@ QdsmcParticleContainer::SetV (int lev, auto const np = pti.numParticles(); auto& attribs = pti.GetStructOfArrays().GetRealData(); - // making box cell centered - // amrex::Box box = pti.tilebox(); - // box.grow(Ux.nGrowVect()); - amrex::Box box = pti.tilebox(); const amrex::XDim3 xyzmin = WarpX::LowerCorner(box, lev, 0._rt); @@ -233,8 +229,8 @@ QdsmcParticleContainer::SetV (int lev, const auto &arrUyfield = Uy[pti].array(); const auto &arrUzfield = Uz[pti].array(); - // Gather drift velocity directly from nodes - // since particles are located at the node positions before PushX + // Gather drift velocity directly from nodes since + // particles are located at the node positions before PushX amrex::ParallelFor( np, [=] AMREX_GPU_DEVICE (long ip) { amrex::Real vxp; diff --git a/Source/Fluids/Qdsmc_K.H b/Source/Fluids/Qdsmc_K.H index 4139f5a291b..dfd19a4b641 100644 --- a/Source/Fluids/Qdsmc_K.H +++ b/Source/Fluids/Qdsmc_K.H @@ -17,10 +17,10 @@ /** - * Qdsmc particles x0,y0,z0 positions should be cell centered + * Qdsmc particles x0,y0,z0 positions should at nodes positions * at the beginning of every step, so the field gathering * does not require any interpolation, it only needs to read - * the values from the cell centered nodes. + * the values from the nodes. */ AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void gather_density_entropy(const amrex::ParticleReal x0p, @@ -48,9 +48,10 @@ void gather_density_entropy(const amrex::ParticleReal x0p, /** - * Qdsmc particles x0,y0,z0 positions should match node positions + * Qdsmc particles x0,y0,z0 positions should at nodes positions * at the beginning of every step, so the field gathering - * does not require any interpolation. + * does not require any interpolation, it only needs to read + * the values from the nodes. */ AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void gather_vector_field_qdsmc(const amrex::ParticleReal x0p, @@ -106,7 +107,9 @@ void push_qdsmc_particle(const amrex::ParticleReal x0p, */ AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void do_deposit_scalar(amrex::Array4 const& scalar_field, - const amrex::ParticleReal xp, const amrex::ParticleReal yp, const amrex::ParticleReal zp, + const amrex::ParticleReal xp, + const amrex::ParticleReal yp, + const amrex::ParticleReal zp, const amrex::XDim3 & xyzmin, const amrex::XDim3 & dinv, const amrex::ParticleReal valuep) diff --git a/Source/Fluids/WarpXFluidContainer.cpp b/Source/Fluids/WarpXFluidContainer.cpp index 556ccdf6bdb..ec25203e2fa 100644 --- a/Source/Fluids/WarpXFluidContainer.cpp +++ b/Source/Fluids/WarpXFluidContainer.cpp @@ -1478,7 +1478,7 @@ void WarpXFluidContainer::HybridInitializeKe (ablastr::fields::MultiFabRegister& { using warpx::fields::FieldType; ablastr::fields::ScalarField rho_fp = fields.get(FieldType::rho_fp, lev); - + #ifdef AMREX_USE_OMP #pragma omp parallel if (amrex::Gpu::notInLaunchRegion()) #endif @@ -1492,7 +1492,8 @@ void WarpXFluidContainer::HybridInitializeKe (ablastr::fields::MultiFabRegister& ParallelFor(tilebox, [=] AMREX_GPU_DEVICE (int i, int j, int k) { if( rho(i, j, k) > 0.0_rt ){ - Ke(i, j, k) = Te(i, j, k)*std::pow((rho(i, j, k)/PhysConst::q_e), 1-gamma); + amrex::Real ne = rho(i, j, k)/PhysConst::q_e; + Ke(i, j, k) = Te(i, j, k)*std::pow(ne, 1-gamma); } }); } @@ -1504,6 +1505,11 @@ void WarpXFluidContainer::HybridUpdateTe (ablastr::fields::MultiFabRegister& fie using warpx::fields::FieldType; ablastr::fields::ScalarField rho_fp = fields.get(FieldType::rho_fp, lev); + WarpX &warpx = WarpX::GetInstance(); + const amrex::Geometry &geom = warpx.Geom(lev); + const auto dx = geom.CellSizeArray(); + const auto cell_volume = dx[0]*dx[1]*dx[2]; + #ifdef AMREX_USE_OMP #pragma omp parallel if (amrex::Gpu::notInLaunchRegion()) #endif @@ -1518,7 +1524,9 @@ void WarpXFluidContainer::HybridUpdateTe (ablastr::fields::MultiFabRegister& fie ParallelFor(tilebox, [=] AMREX_GPU_DEVICE (int i, int j, int k) { Te(i, j, k) = 0; if( rho(i, j, k) > 0.0_rt ){ - Te(i, j, k) = Ke(i, j, k)*std::pow((rho(i, j, k)/PhysConst::q_e), gamma-1)*PhysConst::q_e; + amrex::Real ne = rho(i, j, k)/PhysConst::q_e; + amrex::Real N_cell = ne*cell_volume; + Te(i, j, k) = (Ke(i, j, k)/std::pow(ne, 1-gamma))/N_cell; } }); }