Skip to content

Commit b87251a

Browse files
BrianMarreBrianMarre
committed
fix skipped applyIPDIOnization at end of FLYonPIC loop
1 parent 1465693 commit b87251a

File tree

6 files changed

+95
-49
lines changed

6 files changed

+95
-49
lines changed

‎include/picongpu/particles/atomicPhysics/ionizationPotentialDepression/StewartPyattIPD.hpp‎

Lines changed: 6 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
/* Copyright 2024-2024 Brian Marre
1+
/* Copyright 2024-2025 Brian Marre
22
*
33
* This file is part of PIConGPU.
44
*
@@ -198,13 +198,15 @@ namespace picongpu::particles::atomicPhysics::ionizationPotentialDepression
198198
*
199199
* @attention collective over all ion species
200200
*/
201-
template<typename T_AtomicPhysicsIonSpeciesList>
201+
template<typename T_AtomicPhysicsIonSpeciesList, bool T_SkipFinishedSuperCell>
202202
HINLINE static void applyIPDIonization(picongpu::MappingDesc const mappingDesc, uint32_t const)
203203
{
204204
using ForEachIonSpeciesApplyIPDIonization = pmacc::meta::ForEach<
205205
T_AtomicPhysicsIonSpeciesList,
206-
s_IPD::stage::ApplyIPDIonization<boost::mpl::_1, StewartPyattIPD<T_TemperatureFunctor>>>;
207-
206+
s_IPD::stage::ApplyIPDIonization<
207+
boost::mpl::_1,
208+
StewartPyattIPD<T_TemperatureFunctor>,
209+
std::integral_constant<bool, T_SkipFinishedSuperCell>>>;
208210
ForEachIonSpeciesApplyIPDIonization{}(mappingDesc);
209211
};
210212

‎include/picongpu/particles/atomicPhysics/ionizationPotentialDepression/kernel/ApplyIPDIonization.kernel‎

Lines changed: 17 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -198,21 +198,21 @@ namespace picongpu::particles::atomicPhysics::ionizationPotentialDepression::ker
198198
ipd += BarrierSupressionIonization::getIPD(screenedCharge, eFieldNormAU);
199199
}
200200

201+
201202
bool const stateIsUnbound = ((ionizationEnergy - ipd) <= 0._X);
202-
bool const stateIsGenerallyUnbound = (ionizationEnergy <= 0._X);
203203

204-
/** @details states that are unbound without an IPD contribution must relax via auto-ionization, electronic
205-
* collisional ionization or deexcitation channels in the regular rate solver since IPD lacks the energy
206-
* for ionization.
204+
/** @details states that are unbound without an IPD contribution should have option to relax via
205+
* auto-ionization, electronic collisional ionization or deexcitation channels in the regular rate solver
206+
* since IPD lacks the energy for ionization.
207207
*
208-
* These states are typically low shell hole states, IPD-Ionization will not relax these states until
209-
* very high charge states, causing:
208+
* This is currently not possible since IPD-ionization is instantaneous, this can create problems for low
209+
* shell hole states, IPD-Ionization will not relax these states until very high charge states, causing:
210210
* - numerical energy creation as the IPD actually lacks the energy required for ionization
211211
* - instant ionization cascades, which are not consistent with the IPD equilibrium description
212212
*
213213
* @todo improve by developing a non equlibirium IPD description, Brian Marre, 2025
214214
*/
215-
if(!stateIsGenerallyUnbound && stateIsUnbound)
215+
if(stateIsUnbound)
216216
{
217217
/* we only update the atomic state since IPD-Ionization is not a regular transition and does not use
218218
* shared resources */
@@ -271,10 +271,14 @@ namespace picongpu::particles::atomicPhysics::ionizationPotentialDepression::ker
271271
}
272272
};
273273

274+
///@todo switch to 2 body collision for IPD ionization to ensure energy conservation,
275+
/// even from highly excited states, Brian Marre, 2025
276+
274277
//! moduleConfig of the ApplyIPDIonizationKernel
278+
template<template<typename... T_KernelOptions> typename T_SkipSuperCellFunctor>
275279
using ApplyIPDIonizationModulConfig = particles::creation::ModuleConfig<
276280
ApplyIPDIonizationSanityCheckInputs,
277-
spawnFromSourceSpeciesModules::SkipFinishedSuperCellsAtomicPhysics,
281+
T_SkipSuperCellFunctor,
278282
IPDIonizationPredictor,
279283
spawnFromSourceSpeciesModules::InitAsCoMoving,
280284
KernelState,
@@ -284,11 +288,14 @@ namespace picongpu::particles::atomicPhysics::ionizationPotentialDepression::ker
284288
WriteFoundUnboundToSuperCellField>;
285289

286290
//! actual definition of ApplyIPDIonizationKernel
287-
template<typename T_IPDModel, typename T_fieldIonizationActive>
291+
template<
292+
template<typename... T_KernelOptions> typename T_SkipSuperCellFunctor,
293+
typename T_IPDModel,
294+
typename T_fieldIonizationActive>
288295
using ApplyIPDIonizationKernel = picongpu::particles::creation::SpawnFromSourceSpeciesKernelFramework<
289296
// T_TypeNumber
290297
uint8_t,
291-
ApplyIPDIonizationModulConfig,
298+
ApplyIPDIonizationModulConfig<T_SkipSuperCellFunctor>,
292299
T_IPDModel,
293300
T_fieldIonizationActive>;
294301
} // namespace picongpu::particles::atomicPhysics::ionizationPotentialDepression::kernel

‎include/picongpu/particles/atomicPhysics/ionizationPotentialDepression/stage/ApplyIPDIonization.def‎

Lines changed: 7 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
/* Copyright 2024-2024 Brian Marre
1+
/* Copyright 2024-2025 Brian Marre
22
*
33
* This file is part of PIConGPU.
44
*
@@ -31,9 +31,13 @@ namespace picongpu::particles::atomicPhysics::ionizationPotentialDepression::sta
3131
*
3232
* @todo implement version for non atomicPhysics data species
3333
*
34-
* @tparam ion species with atomic data
34+
* @tparam T_IonSpecies species with atomic data to apply IPD to
35+
* @tparam T_IPDModel ionization potential depression model to use
36+
* @tparam T_SkipFinishedSuperCell std::integral wrapped boolean,
37+
* true =^= skip finished super cells in the applyIPDIonization check
38+
* false =^= **do not** skip finished super cells in the applyIPDIonization check
3539
*/
36-
template<typename T_IonSpecies, typename T_IPDModel>
40+
template<typename T_IonSpecies, typename T_IPDModel, typename T_SkipFinishedSuperCell>
3741
struct ApplyIPDIonization
3842
{
3943
//! call of kernel for every superCell

‎include/picongpu/particles/atomicPhysics/ionizationPotentialDepression/stage/ApplyIPDIonization.hpp‎

Lines changed: 51 additions & 22 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
/* Copyright 2024-2024 Brian Marre
1+
/* Copyright 2024-2025 Brian Marre
22
*
33
* This file is part of PIConGPU.
44
*
@@ -26,6 +26,8 @@
2626
#include "picongpu/particles/atomicPhysics/ionizationPotentialDepression/stage/ApplyIPDIonization.def"
2727
#include "picongpu/particles/atomicPhysics/localHelperFields/FoundUnboundIonField.hpp"
2828
#include "picongpu/particles/atomicPhysics/localHelperFields/TimeRemainingField.hpp"
29+
#include "picongpu/particles/atomicPhysics/spawnFromSourceSpeciesModules/NeverSkipSuperCells.hpp"
30+
#include "picongpu/particles/atomicPhysics/spawnFromSourceSpeciesModules/SkipFinishedSuperCellsAtomicPhysics.hpp"
2931
#include "picongpu/particles/param.hpp"
3032
#include "picongpu/particles/traits/GetAtomicDataType.hpp"
3133
#include "picongpu/particles/traits/GetIonizationElectronSpecies.hpp"
@@ -39,8 +41,8 @@ namespace picongpu::particles::atomicPhysics::ionizationPotentialDepression::sta
3941
//! short hand for IPD namespace
4042
namespace s_IPD = picongpu::particles::atomicPhysics::ionizationPotentialDepression;
4143

42-
template<typename T_IonSpecies, typename T_IPDModel>
43-
HINLINE void ApplyIPDIonization<T_IonSpecies, T_IPDModel>::operator()(
44+
template<typename T_IonSpecies, typename T_IPDModel, typename T_SkipFinishedSuperCell>
45+
HINLINE void ApplyIPDIonization<T_IonSpecies, T_IPDModel, T_SkipFinishedSuperCell>::operator()(
4446
picongpu::MappingDesc const mappingDesc) const
4547
{
4648
// might be alias, from here on out no more
@@ -80,25 +82,52 @@ namespace picongpu::particles::atomicPhysics::ionizationPotentialDepression::sta
8082

8183
auto& fieldE = *dc.get<FieldE>(FieldE::getName());
8284

83-
// macro for call of kernel on every superCell, see pull request #4321
84-
PMACC_LOCKSTEP_KERNEL(
85-
s_IPD::kernel::ApplyIPDIonizationKernel<
86-
T_IPDModel,
87-
std::integral_constant<bool, AtomicDataType::switchFieldIonization>>())
88-
.config(mapper.getGridDim(), ions)(
89-
mapper,
90-
idProvider->getDeviceGenerator(),
91-
ions.getDeviceParticlesBox(),
92-
electrons.getDeviceParticlesBox(),
93-
timeRemainingField.getDeviceDataBox(),
94-
foundUnboundIonField.getDeviceDataBox(),
95-
atomicData.template getChargeStateDataDataBox</*on device*/ false>(),
96-
atomicData.template getAtomicStateDataDataBox</*on device*/ false>(),
97-
atomicData.template getIPDIonizationStateDataBox</*on device*/ false>(),
98-
fieldE.getDeviceDataBox(),
99-
debyeLengthField.getDeviceDataBox(),
100-
temperatureEnergyField.getDeviceDataBox(),
101-
zStarField.getDeviceDataBox());
85+
/** @details must use if-constexpr since "apply" in PMACC_LOCKSTEP_KERNEL macro is only able to handle
86+
* typenames, not typename aliases, ... for some reason. */
87+
if constexpr(T_SkipFinishedSuperCell::value)
88+
{
89+
PMACC_LOCKSTEP_KERNEL(
90+
s_IPD::kernel::ApplyIPDIonizationKernel<
91+
particles::atomicPhysics::spawnFromSourceSpeciesModules::SkipFinishedSuperCellsAtomicPhysics,
92+
T_IPDModel,
93+
std::integral_constant<bool, AtomicDataType::switchFieldIonization>>())
94+
.config(mapper.getGridDim(), ions)(
95+
mapper,
96+
idProvider->getDeviceGenerator(),
97+
ions.getDeviceParticlesBox(),
98+
electrons.getDeviceParticlesBox(),
99+
timeRemainingField.getDeviceDataBox(),
100+
foundUnboundIonField.getDeviceDataBox(),
101+
atomicData.template getChargeStateDataDataBox</*on device*/ false>(),
102+
atomicData.template getAtomicStateDataDataBox</*on device*/ false>(),
103+
atomicData.template getIPDIonizationStateDataBox</*on device*/ false>(),
104+
fieldE.getDeviceDataBox(),
105+
debyeLengthField.getDeviceDataBox(),
106+
temperatureEnergyField.getDeviceDataBox(),
107+
zStarField.getDeviceDataBox());
108+
}
109+
else
110+
{
111+
PMACC_LOCKSTEP_KERNEL(
112+
s_IPD::kernel::ApplyIPDIonizationKernel<
113+
particles::atomicPhysics::spawnFromSourceSpeciesModules::NeverSkipSuperCells,
114+
T_IPDModel,
115+
std::integral_constant<bool, AtomicDataType::switchFieldIonization>>())
116+
.config(mapper.getGridDim(), ions)(
117+
mapper,
118+
idProvider->getDeviceGenerator(),
119+
ions.getDeviceParticlesBox(),
120+
electrons.getDeviceParticlesBox(),
121+
timeRemainingField.getDeviceDataBox(),
122+
foundUnboundIonField.getDeviceDataBox(),
123+
atomicData.template getChargeStateDataDataBox</*on device*/ false>(),
124+
atomicData.template getAtomicStateDataDataBox</*on device*/ false>(),
125+
atomicData.template getIPDIonizationStateDataBox</*on device*/ false>(),
126+
fieldE.getDeviceDataBox(),
127+
debyeLengthField.getDeviceDataBox(),
128+
temperatureEnergyField.getDeviceDataBox(),
129+
zStarField.getDeviceDataBox());
130+
}
102131

103132
// no need to call fillAllGaps, since we do not leave any gaps
104133

‎include/picongpu/particles/creation/SpawnFromSourceSpecies.kernel‎

Lines changed: 1 addition & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -113,7 +113,6 @@ namespace picongpu::particles::creation
113113
using SourceSpeciesFramePtr = typename T_SourceParticleBox::FramePtr;
114114
using ProductSpeciesFramePtr = typename T_ProductParticleBox::FramePtr;
115115

116-
// get frame size
117116
PMACC_CASSERT_MSG(
118117
sourceSpecies_framesize_must_less_or_equal_productSpecies_framesize,
119118
T_SourceParticleBox::frameSize == T_ProductParticleBox::frameSize);
@@ -123,7 +122,7 @@ namespace picongpu::particles::creation
123122
pmacc::DataSpace<picongpu::simDim> const superCellIndex
124123
= areaMapping.getSuperCellIndex(worker.blockDomIdxND());
125124

126-
// renaming necessary for compiler unknown reasons
125+
// renaming necessary for compiler, unknown reasons
127126
using IndexFunctor = typename Modules::template AdditionalDataIndexFunctor<T_KernelConfigOptions...>;
128127
auto const additionalDataIndex = IndexFunctor::getIndex(areaMapping, superCellIndex);
129128

@@ -139,7 +138,6 @@ namespace picongpu::particles::creation
139138
std::forward<T_AdditionalData>(additonalData)...);
140139
bool const noParticlesInSuperCell = !sourceSpeciesFrame.isValid();
141140

142-
// test for skipping superCell due to no sourceSpecies particles or SkipSuperCell == true
143141
if(skipSuperCell || noParticlesInSuperCell)
144142
return;
145143

‎include/picongpu/simulation/stage/AtomicPhysics.x.cpp‎

Lines changed: 13 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -550,7 +550,7 @@ namespace picongpu::simulation::stage
550550
mappingDesc);
551551
}
552552

553-
template<typename T_DeviceReduce>
553+
template<bool T_SkipFinishedSuperCell, typename T_DeviceReduce>
554554
HINLINE static void applyIPDIonization(
555555
picongpu::MappingDesc const& mappingDesc,
556556
uint32_t const currentStep,
@@ -562,15 +562,15 @@ namespace picongpu::simulation::stage
562562
DataSpace<picongpu::simDim> const fieldGridLayoutFoundUnbound
563563
= foundUnboundIonField.getGridLayout().sizeWithoutGuardND();
564564

565-
// pressure ionization loop, ends when no ion in unbound state anymore
565+
// ipd ionization loop, ends when no ion is in unbound state anymore
566566
bool foundUnbound = true;
567567
while(foundUnbound)
568568
{
569569
resetFoundUnboundIon(foundUnboundIonField);
570570
calculateIPDInput(mappingDesc, currentStep);
571-
picongpu::atomicPhysics::IPDModel::template applyIPDIonization<AtomicPhysicsIonSpecies>(
572-
mappingDesc,
573-
currentStep);
571+
picongpu::atomicPhysics::IPDModel::template applyIPDIonization<
572+
AtomicPhysicsIonSpecies,
573+
T_SkipFinishedSuperCell>(mappingDesc, currentStep);
574574

575575
auto linearizedFoundUnboundIonBox = LinearizedBox<FoundUnboundField>(
576576
foundUnboundIonField.getDeviceDataBox(),
@@ -700,9 +700,15 @@ namespace picongpu::simulation::stage
700700
while(!isSubSteppingComplete)
701701
{
702702
debugForceConstantElectronTemperature(currentStep);
703-
applyIPDIonization(mappingDesc, currentStep, deviceLocalReduce);
703+
applyIPDIonization</*skip finished super cells*/ true>(
704+
mappingDesc,
705+
currentStep,
706+
deviceLocalReduce);
707+
704708
resetAcceptedStatus(mappingDesc);
705709
applyInstantTransitions(mappingDesc, currentStep, deviceLocalReduce);
710+
resetAcceptedStatus(mappingDesc);
711+
706712
resetElectronEnergyHistogram();
707713
binElectronsToEnergyHistogram(mappingDesc);
708714
resetTimeStep(mappingDesc);
@@ -752,7 +758,7 @@ namespace picongpu::simulation::stage
752758
} // end atomicPhysics sub-stepping loop
753759

754760
// ensure no unbound states are visible to the rest of the loop
755-
applyIPDIonization(mappingDesc, currentStep, deviceLocalReduce);
761+
applyIPDIonization</*skip finished super cells*/ false>(mappingDesc, currentStep, deviceLocalReduce);
756762
}
757763
};
758764

0 commit comments

Comments
 (0)