From adc6c03614844b698a4c75bc8f2a535fe7a2cc56 Mon Sep 17 00:00:00 2001 From: Emanuele Marchetti Date: Mon, 9 Feb 2026 16:52:35 +0100 Subject: [PATCH 1/4] Temperature deposition shape factor - General implementation --- docs/source/run/parameters.rst | 5 + src/Hipace.H | 2 + src/Hipace.cpp | 1 + src/fields/Fields.cpp | 5 + .../deposition/PlasmaDepositCurrent.cpp | 2 +- .../deposition/TemperatureDeposition.cpp | 96 +++++++++++++------ 6 files changed, 79 insertions(+), 32 deletions(-) diff --git a/docs/source/run/parameters.rst b/docs/source/run/parameters.rst index 501e7b6ea1..6af8fa7d0a 100644 --- a/docs/source/run/parameters.rst +++ b/docs/source/run/parameters.rst @@ -1158,6 +1158,11 @@ Field diagnostics will be deposited into individual fields accessible as ``w``, ``ux_`` or ``ux^2_`` (similarly for ``uy`` and ``uz``) in ``diagnostic.field_data``. +* ``hipace.temp_depos_order`` (`int`) optional (default `2`) + When ``hipace.deposit_temp_individual`` is turned on, + this option specifies the shape order of the deposited fields. + Currently, 0,1,2,3 are implemented. + In-situ diagnostics ^^^^^^^^^^^^^^^^^^^ diff --git a/src/Hipace.H b/src/Hipace.H index 240e357775..8d29d980aa 100644 --- a/src/Hipace.H +++ b/src/Hipace.H @@ -242,6 +242,8 @@ public: inline static bool m_deposit_rho_individual = false; /** Whether to deposit temperature for every individual plasma for diagnostics */ inline static bool m_deposit_temp_individual = false; + /** Order of the temperature deposition shape factor */ + inline static int m_temp_depos_order = 2; /** Whether to interpolate the neutralizing background to MR levels 1 and 2 instead of depositing it */ inline static bool m_interpolate_neutralizing_background = false; /** Whether to use tiling for particle operations */ diff --git a/src/Hipace.cpp b/src/Hipace.cpp index 5d551e7dab..2d205f73c1 100644 --- a/src/Hipace.cpp +++ b/src/Hipace.cpp @@ -130,6 +130,7 @@ Hipace::ReadParameters () queryWithParser(pph, "deposit_rho_individual", m_deposit_rho_individual); m_deposit_temp_individual = m_diags.needsTempIndividual(); queryWithParser(pph, "deposit_temp_individual", m_deposit_temp_individual); + queryWithParser(pph, "temp_depos_order", m_temp_depos_order); queryWithParser(pph, "interpolate_neutralizing_background", m_interpolate_neutralizing_background); bool do_mfi_sync = false; diff --git a/src/fields/Fields.cpp b/src/fields/Fields.cpp index 4eb0208938..2693fcd305 100644 --- a/src/fields/Fields.cpp +++ b/src/fields/Fields.cpp @@ -62,6 +62,11 @@ Fields::AllocData ( // Need 1 extra guard cell transversally for transverse derivative int nguards_xy = (Hipace::m_depos_order_xy + 1) / 2 + 1; + // Check the temperature deposition order, if enabled + if (Hipace::m_deposit_temp_individual && + Hipace::m_temp_depos_order>Hipace::m_depos_order_xy) { + nguards_xy = (Hipace::m_temp_depos_order + 1) / 2 + 1; + } m_slices_nguards = amrex::IntVect{nguards_xy, nguards_xy, 0}; m_explicit = Hipace::m_explicit; diff --git a/src/particles/deposition/PlasmaDepositCurrent.cpp b/src/particles/deposition/PlasmaDepositCurrent.cpp index 35384e4515..e37d94675f 100644 --- a/src/particles/deposition/PlasmaDepositCurrent.cpp +++ b/src/particles/deposition/PlasmaDepositCurrent.cpp @@ -131,7 +131,7 @@ DepositCurrent (PlasmaParticleContainer& plasma, Fields & fields, [=] AMREX_GPU_DEVICE (int ip, auto ptd, auto /*depos_order*/, auto /*can_ionize*/, - auto /*use_laserr*/) + auto /*use_laser*/) { // only deposit plasma currents on or below their according MR level return ptd.id(ip).is_valid() && (lev == 0 || ptd.cpu(ip) >= lev); diff --git a/src/particles/deposition/TemperatureDeposition.cpp b/src/particles/deposition/TemperatureDeposition.cpp index 1324a77397..4b85ce7baf 100644 --- a/src/particles/deposition/TemperatureDeposition.cpp +++ b/src/particles/deposition/TemperatureDeposition.cpp @@ -54,22 +54,53 @@ DepositTemperature (PlasmaParticleContainer& plasma, // extract laser properties and boolean for the presence of the laser and for ionization const PhysConst pc = get_phys_const(); const int aabs = Hipace::m_use_laser ? Comps[WhichSlice::This]["aabs"] : -1; - const bool can_ionize = plasma.m_can_ionize; - const bool use_laser = Hipace::m_use_laser; const amrex::Real laser_norm = (plasma.m_charge/pc.q_e) * (pc.m_e/plasma.m_mass) * (plasma.m_charge/pc.q_e) * (pc.m_e/plasma.m_mass); - + // Loop over particles - SharedMemoryDeposition<1, 1, true>( - int(pti.numParticles()), + amrex::AnyCTO( + // use compile-time options + amrex::TypeList< + amrex::CompileTimeOptions<0, 1, 2, 3>, // depos_order + amrex::CompileTimeOptions, // can_ionize + amrex::CompileTimeOptions // use_laser + >{}, { + Hipace::m_temp_depos_order, + plasma.m_can_ionize, + Hipace::m_use_laser + }, + // call deposition function + // The three functions passed as arguments to this lambda + // are defined below as the next arguments. + [&](auto is_valid, auto get_start_cell, auto deposit){ + constexpr auto ctos = deposit.GetOptions(); + constexpr int depos_order = ctos[0]; + constexpr int use_laser = ctos[2]; + constexpr int stencil_size = depos_order + 1; + + SharedMemoryDeposition( + int(pti.numParticles()), is_valid, get_start_cell, deposit, isl_fab.array(), + isl_fab.box(), pti.GetParticleTile().getParticleTileData(), + amrex::GpuArray{aabs}, + amrex::GpuArray{w, ux, uy, uz, uxsq, uysq, uzsq}); + }, + // is_valid // return whether the particle is valid and should deposit - [=] AMREX_GPU_DEVICE (int ip, auto ptd) + [=] AMREX_GPU_DEVICE (int ip, auto ptd, + auto /*depos_order*/, + auto /*can_ionize*/, + auto /*use_laser*/) { // only deposit on or below their according MR level return ptd.id(ip).is_valid() && (lev == 0 || ptd.cpu(ip) >= lev); }, - [=] AMREX_GPU_DEVICE (int ip, auto ptd) -> amrex::IntVectND<2> + // get_start_cell + // return the lowest cell index that the particle deposits into + [=] AMREX_GPU_DEVICE (int ip, auto ptd, + auto depos_order, + auto /*can_ionize*/, + auto /*use_laser*/) -> amrex::IntVectND<2> { const amrex::Real xp = ptd.pos(0, ip); const amrex::Real yp = ptd.pos(1, ip); @@ -78,17 +109,21 @@ DepositTemperature (PlasmaParticleContainer& plasma, const amrex::Real ymid = (yp - y_pos_offset) * dy_inv; auto [shape_x, i] = - compute_single_shape_factor(xmid, 0); + compute_single_shape_factor(xmid, 0); auto [shape_y, j] = - compute_single_shape_factor(ymid, 0); + compute_single_shape_factor(ymid, 0); return {i, j}; }, + // do_deposit // deposit of weight, momentum (ux, uy, uz) and their squares (uxsq, uysq, uzsq) [=] AMREX_GPU_DEVICE (int ip, auto ptd, - Array3 arr, - auto cache_idx, auto depos_idx) noexcept + Array3 arr, + auto cache_idx, auto depos_idx, + auto depos_order, + auto can_ionize, + auto use_laser) noexcept { const amrex::Real xp = ptd.pos(0, ip); const amrex::Real yp = ptd.pos(1, ip); @@ -101,7 +136,7 @@ DepositTemperature (PlasmaParticleContainer& plasma, ptd.idata(PlasmaIdx::ion_lev)[ip] * ptd.idata(PlasmaIdx::ion_lev)[ip]; } doLaserGatherShapeN<2>(xp, yp, Aabssqp, arr, cache_idx[0], - dx_inv, dy_inv, x_pos_offset, y_pos_offset); + dx_inv, dy_inv, x_pos_offset, y_pos_offset); Aabssqp *= laser_norm_ion; } @@ -116,25 +151,24 @@ DepositTemperature (PlasmaParticleContainer& plasma, const amrex::Real xmid = (xp - x_pos_offset) * dx_inv; const amrex::Real ymid = (yp - y_pos_offset) * dy_inv; - // --- Compute shape factors - // x direction - auto [shape_x, i] = compute_single_shape_factor(xmid, 0); - // y direction - auto [shape_y, j] = compute_single_shape_factor(ymid, 0); - - amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[0]), wp); - amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[1]), wp*uxp); - amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[2]), wp*uyp); - amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[3]), wp*uzp); - amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[4]), wp*uxp*uxp); - amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[5]), wp*uyp*uyp); - amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[6]), wp*uzp*uzp); - }, - isl_fab.array(), - isl_fab.box(), pti.GetParticleTile().getParticleTileData(), - amrex::GpuArray{aabs}, - amrex::GpuArray{w, ux, uy, uz, uxsq, uysq, uzsq} - ); + for (int iy=0; iy <= depos_order; ++iy) { + for (int ix=0; ix <= depos_order; ++ix) { + // --- Compute shape factors + // x direction + auto [shape_x, i] = compute_single_shape_factor(xmid, ix); + // y direction + auto [shape_y, j] = compute_single_shape_factor(ymid, iy); + + amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[0]), shape_x*shape_y*wp); + amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[1]), shape_x*shape_y*wp*uxp); + amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[2]), shape_x*shape_y*wp*uyp); + amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[3]), shape_x*shape_y*wp*uzp); + amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[4]), shape_x*shape_y*wp*uxp*uxp); + amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[5]), shape_x*shape_y*wp*uyp*uyp); + amrex::Gpu::Atomic::Add(arr.ptr(i, j, depos_idx[6]), shape_x*shape_y*wp*uzp*uzp); + } + } + }); Array3 field_arr = isl_fab.array(); // Normalize the components of momentum (ux, uy, uz) and their squares (uxsq, uysq, uzsq) From 7d1f4be90fcef8422740be6a32277ce1841ba63c Mon Sep 17 00:00:00 2001 From: Emanuele Marchetti Date: Thu, 19 Feb 2026 11:07:54 +0100 Subject: [PATCH 2/4] Fix unused variable and remove trailing whitespace --- docs/source/run/parameters.rst | 2 +- src/fields/Fields.cpp | 2 +- .../deposition/TemperatureDeposition.cpp | 25 ++++++++++++------- 3 files changed, 18 insertions(+), 11 deletions(-) diff --git a/docs/source/run/parameters.rst b/docs/source/run/parameters.rst index 6af8fa7d0a..85656506b1 100644 --- a/docs/source/run/parameters.rst +++ b/docs/source/run/parameters.rst @@ -1159,7 +1159,7 @@ Field diagnostics ``ux^2_`` (similarly for ``uy`` and ``uz``) in ``diagnostic.field_data``. * ``hipace.temp_depos_order`` (`int`) optional (default `2`) - When ``hipace.deposit_temp_individual`` is turned on, + When ``hipace.deposit_temp_individual`` is turned on, this option specifies the shape order of the deposited fields. Currently, 0,1,2,3 are implemented. diff --git a/src/fields/Fields.cpp b/src/fields/Fields.cpp index 2693fcd305..0735d4c7a5 100644 --- a/src/fields/Fields.cpp +++ b/src/fields/Fields.cpp @@ -63,7 +63,7 @@ Fields::AllocData ( // Need 1 extra guard cell transversally for transverse derivative int nguards_xy = (Hipace::m_depos_order_xy + 1) / 2 + 1; // Check the temperature deposition order, if enabled - if (Hipace::m_deposit_temp_individual && + if (Hipace::m_deposit_temp_individual && Hipace::m_temp_depos_order>Hipace::m_depos_order_xy) { nguards_xy = (Hipace::m_temp_depos_order + 1) / 2 + 1; } diff --git a/src/particles/deposition/TemperatureDeposition.cpp b/src/particles/deposition/TemperatureDeposition.cpp index cf46cfb27e..50bbc3f650 100644 --- a/src/particles/deposition/TemperatureDeposition.cpp +++ b/src/particles/deposition/TemperatureDeposition.cpp @@ -56,7 +56,7 @@ DepositTemperature (PlasmaParticleContainer& plasma, const int aabs = Hipace::m_use_laser ? Comps[WhichSlice::This]["aabs"] : -1; const amrex::Real laser_norm = (plasma.m_charge/pc.q_e) * (pc.m_e/plasma.m_mass) * (plasma.m_charge/pc.q_e) * (pc.m_e/plasma.m_mass); - + // Loop over particles amrex::AnyCTO( // use compile-time options @@ -78,16 +78,23 @@ DepositTemperature (PlasmaParticleContainer& plasma, constexpr int use_laser = ctos[2]; constexpr int stencil_size = depos_order + 1; - SharedMemoryDeposition( - int(pti.numParticles()), is_valid, get_start_cell, deposit, isl_fab.array(), - isl_fab.box(), pti.GetParticleTile().getParticleTileData(), - amrex::GpuArray{aabs}, - amrex::GpuArray{w, ux, uy, uz, uxsq, uysq, uzsq}); + if constexpr (use_laser) { + SharedMemoryDeposition( + int(pti.numParticles()), is_valid, get_start_cell, deposit, isl_fab.array(), + isl_fab.box(), pti.GetParticleTile().getParticleTileData(), + amrex::GpuArray{aabs}, + amrex::GpuArray{w, ux, uy, uz, uxsq, uysq, uzsq}); + } else { + SharedMemoryDeposition( + int(pti.numParticles()), is_valid, get_start_cell, deposit, isl_fab.array(), + isl_fab.box(), pti.GetParticleTile().getParticleTileData(), + amrex::GpuArray{}, + amrex::GpuArray{w, ux, uy, uz, uxsq, uysq, uzsq}); + } }, - // is_valid // return whether the particle is valid and should deposit - [=] AMREX_GPU_DEVICE (int ip, auto ptd, + [=] AMREX_GPU_DEVICE (int ip, auto ptd, auto /*depos_order*/, auto /*can_ionize*/, auto /*use_laser*/) @@ -97,7 +104,7 @@ DepositTemperature (PlasmaParticleContainer& plasma, }, // get_start_cell // return the lowest cell index that the particle deposits into - [=] AMREX_GPU_DEVICE (int ip, auto ptd, + [=] AMREX_GPU_DEVICE (int ip, auto ptd, auto depos_order, auto /*can_ionize*/, auto /*use_laser*/) -> amrex::IntVectND<2> From aea912501b280518649ec7ddc6cdc288ecad2057 Mon Sep 17 00:00:00 2001 From: Emanuele Marchetti Date: Thu, 19 Feb 2026 18:02:49 +0100 Subject: [PATCH 3/4] Rename temp -> temperature and fix CI --- docs/source/run/parameters.rst | 2 +- src/Hipace.H | 2 +- src/Hipace.cpp | 2 +- src/fields/Fields.cpp | 4 ++-- src/particles/deposition/TemperatureDeposition.cpp | 2 +- tests/laser_ionization.1Rank.sh | 2 ++ 6 files changed, 8 insertions(+), 6 deletions(-) diff --git a/docs/source/run/parameters.rst b/docs/source/run/parameters.rst index 85656506b1..1b7333684e 100644 --- a/docs/source/run/parameters.rst +++ b/docs/source/run/parameters.rst @@ -1158,7 +1158,7 @@ Field diagnostics will be deposited into individual fields accessible as ``w``, ``ux_`` or ``ux^2_`` (similarly for ``uy`` and ``uz``) in ``diagnostic.field_data``. -* ``hipace.temp_depos_order`` (`int`) optional (default `2`) +* ``hipace.temperature_depos_order`` (`int`) optional (default `2`) When ``hipace.deposit_temp_individual`` is turned on, this option specifies the shape order of the deposited fields. Currently, 0,1,2,3 are implemented. diff --git a/src/Hipace.H b/src/Hipace.H index 8d29d980aa..0a6c499148 100644 --- a/src/Hipace.H +++ b/src/Hipace.H @@ -243,7 +243,7 @@ public: /** Whether to deposit temperature for every individual plasma for diagnostics */ inline static bool m_deposit_temp_individual = false; /** Order of the temperature deposition shape factor */ - inline static int m_temp_depos_order = 2; + inline static int m_temperature_depos_order = 2; /** Whether to interpolate the neutralizing background to MR levels 1 and 2 instead of depositing it */ inline static bool m_interpolate_neutralizing_background = false; /** Whether to use tiling for particle operations */ diff --git a/src/Hipace.cpp b/src/Hipace.cpp index 2d205f73c1..c7ffaa63e3 100644 --- a/src/Hipace.cpp +++ b/src/Hipace.cpp @@ -130,7 +130,7 @@ Hipace::ReadParameters () queryWithParser(pph, "deposit_rho_individual", m_deposit_rho_individual); m_deposit_temp_individual = m_diags.needsTempIndividual(); queryWithParser(pph, "deposit_temp_individual", m_deposit_temp_individual); - queryWithParser(pph, "temp_depos_order", m_temp_depos_order); + queryWithParser(pph, "temperature_depos_order", m_temperature_depos_order); queryWithParser(pph, "interpolate_neutralizing_background", m_interpolate_neutralizing_background); bool do_mfi_sync = false; diff --git a/src/fields/Fields.cpp b/src/fields/Fields.cpp index 0735d4c7a5..069e0a6872 100644 --- a/src/fields/Fields.cpp +++ b/src/fields/Fields.cpp @@ -64,8 +64,8 @@ Fields::AllocData ( int nguards_xy = (Hipace::m_depos_order_xy + 1) / 2 + 1; // Check the temperature deposition order, if enabled if (Hipace::m_deposit_temp_individual && - Hipace::m_temp_depos_order>Hipace::m_depos_order_xy) { - nguards_xy = (Hipace::m_temp_depos_order + 1) / 2 + 1; + Hipace::m_temperature_depos_order>Hipace::m_depos_order_xy) { + nguards_xy = (Hipace::m_temperature_depos_order + 1) / 2 + 1; } m_slices_nguards = amrex::IntVect{nguards_xy, nguards_xy, 0}; diff --git a/src/particles/deposition/TemperatureDeposition.cpp b/src/particles/deposition/TemperatureDeposition.cpp index 50bbc3f650..5f7ecbfb9d 100644 --- a/src/particles/deposition/TemperatureDeposition.cpp +++ b/src/particles/deposition/TemperatureDeposition.cpp @@ -65,7 +65,7 @@ DepositTemperature (PlasmaParticleContainer& plasma, amrex::CompileTimeOptions, // can_ionize amrex::CompileTimeOptions // use_laser >{}, { - Hipace::m_temp_depos_order, + Hipace::m_temperature_depos_order, plasma.m_can_ionize, Hipace::m_use_laser }, diff --git a/tests/laser_ionization.1Rank.sh b/tests/laser_ionization.1Rank.sh index eab7146909..83499c80ce 100755 --- a/tests/laser_ionization.1Rank.sh +++ b/tests/laser_ionization.1Rank.sh @@ -30,6 +30,7 @@ mpiexec -n 1 $HIPACE_EXECUTABLE $HIPACE_EXAMPLE_DIR/inputs_laser_ionization \ my_constants.a0 = 0.00885126 \ plasmas.do_push = 0 \ hipace.file_prefix=$TEST_NAME/linear \ + hipace.temperature_depos_order = 0 \ plasmas.insitu_file_prefix = $TEST_NAME/insitu_linear mpiexec -n 1 $HIPACE_EXECUTABLE $HIPACE_EXAMPLE_DIR/inputs_laser_ionization \ @@ -37,6 +38,7 @@ mpiexec -n 1 $HIPACE_EXECUTABLE $HIPACE_EXAMPLE_DIR/inputs_laser_ionization \ lasers.polarization = circular \ plasmas.do_push = 0 \ hipace.file_prefix=$TEST_NAME/circular \ + hipace.temperature_depos_order = 0 \ plasmas.insitu_file_prefix = $TEST_NAME/insitu_circular # Compare the result with theory From 90a74ffa85534e3828f128b46413d1680b14156a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Maxence=20Th=C3=A9venet?= Date: Tue, 17 Mar 2026 10:38:56 +0100 Subject: [PATCH 4/4] Apply suggestion from @MaxThevenet --- src/fields/Fields.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/fields/Fields.cpp b/src/fields/Fields.cpp index e62b4e631f..78f69dbc51 100644 --- a/src/fields/Fields.cpp +++ b/src/fields/Fields.cpp @@ -65,7 +65,7 @@ Fields::AllocData ( int nguards_xy = (Hipace::m_depos_order_xy + 1) / 2 + 1; // Check the temperature deposition order, if enabled if (Hipace::m_deposit_temp_individual && - Hipace::m_temperature_depos_order>Hipace::m_depos_order_xy) { + Hipace::m_temperature_depos_order > Hipace::m_depos_order_xy) { nguards_xy = (Hipace::m_temperature_depos_order + 1) / 2 + 1; } m_slices_nguards = amrex::IntVect{nguards_xy, nguards_xy, 0};