Geant4 11.4.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
MCGIDI::ProtareComposite Class Reference

#include <MCGIDI.hpp>

Inheritance diagram for MCGIDI::ProtareComposite:

Public Member Functions

LUPI_HOST_DEVICE ProtareComposite ()
LUPI_HOST ProtareComposite (LUPI::StatusMessageReporting &a_smr, GIDI::ProtareComposite const &a_protare, PoPI::Database const &a_pops, Transporting::MC &a_settings, GIDI::Transporting::Particles const &a_particles, DomainHash const &a_domainHash, GIDI::Styles::TemperatureInfos const &a_temperatureInfos, GIDI::ExcludeReactionsSet const &a_reactionsToExclude, std::size_t a_reactionsToExcludeOffset=0, bool a_allowFixedGrid=true)
LUPI_HOST_DEVICE ~ProtareComposite ()
Vector< ProtareSingle * > protares () const
LUPI_HOST void setUserParticleIndex2 (int a_particleIndex, int a_userParticleIndex)
LUPI_HOST void setUserParticleIndexViaIntid2 (int a_particleIntid, int a_userParticleIndex)
LUPI_HOST_DEVICE std::size_t numberOfProtares () const
LUPI_HOST_DEVICE ProtareSingle const * protare (std::size_t a_index) const
LUPI_HOST_DEVICE ProtareSingleprotare (std::size_t a_index)
LUPI_HOST_DEVICE ProtareSingle const * protareWithReaction (std::size_t a_index) const
LUPI_HOST_DEVICE double minimumEnergy () const
LUPI_HOST_DEVICE double maximumEnergy () const
LUPI_HOST_DEVICE Vector< double > temperatures (std::size_t a_index=0) const
LUPI_HOST Vector< double > const & projectileMultiGroupBoundaries () const
LUPI_HOST Vector< double > const & projectileMultiGroupBoundariesCollapsed () const
LUPI_HOST_DEVICE std::size_t numberOfReactions () const
LUPI_HOST_DEVICE Reaction const * reaction (std::size_t a_index) const
LUPI_HOST_DEVICE std::size_t numberOfOrphanProducts () const
LUPI_HOST_DEVICE Reaction const * orphanProduct (std::size_t a_index) const
LUPI_HOST_DEVICE bool hasFission () const
LUPI_HOST_DEVICE bool hasIncoherentDoppler () const
LUPI_HOST_DEVICE int URR_index () const
LUPI_HOST_DEVICE bool hasURR_probabilityTables () const
LUPI_HOST_DEVICE double URR_domainMin () const
LUPI_HOST_DEVICE double URR_domainMax () const
LUPI_HOST_DEVICE bool reactionHasURR_probabilityTables (std::size_t a_index) const
LUPI_HOST_DEVICE double threshold (std::size_t a_index) const
LUPI_HOST_DEVICE double crossSection (URR_protareInfos const &a_URR_protareInfos, std::size_t a_hashIndex, double a_temperature, double a_energy, bool a_sampling=false) const
LUPI_HOST_DEVICE void crossSectionVector (double a_temperature, double a_userFactor, std::size_t a_numberAllocated, double *a_crossSectionVector) const
LUPI_HOST_DEVICE double reactionCrossSection (std::size_t a_reactionIndex, URR_protareInfos const &a_URR_protareInfos, std::size_t a_hashIndex, double a_temperature, double a_energy, bool a_sampling=false) const
LUPI_HOST_DEVICE double reactionCrossSection (std::size_t a_reactionIndex, URR_protareInfos const &a_URR_protareInfos, double a_temperature, double a_energy) const
template<typename RNG>
LUPI_HOST_DEVICE std::size_t sampleReaction (Sampling::Input &a_input, URR_protareInfos const &a_URR_protareInfos, std::size_t a_hashIndex, double a_crossSection, RNG &&a_rng) const
LUPI_HOST_DEVICE double depositionEnergy (std::size_t a_hashIndex, double a_temperature, double a_energy) const
LUPI_HOST_DEVICE double depositionMomentum (std::size_t a_hashIndex, double a_temperature, double a_energy) const
LUPI_HOST_DEVICE double productionEnergy (std::size_t a_hashIndex, double a_temperature, double a_energy) const
LUPI_HOST_DEVICE double gain (std::size_t a_hashIndex, double a_temperature, double a_energy, int a_particleIndex) const
LUPI_HOST_DEVICE double gainViaIntid (std::size_t a_hashIndex, double a_temperature, double a_energy, int a_particleIntid) const
LUPI_HOST_DEVICE Vector< double > const & upscatterModelAGroupVelocities () const
LUPI_HOST_DEVICE void serialize2 (LUPI::DataBuffer &a_buffer, LUPI::DataBuffer::Mode a_mode)
LUPI_HOST_DEVICE long sizeOf2 () const
Public Member Functions inherited from MCGIDI::Protare
LUPI_HOST_DEVICE Protare (ProtareType a_protareType)
LUPI_HOST Protare (ProtareType a_protareType, GIDI::Protare const &a_protare, Transporting::MC const &a_settings, PoPI::Database const &a_pops)
virtual LUPI_HOST_DEVICE ~Protare ()
LUPI_HOST_DEVICE ProtareType protareType () const
LUPI_HOST_DEVICE String const & projectileID () const
LUPI_HOST_DEVICE int projectileIntid () const
LUPI_HOST_DEVICE int projectileIndex () const
LUPI_HOST_DEVICE int projectileUserIndex () const
LUPI_HOST_DEVICE double projectileMass () const
LUPI_HOST_DEVICE double projectileExcitationEnergy () const
LUPI_HOST_DEVICE String const & targetID () const
LUPI_HOST_DEVICE int targetIntid () const
LUPI_HOST_DEVICE int targetIndex () const
LUPI_HOST_DEVICE int targetUserIndex () const
LUPI_HOST_DEVICE double targetMass () const
LUPI_HOST_DEVICE double targetExcitationEnergy () const
LUPI_HOST_DEVICE int photonIndex () const
LUPI_HOST_DEVICE int userPhotonIndex () const
LUPI_HOST_DEVICE String evaluation () const
LUPI_HOST GIDI::Frame projectileFrame () const
LUPI_HOST Vector< int > const & productIntids (bool a_transportablesOnly) const
LUPI_HOST Vector< int > const & productIndices (bool a_transportablesOnly) const
LUPI_HOST Vector< int > const & userProductIndices (bool a_transportablesOnly) const
LUPI_HOST void setUserParticleIndex (int a_particleIndex, int a_userParticleIndex)
LUPI_HOST void setUserParticleIndexViaIntid (int a_particleIntid, int a_userParticleIndex)
LUPI_HOST_DEVICE bool isTNSL_ProtareSingle () const
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE std::size_t numberOfProtares () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE ProtareSingle const * protare (std::size_t a_index) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE ProtareSingleprotare (std::size_t a_index) MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE ProtareSingle const * protareWithReaction (std::size_t a_index) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double minimumEnergy () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double maximumEnergy () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE Vector< double > temperatures (std::size_t a_index=0) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST Vector< double > const & projectileMultiGroupBoundaries () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST Vector< double > const & projectileMultiGroupBoundariesCollapsed () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE std::size_t numberOfReactions () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE Reaction const * reaction (std::size_t a_index) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE std::size_t numberOfOrphanProducts () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE Reaction const * orphanProduct (std::size_t a_index) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE bool hasFission () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE bool hasIncoherentDoppler () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE int URR_index () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE bool hasURR_probabilityTables () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double URR_domainMin () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double URR_domainMax () const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE bool reactionHasURR_probabilityTables (std::size_t a_index) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double threshold (std::size_t a_index) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double crossSection (URR_protareInfos const &a_URR_protareInfos, std::size_t a_hashIndex, double a_temperature, double a_energy, bool a_sampling=false) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE void crossSectionVector (double a_temperature, double a_userFactor, std::size_t a_numberAllocated, double *a_crossSectionVector) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double reactionCrossSection (std::size_t a_reactionIndex, URR_protareInfos const &a_URR_protareInfos, std::size_t a_hashIndex, double a_temperature, double a_energy, bool a_sampling=false) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double reactionCrossSection (std::size_t a_reactionIndex, URR_protareInfos const &a_URR_protareInfos, double a_temperature, double a_energy) const MCGIDI_TRUE_VIRTUAL
template<typename RNG>
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE std::size_t sampleReaction (Sampling::Input &a_input, URR_protareInfos const &a_URR_protareInfos, std::size_t a_hashIndex, double a_crossSection, RNG &&a_rng) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double depositionEnergy (std::size_t a_hashIndex, double a_temperature, double a_energy) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double depositionMomentum (std::size_t a_hashIndex, double a_temperature, double a_energy) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double productionEnergy (std::size_t a_hashIndex, double a_temperature, double a_energy) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double gain (std::size_t a_hashIndex, double a_temperature, double a_energy, int a_particleIndex) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE double gainViaIntid (std::size_t a_hashIndex, double a_temperature, double a_energy, int a_particleIntid) const MCGIDI_TRUE_VIRTUAL
MCGIDI_VIRTUAL_FUNCTION LUPI_HOST_DEVICE Vector< double > const & upscatterModelAGroupVelocities () const MCGIDI_TRUE_VIRTUAL
LUPI_HOST_DEVICE void serialize (LUPI::DataBuffer &a_buffer, LUPI::DataBuffer::Mode a_mode)
LUPI_HOST_DEVICE void serialize2 (LUPI::DataBuffer &a_buffer, LUPI::DataBuffer::Mode a_mode)
LUPI_HOST_DEVICE void serializeCommon (LUPI::DataBuffer &a_buffer, LUPI::DataBuffer::Mode a_mode)
LUPI_HOST_DEVICE long sizeOf () const
LUPI_HOST_DEVICE long memorySize ()
LUPI_HOST_DEVICE void incrementMemorySize (long &a_totalMemory, long &a_sharedMemory)
template<typename RNG>
LUPI_HOST_DEVICE std::size_t sampleReaction (Sampling::Input &a_input, URR_protareInfos const &a_URR_protareInfos, std::size_t a_hashIndex, double a_crossSection, RNG &&a_rng) const

Additional Inherited Members

Public Attributes inherited from MCGIDI::Protare
friend ProtareSingle
friend ProtareComposite
friend ProtareTNSL

Detailed Description

Class representing a GNDS <reactionSuite> node with only data needed for Monte Carlo transport. The data are also stored in a way that is better suited for Monte Carlo transport. For example, cross section data for each reaction are not stored with its reaction, but within the HeatedCrossSections member of the Protare.

Definition at line 1750 of file MCGIDI.hpp.

Constructor & Destructor Documentation

◆ ProtareComposite() [1/2]

LUPI_HOST_DEVICE MCGIDI::ProtareComposite::ProtareComposite ( )

Default constructor used when broadcasting a Protare as needed by MPI or GPUs.

Definition at line 26 of file MCGIDI_protareComposite.cc.

26 :
28 m_numberOfReactions( 0 ),
29 m_numberOfOrphanProducts( 0 ),
30 m_minimumEnergy( 0.0 ),
31 m_maximumEnergy( 0.0 ) {
32
33}
LUPI_HOST_DEVICE Protare(ProtareType a_protareType)

Referenced by ProtareComposite().

◆ ProtareComposite() [2/2]

LUPI_HOST MCGIDI::ProtareComposite::ProtareComposite ( LUPI::StatusMessageReporting & a_smr,
GIDI::ProtareComposite const & a_protare,
PoPI::Database const & a_pops,
Transporting::MC & a_settings,
GIDI::Transporting::Particles const & a_particles,
DomainHash const & a_domainHash,
GIDI::Styles::TemperatureInfos const & a_temperatureInfos,
GIDI::ExcludeReactionsSet const & a_reactionsToExclude,
std::size_t a_reactionsToExcludeOffset = 0,
bool a_allowFixedGrid = true )

◆ ~ProtareComposite()

LUPI_HOST_DEVICE MCGIDI::ProtareComposite::~ProtareComposite ( )

Definition at line 94 of file MCGIDI_protareComposite.cc.

94 {
95
96 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
97
98 for( std::size_t i1 = 0; i1 < length; ++i1 ) delete m_protares[i1];
99}

Member Function Documentation

◆ crossSection()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::crossSection ( URR_protareInfos const & a_URR_protareInfos,
std::size_t a_hashIndex,
double a_temperature,
double a_energy,
bool a_sampling = false ) const

Returns the total cross section.

Parameters
a_URR_protareInfos[in] URR information.
a_hashIndex[in] The cross section hash index.
a_temperature[in] The target temperature.
a_energy[in] The projectile energy.
a_sampling[in] Only used for multi-group cross sections. When sampling, the cross section in the group where threshold is present the cross section is augmented.
Returns
The total cross section.

Definition at line 406 of file MCGIDI_protareComposite.cc.

406 {
407
408 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
409 double cross_section = 0.0;
410
411 for( std::size_t i1 = 0; i1 < length; ++i1 ) cross_section += m_protares[i1]->crossSection( a_URR_protareInfos, a_hashIndex, a_temperature, a_energy, a_sampling );
412
413 return( cross_section );
414}
LUPI_HOST_DEVICE double crossSection(URR_protareInfos const &a_URR_protareInfos, std::size_t a_hashIndex, double a_temperature, double a_energy, bool a_sampling=false) const

Referenced by crossSection().

◆ crossSectionVector()

LUPI_HOST_DEVICE void MCGIDI::ProtareComposite::crossSectionVector ( double a_temperature,
double a_userFactor,
std::size_t a_numberAllocated,
double * a_crossSectionVector ) const

Adds the energy dependent, total cross section corresponding to the temperature a_temperature multiplied by a_userFactor to a_crossSectionVector.

Parameters
a_temperature[in] Specifies the temperature of the material.
a_userFactor[in] User factor which all cross sections are multiplied by.
a_numberAllocated[in] The length of memory allocated for a_crossSectionVector.
a_crossSectionVector[in/out] The energy dependent, total cross section to add cross section data to.

Definition at line 425 of file MCGIDI_protareComposite.cc.

425 {
426
427 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
428
429 for( std::size_t i1 = 0; i1 < length; ++i1 ) m_protares[i1]->crossSectionVector( a_temperature, a_userFactor, a_numberAllocated, a_crossSectionVector );
430}
LUPI_HOST_DEVICE void crossSectionVector(double a_temperature, double a_userFactor, std::size_t a_numberAllocated, double *a_crossSectionVector) const

Referenced by crossSectionVector().

◆ depositionEnergy()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::depositionEnergy ( std::size_t a_hashIndex,
double a_temperature,
double a_energy ) const

Returns the total deposition energy.

Parameters
a_hashIndex[in] The cross section hash index.
a_temperature[in] The target temperature.
a_energy[in] The projectile energy.
Returns
The total deposition energy.

Definition at line 503 of file MCGIDI_protareComposite.cc.

503 {
504
505 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
506 double deposition_energy = 0.0;
507
508 for( std::size_t i1 = 0; i1 < length; ++i1 ) deposition_energy += m_protares[i1]->depositionEnergy( a_hashIndex, a_temperature, a_energy );
509
510 return( deposition_energy );
511}
LUPI_HOST_DEVICE double depositionEnergy(std::size_t a_hashIndex, double a_temperature, double a_energy) const

Referenced by depositionEnergy().

◆ depositionMomentum()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::depositionMomentum ( std::size_t a_hashIndex,
double a_temperature,
double a_energy ) const

Returns the total deposition momentum.

Parameters
a_hashIndex[in] The cross section hash index.
a_temperature[in] The target temperature.
a_energy[in] The projectile energy.
Returns
The total deposition momentum.

Definition at line 523 of file MCGIDI_protareComposite.cc.

523 {
524
525 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
526 double deposition_momentum = 0.0;
527
528 for( std::size_t i1 = 0; i1 < length; ++i1 ) deposition_momentum += m_protares[i1]->depositionMomentum( a_hashIndex, a_temperature, a_energy );
529
530 return( deposition_momentum );
531}
LUPI_HOST_DEVICE double depositionMomentum(std::size_t a_hashIndex, double a_temperature, double a_energy) const

Referenced by depositionMomentum().

◆ gain()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::gain ( std::size_t a_hashIndex,
double a_temperature,
double a_energy,
int a_particleIndex ) const

Returns the gain for particle with index a_particleIndex.

Parameters
a_hashIndex[in] The continuous energy hash or multi-group index.
a_temperature[in] The temperature of the target.
a_energy[in] The projectile energy.
a_particleIndex[in] The index of the particle whose gain is to be returned.
Returns
[in] A vector of the length of the number of multi-group groups.

Definition at line 564 of file MCGIDI_protareComposite.cc.

564 {
565
566 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
567 double gain1 = m_protares[0]->gain( a_hashIndex, a_temperature, a_energy, a_particleIndex );
568
569 for( std::size_t i1 = 1; i1 < length; ++i1 ) gain1 += m_protares[i1]->gain( a_hashIndex, a_temperature, a_energy, a_particleIndex );
570
571 return( gain1 );
572}
LUPI_HOST_DEVICE double gain(std::size_t a_hashIndex, double a_temperature, double a_energy, int a_particleIndex) const

Referenced by gain().

◆ gainViaIntid()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::gainViaIntid ( std::size_t a_hashIndex,
double a_temperature,
double a_energy,
int a_particleIntid ) const

Returns the gain for particle with intid a_particleIntid.

Parameters
a_hashIndex[in] The continuous energy hash or multi-group index.
a_temperature[in] The temperature of the target.
a_energy[in] The projectile energy.
a_particleIntid[in] The intid of the particle whose gain is to be returned.
Returns
[in] A vector of the length of the number of multi-group groups.

Definition at line 585 of file MCGIDI_protareComposite.cc.

585 {
586
587 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
588 double gain1 = m_protares[0]->gainViaIntid( a_hashIndex, a_temperature, a_energy, a_particleIntid );
589
590 for( std::size_t i1 = 1; i1 < length; ++i1 ) gain1 += m_protares[i1]->gainViaIntid( a_hashIndex, a_temperature, a_energy, a_particleIntid );
591
592 return( gain1 );
593}
LUPI_HOST_DEVICE double gainViaIntid(std::size_t a_hashIndex, double a_temperature, double a_energy, int a_particleIntid) const

Referenced by gainViaIntid().

◆ hasFission()

LUPI_HOST_DEVICE bool MCGIDI::ProtareComposite::hasFission ( ) const

Returns true is one of the protares has a fission channel and false otherwise.

Returns
true is one of the protares has a fission channel and false otherwise.

Definition at line 260 of file MCGIDI_protareComposite.cc.

260 {
261
262 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
263
264 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
265 if( m_protares[i1]->hasFission( ) ) return( true );
266 }
267
268 return( false );
269}
LUPI_HOST_DEVICE bool hasFission() const

Referenced by hasFission().

◆ hasIncoherentDoppler()

LUPI_HOST_DEVICE bool MCGIDI::ProtareComposite::hasIncoherentDoppler ( ) const

Returns true if this has a photoatomic incoherent doppler broadened reaction and false otherwise.

Returns
true is one of the protares has a fission channel and false otherwise.

Definition at line 277 of file MCGIDI_protareComposite.cc.

277 {
278
279 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
280
281 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
282 if( m_protares[i1]->hasIncoherentDoppler( ) ) return( true );
283 }
284
285 return( false );
286}
LUPI_HOST_DEVICE bool hasIncoherentDoppler() const

Referenced by hasIncoherentDoppler().

◆ hasURR_probabilityTables()

LUPI_HOST_DEVICE bool MCGIDI::ProtareComposite::hasURR_probabilityTables ( ) const

Returns true if one the the sub-protares of this has a unresolved resonance region (URR) data and false otherwise.

Returns
true is if this has a URR data.

Definition at line 294 of file MCGIDI_protareComposite.cc.

294 {
295
296 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
297
298 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
299 if( m_protares[i1]->hasURR_probabilityTables( ) ) return( true );
300 }
301
302 return( false );
303}
LUPI_HOST_DEVICE bool hasURR_probabilityTables() const

Referenced by hasURR_probabilityTables(), URR_domainMax(), and URR_domainMin().

◆ maximumEnergy()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::maximumEnergy ( ) const
inline

Returns the value of the m_maximumEnergy member.

Definition at line 1778 of file MCGIDI.hpp.

◆ minimumEnergy()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::minimumEnergy ( ) const
inline

Returns the value of the m_minimumEnergy member.

Definition at line 1777 of file MCGIDI.hpp.

◆ numberOfOrphanProducts()

LUPI_HOST_DEVICE std::size_t MCGIDI::ProtareComposite::numberOfOrphanProducts ( ) const
inline

Returns the value of the m_numberOfOrphanProducts member.

Definition at line 1789 of file MCGIDI.hpp.

◆ numberOfProtares()

LUPI_HOST_DEVICE std::size_t MCGIDI::ProtareComposite::numberOfProtares ( ) const
inline

Returns the number of protares contained in this.

Definition at line 1772 of file MCGIDI.hpp.

◆ numberOfReactions()

LUPI_HOST_DEVICE std::size_t MCGIDI::ProtareComposite::numberOfReactions ( ) const
inline

Returns the value of the m_numberOfReactions member.

Definition at line 1786 of file MCGIDI.hpp.

Referenced by orphanProduct(), protareWithReaction(), reaction(), reactionCrossSection(), reactionCrossSection(), reactionHasURR_probabilityTables(), and threshold().

◆ orphanProduct()

LUPI_HOST_DEVICE Reaction const * MCGIDI::ProtareComposite::orphanProduct ( std::size_t a_index) const

Returns the reaction at index a_index. If a_index is negative, the reaction of the TNSL protare at index -*a_index* is returned; otherwise, the reaction from the regular protare at index a_index is returned.

Parameters
a_index[in] The index of the reaction to return.
Returns
The reaction at index a_index.

Definition at line 240 of file MCGIDI_protareComposite.cc.

240 {
241
242 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
243
244 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
245 std::size_t numberOfReactions = m_protares[i1]->numberOfOrphanProducts( );
246
247 if( a_index < numberOfReactions ) return( m_protares[i1]->orphanProduct( a_index ) );
248 a_index -= numberOfReactions;
249 }
250
251 return( nullptr );
252}
LUPI_HOST_DEVICE Reaction const * orphanProduct(std::size_t a_index) const
LUPI_HOST_DEVICE std::size_t numberOfReactions() const
Definition MCGIDI.hpp:1786

Referenced by orphanProduct().

◆ productionEnergy()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::productionEnergy ( std::size_t a_hashIndex,
double a_temperature,
double a_energy ) const

Returns the total production energy.

Parameters
a_hashIndex[in] The cross section hash index.
a_temperature[in] The target temperature.
a_energy[in] The projectile energy.
Returns
The total production energy.

Definition at line 543 of file MCGIDI_protareComposite.cc.

543 {
544
545 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
546 double production_energy = 0.0;
547
548 for( std::size_t i1 = 0; i1 < length; ++i1 ) production_energy += m_protares[i1]->productionEnergy( a_hashIndex, a_temperature, a_energy );
549
550 return( production_energy );
551}
LUPI_HOST_DEVICE double productionEnergy(std::size_t a_hashIndex, double a_temperature, double a_energy) const

Referenced by productionEnergy().

◆ projectileMultiGroupBoundaries()

LUPI_HOST Vector< double > const & MCGIDI::ProtareComposite::projectileMultiGroupBoundaries ( ) const
inline

Returns the value of the m_projectileMultiGroupBoundaries member.

Definition at line 1781 of file MCGIDI.hpp.

Referenced by projectileMultiGroupBoundaries().

◆ projectileMultiGroupBoundariesCollapsed()

LUPI_HOST Vector< double > const & MCGIDI::ProtareComposite::projectileMultiGroupBoundariesCollapsed ( ) const
inline

Returns the value of the m_projectileMultiGroupBoundariesCollapsed member.

Definition at line 1783 of file MCGIDI.hpp.

Referenced by projectileMultiGroupBoundariesCollapsed().

◆ protare() [1/2]

LUPI_HOST_DEVICE ProtareSingle * MCGIDI::ProtareComposite::protare ( std::size_t a_index)

Returns the pointer representing the (a_index - 1)th ProtareSingle.

Parameters
a_index[in] Index of the ProtareSingle to return.
Returns
Pointer to the requested protare or nullptr if invalid a_index..

Definition at line 153 of file MCGIDI_protareComposite.cc.

153 {
154
155 for( std::size_t i1 = 0; i1 < m_protares.size( ); ++i1 ) {
156 std::size_t number = m_protares[i1]->numberOfProtares( );
157
158 if( number > a_index ) return( m_protares[i1]->protare( a_index ) );
159 a_index -= number;
160 }
161
162 return( nullptr );
163}
LUPI_HOST_DEVICE ProtareSingle const * protare(std::size_t a_index) const

◆ protare() [2/2]

LUPI_HOST_DEVICE ProtareSingle const * MCGIDI::ProtareComposite::protare ( std::size_t a_index) const

Returns the pointer representing the (a_index - 1)th ProtareSingle.

Parameters
a_index[in] Index of the ProtareSingle to return.
Returns
Pointer to the requested protare or nullptr if invalid a_index..

Definition at line 133 of file MCGIDI_protareComposite.cc.

133 {
134
135 for( std::size_t i1 = 0; i1 < m_protares.size( ); ++i1 ) {
136 std::size_t number = m_protares[i1]->numberOfProtares( );
137
138 if( number > a_index ) return( m_protares[i1]->protare( a_index ) );
139 a_index -= number;
140 }
141
142 return( nullptr );
143}

Referenced by protare(), and protare().

◆ protares()

Vector< ProtareSingle * > MCGIDI::ProtareComposite::protares ( ) const
inline

Returns the value of the m_protares member.

Definition at line 1766 of file MCGIDI.hpp.

◆ protareWithReaction()

LUPI_HOST_DEVICE ProtareSingle const * MCGIDI::ProtareComposite::protareWithReaction ( std::size_t a_index) const

Returns the pointer to the ProtareSingle that contains the (a_index - 1)th reaction.

Parameters
a_index[in] Index of the reaction.
Returns
Pointer to the requested protare or nullptr if invalid a_index..

Definition at line 173 of file MCGIDI_protareComposite.cc.

173 {
174
175 for( std::size_t i1 = 0; i1 < m_protares.size( ); ++i1 ) {
176 std::size_t numberOfReactions = m_protares[i1]->numberOfReactions( );
177
178 if( a_index < numberOfReactions ) return( m_protares[i1] );
179 a_index -= numberOfReactions;
180 }
181
182 return( nullptr );
183}

◆ reaction()

LUPI_HOST_DEVICE Reaction const * MCGIDI::ProtareComposite::reaction ( std::size_t a_index) const

Returns the reaction at index a_index. If a_index is negative, the reaction of the TNSL protare at index -*a_index* is returned; otherwise, the reaction from the regular protare at index a_index is returned.

Parameters
a_index[in] The index of the reaction to return.
Returns
The reaction at index a_index.

Definition at line 217 of file MCGIDI_protareComposite.cc.

217 {
218
219 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
220
221 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
222 std::size_t numberOfReactions = m_protares[i1]->numberOfReactions( );
223
224 if( a_index < numberOfReactions ) return( m_protares[i1]->reaction( a_index ) );
225 a_index -= numberOfReactions;
226 }
227
228 return( nullptr );
229}
LUPI_HOST_DEVICE Reaction const * reaction(std::size_t a_index) const

Referenced by reaction().

◆ reactionCrossSection() [1/2]

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::reactionCrossSection ( std::size_t a_reactionIndex,
URR_protareInfos const & a_URR_protareInfos,
double a_temperature,
double a_energy ) const

Returns the cross section for reaction at index a_reactionIndex.

Parameters
a_reactionIndex[in] The index of the reaction.
a_URR_protareInfos[in] URR information.
a_temperature[in] The target temperature.
a_energy[in] The projectile energy.
Returns
The total cross section.

Definition at line 475 of file MCGIDI_protareComposite.cc.

475 {
476
477 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
478 double cross_section = 0.0;
479
480 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
481 std::size_t numberOfReactions = m_protares[i1]->numberOfReactions( );
482
483 if( a_reactionIndex < numberOfReactions ) {
484 cross_section = m_protares[i1]->reactionCrossSection( a_reactionIndex, a_URR_protareInfos, a_temperature, a_energy );
485 break;
486 }
487 a_reactionIndex -= numberOfReactions;
488 }
489
490 return( cross_section );
491}

◆ reactionCrossSection() [2/2]

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::reactionCrossSection ( std::size_t a_reactionIndex,
URR_protareInfos const & a_URR_protareInfos,
std::size_t a_hashIndex,
double a_temperature,
double a_energy,
bool a_sampling = false ) const

Returns the cross section for reaction at index a_reactionIndex.

Parameters
a_reactionIndex[in] The index of the reaction.
a_URR_protareInfos[in] URR information.
a_hashIndex[in] The cross section hash index.
a_temperature[in] The target temperature.
a_energy[in] The projectile energy.
a_sampling[in] Only used for multi-group cross sections. When sampling, the cross section in the group where threshold is present the cross section is augmented.
Returns
The total cross section.

Definition at line 446 of file MCGIDI_protareComposite.cc.

446 {
447
448 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
449 double cross_section = 0.0;
450
451 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
452 std::size_t numberOfReactions = m_protares[i1]->numberOfReactions( );
453
454 if( a_reactionIndex < numberOfReactions ) {
455 cross_section = m_protares[i1]->reactionCrossSection( a_reactionIndex, a_URR_protareInfos, a_hashIndex, a_temperature, a_energy, a_sampling );
456 break;
457 }
458 a_reactionIndex -= numberOfReactions;
459 }
460
461 return( cross_section );
462}

◆ reactionHasURR_probabilityTables()

LUPI_HOST_DEVICE bool MCGIDI::ProtareComposite::reactionHasURR_probabilityTables ( std::size_t a_index) const

Returns true if the reaction at index a_index has URR robability tables and false otherwise.

Parameters
a_index[in] The index of the reaction.
Returns
true if the reaction has URR robability tables and false otherwise.

Definition at line 356 of file MCGIDI_protareComposite.cc.

356 {
357
358 std::size_t length = m_protares.size( );
359
360 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
361 std::size_t numberOfReactions = m_protares[i1]->numberOfReactions( );
362
363 if( a_index < numberOfReactions ) return( m_protares[i1]->reactionHasURR_probabilityTables( a_index ) );
364 a_index -= numberOfReactions;
365 }
366
367 return( false );
368}
LUPI_HOST_DEVICE bool reactionHasURR_probabilityTables(std::size_t a_index) const

Referenced by reactionHasURR_probabilityTables().

◆ sampleReaction()

template<typename RNG>
LUPI_HOST_DEVICE std::size_t MCGIDI::ProtareComposite::sampleReaction ( Sampling::Input & a_input,
URR_protareInfos const & a_URR_protareInfos,
std::size_t a_hashIndex,
double a_crossSection,
RNG && a_rng ) const
inline

Samples a reaction of this and returns its index.

Parameters
a_input[in] Sample options requested by user.
a_URR_protareInfos[in] URR information.
a_hashIndex[in] Specifies the continuous energy hash index or multi-group index.
a_crossSection[in] The total cross section for the protare at a_input.temperature() and a_input.energy().
a_rng[in] The random number generator function that returns a double in the range [0, 1.0).
Returns
The index of the sampled reaction.

Definition at line 2888 of file MCGIDI_headerSource.hpp.

2889 {
2890
2891 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
2892 std::size_t reaction_index = 0;
2893 double cross_section_sum = 0.0;
2894 double cross_section_rng = a_rng( ) * a_crossSection;
2895
2896 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
2897 double cross_section = m_protares[i1]->crossSection( a_URR_protareInfos, a_hashIndex, a_input.temperature( ), a_input.energy( ), true );
2898
2899 cross_section_sum += cross_section;
2900 if( cross_section_sum > cross_section_rng ) {
2901 std::size_t reaction_index2 = m_protares[i1]->sampleReaction( a_input, a_URR_protareInfos, a_hashIndex, cross_section, a_rng );
2902
2903 reaction_index += reaction_index2;
2904 if( reaction_index2 == MCGIDI_nullReaction ) reaction_index = MCGIDI_nullReaction;
2905 break;
2906 }
2907 reaction_index += m_protares[i1]->numberOfReactions( );
2908 }
2909
2910 return( reaction_index );
2911}
#define MCGIDI_nullReaction
Definition MCGIDI.hpp:64

◆ serialize2()

LUPI_HOST_DEVICE void MCGIDI::ProtareComposite::serialize2 ( LUPI::DataBuffer & a_buffer,
LUPI::DataBuffer::Mode a_mode )

This method serializes this for broadcasting as needed for MPI and GPUs. The method can count the number of required bytes, pack this or unpack this depending on a_mode.

Parameters
a_buffer[in] The buffer to read or write data to depending on a_mode.
a_mode[in] Specifies the action of this method.

Definition at line 603 of file MCGIDI_protareComposite.cc.

603 {
604
605 std::size_t vectorSize = m_protares.size( );
606 int vectorSizeInt = static_cast<int>( vectorSize );
607 LUPI::DataBuffer *workingBuffer = &a_buffer;
608
609 DATA_MEMBER_SIZE_T( m_numberOfReactions, a_buffer, a_mode );
610 DATA_MEMBER_SIZE_T( m_numberOfOrphanProducts, a_buffer, a_mode );
611 DATA_MEMBER_DOUBLE( m_minimumEnergy, a_buffer, a_mode );
612 DATA_MEMBER_DOUBLE( m_maximumEnergy, a_buffer, a_mode );
613
614 DATA_MEMBER_INT( vectorSizeInt, *workingBuffer, a_mode );
615 vectorSize = static_cast<std::size_t>( vectorSizeInt );
616
617 if( a_mode == LUPI::DataBuffer::Mode::Unpack ) {
618 m_protares.resize( vectorSize, &(workingBuffer->m_placement) );
619 for( std::size_t vectorIndex = 0; vectorIndex < vectorSize; ++vectorIndex ) {
620 if( workingBuffer->m_placement != nullptr ) {
621 m_protares[vectorIndex] = new(workingBuffer->m_placement) ProtareSingle;
622 workingBuffer->incrementPlacement( sizeof( ProtareSingle ) ); }
623 else {
624 m_protares[vectorIndex] = new ProtareSingle;
625 }
626 }
627 }
628 if( a_mode == LUPI::DataBuffer::Mode::Memory ) {
629 a_buffer.m_placement += m_protares.internalSize();
630 a_buffer.incrementPlacement( sizeof( ProtareSingle ) * vectorSize );
631 }
632
633 for( std::size_t i1 = 0; i1 < vectorSize; ++i1 ) {
634 m_protares[i1]->serializeCommon( a_buffer, a_mode );
635 m_protares[i1]->serialize2( a_buffer, a_mode );
636 }
637}
#define DATA_MEMBER_SIZE_T(member, buf, mode)
#define DATA_MEMBER_DOUBLE(member, buf, mode)
#define DATA_MEMBER_INT( member, buf, mode)
LUPI_HOST_DEVICE void incrementPlacement(std::size_t a_delta)
friend ProtareSingle
Definition MCGIDI.hpp:1599

◆ setUserParticleIndex2()

LUPI_HOST void MCGIDI::ProtareComposite::setUserParticleIndex2 ( int a_particleIndex,
int a_userParticleIndex )

Updates the m_userParticleIndex to a_userParticleIndex for all particles with PoPs index a_particleIndex.

Parameters
a_particleIndex[in] The PoPs index of the particle whose user index is to be set.
a_userParticleIndex[in] The particle index specified by the user.

Definition at line 108 of file MCGIDI_protareComposite.cc.

108 {
109
110 for( auto iter = m_protares.begin( ); iter != m_protares.end( ); ++iter ) (*iter)->setUserParticleIndex( a_particleIndex, a_userParticleIndex );
111}

◆ setUserParticleIndexViaIntid2()

LUPI_HOST void MCGIDI::ProtareComposite::setUserParticleIndexViaIntid2 ( int a_particleIntid,
int a_userParticleIndex )

Updates the m_userParticleIndex to a_userParticleIndex for all particles with PoPs intid a_particleIntid.

Parameters
a_particleIntid[in] The PoPs intid of the particle whose user index is to be set.
a_userParticleIndex[in] The particle index specified by the user.

Definition at line 120 of file MCGIDI_protareComposite.cc.

120 {
121
122 for( auto iter = m_protares.begin( ); iter != m_protares.end( ); ++iter ) (*iter)->setUserParticleIndexViaIntid( a_particleIntid, a_userParticleIndex );
123}

◆ sizeOf2()

LUPI_HOST_DEVICE long MCGIDI::ProtareComposite::sizeOf2 ( ) const
inline

Definition at line 1825 of file MCGIDI.hpp.

1825{ return sizeof(*this); }

◆ temperatures()

LUPI_HOST_DEVICE Vector< double > MCGIDI::ProtareComposite::temperatures ( std::size_t a_index = 0) const

Returns the list of temperatures for the requested ProtareSingle.

Parameters
a_index[in] Index of the reqested ProtareSingle.
Returns
Vector of doubles.

Definition at line 193 of file MCGIDI_protareComposite.cc.

193 {
194
195 for( std::size_t i1 = 0; i1 < m_protares.size( ); ++i1 ) {
196 std::size_t number = m_protares[i1]->numberOfProtares( );
197
198 if( number > a_index ) return( m_protares[i1]->temperatures( a_index ) );
199 a_index -= number;
200 }
201
202 LUPI_THROW( "ProtareSingle::temperatures: a_index not in range." );
203
204 Vector<double> temps; // Only to stop compilers from complaining.
205 return( temps );
206}
#define LUPI_THROW(arg)
LUPI_HOST_DEVICE Vector< double > temperatures(std::size_t a_index=0) const

Referenced by temperatures().

◆ threshold()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::threshold ( std::size_t a_index) const

Returns the threshold for the reaction at index a_index. If a_index is negative, it is set to 0 before the threshold in the regular protare is returned.

Parameters
a_index[in] The index of the reaction.
Returns
The threshold for reaction at index a_index.

Definition at line 379 of file MCGIDI_protareComposite.cc.

379 {
380
381 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
382
383 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
384 std::size_t numberOfReactions = m_protares[i1]->numberOfReactions( );
385
386 if( a_index < numberOfReactions ) return( m_protares[i1]->threshold( a_index ) );
387 a_index -= numberOfReactions;
388 }
389
390 return( 0.0 );
391}
LUPI_HOST_DEVICE double threshold(std::size_t a_index) const

Referenced by threshold().

◆ upscatterModelAGroupVelocities()

LUPI_HOST_DEVICE Vector< double > const & MCGIDI::ProtareComposite::upscatterModelAGroupVelocities ( ) const
inline

Returns a reference to the m_upscatterModelAGroupVelocities member.

Definition at line 1821 of file MCGIDI.hpp.

Referenced by upscatterModelAGroupVelocities().

◆ URR_domainMax()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::URR_domainMax ( ) const

Returns the maximum energy for the unresolved resonance region (URR) domain. If no URR data present, returns -1.

Returns
true is if this has a URR data.

Definition at line 333 of file MCGIDI_protareComposite.cc.

333 {
334
335 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
336 double URR_domain_max = -1.0;
337
338 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
339 if( m_protares[i1]->hasURR_probabilityTables( ) ) {
340 if( URR_domain_max < m_protares[i1]->URR_domainMax( ) ) URR_domain_max = m_protares[i1]->URR_domainMax( );
341
342 }
343 }
344
345 return( URR_domain_max );
346}
LUPI_HOST_DEVICE double URR_domainMax() const

Referenced by URR_domainMax().

◆ URR_domainMin()

LUPI_HOST_DEVICE double MCGIDI::ProtareComposite::URR_domainMin ( ) const

Returns the minimum energy for the unresolved resonance region (URR) domain. If no URR data present, returns -1.

Returns
The energy or -1 if not URR data present.

Definition at line 311 of file MCGIDI_protareComposite.cc.

311 {
312
313 std::size_t length = static_cast<std::size_t>( m_protares.size( ) );
314 double URR_domain_min = 1e32;
315
316 for( std::size_t i1 = 0; i1 < length; ++i1 ) {
317 if( m_protares[i1]->hasURR_probabilityTables( ) ) {
318 if( URR_domain_min > m_protares[i1]->URR_domainMin( ) ) URR_domain_min = m_protares[i1]->URR_domainMin( );
319 }
320 }
321
322 if( URR_domain_min == 1e32 ) URR_domain_min = -1.0;
323
324 return( URR_domain_min );
325}
LUPI_HOST_DEVICE double URR_domainMin() const

Referenced by URR_domainMin().

◆ URR_index()

LUPI_HOST_DEVICE int MCGIDI::ProtareComposite::URR_index ( ) const
inline

Definition at line 1796 of file MCGIDI.hpp.

1796{ return( -1 ); }

The documentation for this class was generated from the following files: