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

#include <GIDI.hpp>

Inheritance diagram for GIDI::ProtareSingle:

Public Member Functions

 ProtareSingle (PoPI::Database const &a_pops, std::string const &a_projectileID, std::string const &a_targetID, std::string const &a_evaluation, std::string const &a_interaction, std::string const &a_formatVersion=GNDS_formatVersion_1_10Chars)
 ProtareSingle (Construction::Settings const &a_construction, std::string const &a_fileName, FileType a_fileType, PoPI::Database const &a_pops, ParticleSubstitution const &a_particleSubstitution, std::vector< std::string > const &a_libraries, std::string const &a_interaction, bool a_targetRequiredInGlobalPoPs=true, bool a_requiredInPoPs=true)
 ProtareSingle (Construction::Settings const &a_construction, HAPI::Node const &a_protare, PoPI::Database const &a_pops, ParticleSubstitution const &a_particleSubstitution, std::vector< std::string > const &a_libraries, std::string const &a_interaction, bool a_targetRequiredInGlobalPoPs=true, bool a_requiredInPoPs=true)
 ~ProtareSingle ()
PoPI::NuclideGammaBranchStateInfos const & nuclideGammaBranchStateInfos () const
HAPI::DataManagerdataManager ()
void setDataManager (HAPI::DataManager *a_dataManager)
void incrementNumberOfLazyParsingHelperForms ()
void incrementNumberOfLazyParsingHelperFormsReplaced ()
TargetInfo::TargetInfo const & targetInfo () const
double projectileEnergyMin () const
double projectileEnergyMax () const
bool isTNSL_ProtareSingle () const
bool isPhotoAtomic () const
Suitereactions ()
Suite const & reactions () const
SuiteorphanProducts ()
Suite const & orphanProducts () const
SuiteincompleteReactions ()
Suite const & incompleteReactions () const
Sums::Sumssums ()
Sums::Sums const & sums () const
SuitefissionComponents ()
bool RutherfordScatteringPresent () const
bool onlyRutherfordScatteringPresent () const
Reaction const * nuclearPlusCoulombInterferenceOnlyReaction () const
Reaction const * checkIf_nuclearPlusCoulombInterferenceWanted (Transporting::MG const &a_settings, Reaction const *a_reaction) const
Reaction const * reactionToMultiGroup (Transporting::MG const &a_settings, std::size_t a_index, ExcludeReactionsSet const &a_reactionsToExclude) const
Reaction const * multiGroupSummedReaction () const
OutputChannel const * multiGroupSummedDelayedNeutrons () const
Suite const & ACE_URR_probabilityTables () const
Suite const & photoAtomicIncoherentDoppler () const
GRIN::GRIN_continuumGammas const * GRIN_continuumGammas2 () const
ProtareType protareType () const
std::size_t numberOfProtares () const
ProtareSingleprotare (std::size_t a_index)
ProtareSingle const * protare (std::size_t a_index) const
LUPI::FormatVersion const & formatVersion (LUPI_maybeUnused std::size_t a_index=0) const
std::string const & fileName (LUPI_maybeUnused std::size_t a_index=0) const
std::string const & realFileName (LUPI_maybeUnused std::size_t a_index=0) const
std::vector< std::string > libraries (LUPI_maybeUnused std::size_t a_index=0) const
std::string const & evaluation (LUPI_maybeUnused std::size_t a_index=0) const
std::string const & interaction (LUPI_maybeUnused std::size_t a_index=0) const
Frame projectileFrame (LUPI_maybeUnused std::size_t a_index=0) const
int numberOfLazyParsingHelperForms () const
int numberOfLazyParsingHelperFormsReplaced () const
double thresholdFactor () const
Documentation_1_10::Suitedocumentations ()
ExternalFile const & externalFile (std::string const a_label) const
ExternalFiles::Suite const & externalFiles () const
Styles::Basestyle (std::string const &a_label)
Styles::Base const & style (std::string const &a_label) const
Styles::Suitestyles ()
Styles::Suite const & styles () const
PoPI::Database const & internalPoPs () const
PoPI::DatabaseinternalPoPs ()
int intid (std::string const &a_id) const
void productIDs (std::set< std::string > &a_ids, Transporting::Particles const &a_particles, bool a_transportablesOnly) const
int maximumLegendreOrder (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, std::string const &a_productID) const
Styles::TemperatureInfos temperatures () const
std::size_t numberOfReactions () const
Reactionreaction (std::size_t a_index)
Reaction const * reaction (std::size_t a_index) const
Reaction const * reaction (std::size_t a_index, Transporting::MG const &a_settings, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
std::size_t numberOfInactiveReactions () const
std::size_t numberOfOrphanProducts () const
ReactionorphanProduct (std::size_t a_index)
Reaction const * orphanProduct (std::size_t a_index) const
std::size_t numberOfIncompleteReactions () const
ReactionincompleteReaction (std::size_t a_index)
Reaction const * incompleteReaction (std::size_t a_index) const
void updateReactionIndices (std::size_t a_offset) const
bool hasFission () const
bool isDelayedFissionNeutronComplete () const
GUPI::AncestryfindInAncestry3 (std::string const &a_item)
GUPI::Ancestry const * findInAncestry3 (std::string const &a_item) const
std::vector< double > groupBoundaries (Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, std::string const &a_productID) const
Vector multiGroupInverseSpeed (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo) const
Vector multiGroupCrossSection (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}, std::string const &a_label="") const
Vector multiGroupQ (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, bool a_final, bool a_effectivePhotoAtomic=true, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupMultiplicity (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, std::string const &a_productID, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupFissionNeutronMultiplicity (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupFissionGammaMultiplicity (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Matrix multiGroupProductMatrix (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, Transporting::Particles const &a_particles, std::string const &a_productID, std::size_t a_order, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Matrix multiGroupFissionMatrix (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, Transporting::Particles const &a_particles, std::size_t a_order, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupTransportCorrection (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, Transporting::Particles const &a_particles, std::size_t a_order, TransportCorrectionType a_transportCorrectionType, double a_temperature, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupAvailableEnergy (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupAverageEnergy (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, std::string const &a_productID, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupDepositionEnergy (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, Transporting::Particles const &a_particles, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupAvailableMomentum (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupAverageMomentum (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, std::string const &a_productID, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupDepositionMomentum (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, Transporting::Particles const &a_particles, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupGain (LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, std::string const &a_productID, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
stringAndDoublePairs muCutoffForCoulombPlusNuclearElastic () const
DelayedNeutronProducts delayedNeutronProducts () const
void incompleteParticles (Transporting::Settings const &a_settings, std::set< std::string > &a_incompleteParticles) const
void saveAs (std::string const &a_fileName) const
void toXMLList (GUPI::WriteInfo &a_writeInfo, std::string const &a_indent="") const
void parseEvaluatedTargetInfo (HAPI::Node const &a_node)
Public Member Functions inherited from GIDI::Protare
 Protare ()
 ~Protare ()
ParticleInfo const & projectile () const
void setProjectile (ParticleInfo const &a_projectile)
ParticleInfo const & target () const
void setTarget (ParticleInfo const &a_target)
ParticleInfo const & GNDS_target () const
virtual LUPI::FormatVersion const & formatVersion (std::size_t a_index=0) const =0
virtual std::string const & fileName (std::size_t a_index=0) const =0
virtual std::string const & realFileName (std::size_t a_index=0) const =0
virtual std::vector< std::string > libraries (std::size_t a_index=0) const =0
virtual std::string const & evaluation (std::size_t a_index=0) const =0
virtual Frame projectileFrame (std::size_t a_index=0) const =0
virtual void TNSL_crossSectionSumCorrection (std::string const &a_label, Functions::XYs1d &a_crossSectionSum)
virtual void TNSL_crossSectionSumCorrection (std::string const &a_label, Functions::Ys1d &a_crossSectionSum)
virtual void TNSL_crossSectionSumCorrection (std::string const &a_label, Vector &a_crossSectionSum)
ExcludeReactionsSet reactionIndicesMatchingENDLCValues (std::set< int > const &a_CValues, bool a_checkActiveState=true)
Public Member Functions inherited from GUPI::Ancestry
 Ancestry (std::string const &a_moniker, std::string const &a_attribute="")
virtual ~Ancestry ()
Ancestryoperator= (Ancestry const &a_ancestry)
std::string const & moniker () const
void setMoniker (std::string const &a_moniker)
Ancestryancestor ()
Ancestry const * ancestor () const
void setAncestor (Ancestry *a_ancestor)
std::string attribute () const
Ancestryroot ()
Ancestry const * root () const
bool isChild (Ancestry *a_instance)
bool isParent (Ancestry *a_parent)
bool isRoot () const
AncestryfindInAncestry (std::string const &a_href)
Ancestry const * findInAncestry (std::string const &a_href) const
virtual LUPI_HOST void serialize (LUPI::DataBuffer &a_buffer, LUPI::DataBuffer::Mode a_mode)
virtual std::string xlinkItemKey () const
std::string toXLink () const
void printXML () const

Additional Inherited Members

Static Public Member Functions inherited from GUPI::Ancestry
static std::string buildXLinkItemKey (std::string const &a_name, std::string const &a_key)
Protected Member Functions inherited from GIDI::Protare
void initialize (HAPI::Node const &a_node, SetupInfo &a_setupInfo, PoPI::Database const &a_pops, PoPI::Database const &a_internalPoPs, bool a_targetRequiredInGlobalPoPs, bool a_requiredInPoPs=true)

Detailed Description

Class to store a GNDS <reactionSuite> node.

Definition at line 4664 of file GIDI.hpp.

Constructor & Destructor Documentation

◆ ProtareSingle() [1/3]

GIDI::ProtareSingle::ProtareSingle ( PoPI::Database const & a_pops,
std::string const & a_projectileID,
std::string const & a_targetID,
std::string const & a_evaluation,
std::string const & a_interaction,
std::string const & a_formatVersion = GNDS_formatVersion_1_10Chars )

Parses a GNDS file to construct the Protare instance. Calls the initialize method which does most of the work.

Parameters
a_pops[in] A PoPs Database instance used to get particle indices and possibly other particle information.
a_projectileID[in] The PoPs id of the projectile.
a_targetID[in] The PoPs id of the target.
a_evaluation[in] The evaluation string.
a_interaction[in] The interaction flag for the protare.
a_formatVersion[in] The GNDS format to use.

Definition at line 155 of file GIDI_protare.cc.

156 :
157 m_doc( nullptr ),
158 m_dataManager( nullptr ),
159 m_numberOfLazyParsingHelperForms( 0 ),
160 m_numberOfLazyParsingHelperFormsReplaced( 0 ),
161 m_formatVersion( a_formatVersion ),
162 m_evaluation( a_evaluation ),
163 m_interaction( a_interaction ),
164 m_projectileFrame( Frame::lab ),
165 m_decayPositronium( false ),
166 m_thresholdFactor( 0.0 ),
167 m_nuclearPlusCoulombInterferenceOnlyReaction( nullptr ),
168 m_pointwiseAverageProductEnergy( GIDI_averageEnergyChars, GIDI_labelChars ) {
169
171 initialize( );
172
173 setProjectile( ParticleInfo( a_projectileID, a_pops, a_pops, true ) );
174 setTarget( ParticleInfo( a_targetID, a_pops, a_pops, true ) );
175}
#define GIDI_labelChars
Definition GIDI.hpp:438
#define GIDI_averageEnergyChars
Definition GIDI.hpp:225
#define GIDI_topLevelChars
Definition GIDI.hpp:169
void setTarget(ParticleInfo const &a_target)
Definition GIDI.hpp:4544
void setProjectile(ParticleInfo const &a_projectile)
Definition GIDI.hpp:4542
void setMoniker(std::string const &a_moniker)
Definition GUPI.hpp:103

Referenced by protare(), and protare().

◆ ProtareSingle() [2/3]

GIDI::ProtareSingle::ProtareSingle ( Construction::Settings const & a_construction,
std::string const & a_fileName,
FileType a_fileType,
PoPI::Database const & a_pops,
ParticleSubstitution const & a_particleSubstitution,
std::vector< std::string > const & a_libraries,
std::string const & a_interaction,
bool a_targetRequiredInGlobalPoPs = true,
bool a_requiredInPoPs = true )

Parses a GNDS file to construct the Protare instance. Calls the initialize method which does most of the work.

Parameters
a_construction[in] Used to pass user options to the constructor.
a_fileName[in] File containing a protare (i.e., reactionSuite) node that is parsed and used to construct the Protare.
a_fileType[in] File type of a_fileType. Currently, only GIDI::FileType::XML and GIDI::FileType::HDF are supported.
a_pops[in] A PoPs Database instance used to get particle indices and possibly other particle information.
a_particleSubstitution[in] Map of particles to substitute with another particles.
a_libraries[in] The list of libraries that were searched to find this.
a_interaction[in] The interaction flag for the protare.
a_targetRequiredInGlobalPoPs[in] If true, the target is required to be in a_pops.
a_requiredInPoPs[in] If true, no particle is required to be in a_pops.

Definition at line 191 of file GIDI_protare.cc.

193 :
194 Protare( ),
195 m_doc( nullptr ),
196 m_dataManager( nullptr ),
197 m_numberOfLazyParsingHelperForms( 0 ),
198 m_numberOfLazyParsingHelperFormsReplaced( 0 ),
199 m_libraries( a_libraries ),
200 m_interaction( a_interaction ),
201 m_fileName( a_fileName ),
202 m_realFileName( LUPI::FileInfo::realPath( a_fileName ) ),
203 m_decayPositronium( a_construction.decayPositronium( ) ),
204 m_pointwiseAverageProductEnergy( GIDI_averageEnergyChars, GIDI_labelChars ) {
205
206#ifdef HAPI_USE_PUGIXML
207 if( a_fileType == GIDI::FileType::XML ) {
208 m_doc = new HAPI::PugiXMLFile( a_fileName.c_str( ), "ProtareSingle::ProtareSingle" );
209 }
210#endif
211#ifdef HAPI_USE_HDF5
212 if( a_fileType == GIDI::FileType::HDF ) {
213 m_doc = new HAPI::HDFFile( a_fileName.c_str( ) );
214 }
215#endif
216 if( m_doc == nullptr ) {
217 throw std::runtime_error( "Only XML/HDF file types supported." );
218 }
219
220#ifdef GIDIP_TEST_PARSING
221 if( a_construction.parseMode( ) != Construction::ParseMode::noParsing ) {
222#endif
223 HAPI::Node protare = m_doc->first_child( );
224
225 SetupInfo setupInfo( this );
226 ParticleSubstitution particleSubstitution( a_particleSubstitution );
227 setupInfo.m_particleSubstitution = &particleSubstitution;
228
229 initialize( a_construction, protare, setupInfo, a_pops, a_targetRequiredInGlobalPoPs, a_requiredInPoPs );
230#ifdef GIDIP_TEST_PARSING
231 }
232#endif
233}
ProtareSingle * protare(std::size_t a_index)
std::map< std::string, ParticleInfo > ParticleSubstitution
Definition GIDI.hpp:487
std::string realPath(std::string const &a_path)
Definition LUPI_file.cc:63

◆ ProtareSingle() [3/3]

GIDI::ProtareSingle::ProtareSingle ( Construction::Settings const & a_construction,
HAPI::Node const & a_protare,
PoPI::Database const & a_pops,
ParticleSubstitution const & a_particleSubstitution,
std::vector< std::string > const & a_libraries,
std::string const & a_interaction,
bool a_targetRequiredInGlobalPoPs = true,
bool a_requiredInPoPs = true )

◆ ~ProtareSingle()

GIDI::ProtareSingle::~ProtareSingle ( )

Definition at line 495 of file GIDI_protare.cc.

495 {
496
497 delete m_doc;
498 delete m_dataManager;
499 delete m_nuclearPlusCoulombInterferenceOnlyReaction;
500 delete m_multiGroupSummedReaction;
501 delete m_multiGroupSummedDelayedNeutrons;
502}

Member Function Documentation

◆ ACE_URR_probabilityTables()

Suite const & GIDI::ProtareSingle::ACE_URR_probabilityTables ( ) const
inline

Returns a const reference to the m_ACE_URR_probabilityTables member.

Definition at line 4770 of file GIDI.hpp.

Referenced by MCGIDI::convertACE_URR_probabilityTablesFromGIDI().

◆ checkIf_nuclearPlusCoulombInterferenceWanted()

Reaction const * GIDI::ProtareSingle::checkIf_nuclearPlusCoulombInterferenceWanted ( Transporting::MG const & a_settings,
Reaction const * a_reaction ) const

Returns a_reaction except when Rutherford scattering present and only nuclear + Coulomb interference is wanted. Otherwise returns m_nuclearPlusCoulombInterferenceOnlyReaction which may be a nullptr.

Parameters
a_settings[in] Specifies user options.
a_reaction[in] Reaction pointer to return if check passes.
Returns
Const reaction pointer or m_nuclearPlusCoulombInterferenceOnlyReaction.

Definition at line 545 of file GIDI_protare.cc.

545 {
546
547 if( a_reaction->RutherfordScatteringPresent( ) && a_settings.nuclearPlusCoulombInterferenceOnly( ) ) {
548 return( m_nuclearPlusCoulombInterferenceOnlyReaction );
549 }
550
551 return( a_reaction );
552}

Referenced by reactionToMultiGroup().

◆ dataManager()

HAPI::DataManager * GIDI::ProtareSingle::dataManager ( )
inline

Returns the value of the m_dataManager member.

Definition at line 4734 of file GIDI.hpp.

Referenced by GIDI::parseValuesOfDoubles(), and GIDI::parseValuesOfInts().

◆ delayedNeutronProducts()

DelayedNeutronProducts GIDI::ProtareSingle::delayedNeutronProducts ( ) const
virtual

Returns the list of DelayedNeutronProduct instances.

Returns
a_delayedNeutronProducts The list of delayed neutrons.

Implements GIDI::Protare.

Definition at line 1490 of file GIDI_protare.cc.

1490 {
1491
1492 DelayedNeutronProducts delayedNeutronProducts1;
1493
1494 if( hasFission( ) ) {
1495 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1496 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
1497
1498 if( !reaction1->active( ) ) continue;
1499 if( reaction1->hasFission( ) ) reaction1->delayedNeutronProducts( delayedNeutronProducts1 );
1500 }
1501 }
1502
1503 return( delayedNeutronProducts1 );
1504}
bool hasFission() const
std::vector< DelayedNeutronProduct > DelayedNeutronProducts
Definition GIDI.hpp:4031

◆ documentations()

Documentation_1_10::Suite & GIDI::ProtareSingle::documentations ( )
inlinevirtual

Returns the value of the m_documentations member.

Implements GIDI::Protare.

Definition at line 4795 of file GIDI.hpp.

◆ evaluation()

std::string const & GIDI::ProtareSingle::evaluation ( LUPI_maybeUnused std::size_t a_index = 0) const
inline

Returns the value of the m_evaluation member.

Definition at line 4786 of file GIDI.hpp.

Referenced by GIDI::Transporting::MG::form(), and toXMLList().

◆ externalFile()

ExternalFile const & GIDI::ProtareSingle::externalFile ( std::string const a_label) const
inline

Returns the external file with label a_label.

Definition at line 4797 of file GIDI.hpp.

◆ externalFiles()

ExternalFiles::Suite const & GIDI::ProtareSingle::externalFiles ( ) const
inline

Returns the value of the m_externalFiles member.

Definition at line 4798 of file GIDI.hpp.

◆ fileName()

std::string const & GIDI::ProtareSingle::fileName ( LUPI_maybeUnused std::size_t a_index = 0) const
inline

Returns the value of the m_fileName member.

Definition at line 4782 of file GIDI.hpp.

◆ findInAncestry3() [1/2]

GUPI::Ancestry * GIDI::ProtareSingle::findInAncestry3 ( std::string const & a_item)
virtual

Used by GUPI::Ancestry to tranverse GNDS nodes. This method returns a pointer to a derived class' a_item member or nullptr if none exists.

Parameters
a_item[in] The name of the class member whose pointer is to be return.
Returns
The pointer to the class member or nullptr if class does not have a member named a_item.

Implements GIDI::Protare.

Definition at line 897 of file GIDI_protare.cc.

897 {
898
899 if( a_item == GIDI_stylesChars ) return( &m_styles );
900 if( a_item == GIDI_reactionsChars ) return( &m_reactions );
901 if( a_item == GIDI_orphanProductsChars ) return( &m_orphanProducts );
902 if( a_item == GIDI_sumsChars ) return( &m_sums );
903 if( a_item == GIDI_fissionComponentsChars ) return( &m_fissionComponents );
904 if( a_item == GIDI_ACE_URR_probabilityTablesChars ) return( &m_ACE_URR_probabilityTables );
905
906 return( nullptr );
907}
#define GIDI_stylesChars
Definition GIDI.hpp:177
#define GIDI_reactionsChars
Definition GIDI.hpp:179
#define GIDI_sumsChars
Definition GIDI.hpp:204
#define GIDI_orphanProductsChars
Definition GIDI.hpp:181
#define GIDI_fissionComponentsChars
Definition GIDI.hpp:184
#define GIDI_ACE_URR_probabilityTablesChars
Definition GIDI.hpp:186

◆ findInAncestry3() [2/2]

GUPI::Ancestry const * GIDI::ProtareSingle::findInAncestry3 ( std::string const & a_item) const
virtual

Used by GUPI::Ancestry to tranverse GNDS nodes. This method returns a pointer to a derived class' a_item member or nullptr if none exists.

Parameters
a_item[in] The name of the class member whose pointer is to be return.
Returns
The pointer to the class member or nullptr if class does not have a member named a_item.

Implements GIDI::Protare.

Definition at line 916 of file GIDI_protare.cc.

916 {
917
918 if( a_item == GIDI_stylesChars ) return( &m_styles );
919 if( a_item == GIDI_reactionsChars ) return( &m_reactions );
920 if( a_item == GIDI_orphanProductsChars ) return( &m_orphanProducts );
921 if( a_item == GIDI_sumsChars ) return( &m_sums );
922 if( a_item == GIDI_fissionComponentsChars ) return( &m_fissionComponents );
923 if( a_item == GIDI_ACE_URR_probabilityTablesChars ) return( &m_ACE_URR_probabilityTables );
924
925 return( nullptr );
926}

◆ fissionComponents()

Suite & GIDI::ProtareSingle::fissionComponents ( )
inline

Returns a reference to the m_fissionComponents member.

Definition at line 4757 of file GIDI.hpp.

◆ formatVersion()

LUPI::FormatVersion const & GIDI::ProtareSingle::formatVersion ( LUPI_maybeUnused std::size_t a_index = 0) const
inline

Returns the value of the m_formatVersion member.

Definition at line 4781 of file GIDI.hpp.

◆ GRIN_continuumGammas2()

GRIN::GRIN_continuumGammas const * GIDI::ProtareSingle::GRIN_continuumGammas2 ( ) const
inline

Returns a const pointer to the m_GRIN_continuumGammas member.

Definition at line 4772 of file GIDI.hpp.

◆ groupBoundaries()

std::vector< double > GIDI::ProtareSingle::groupBoundaries ( Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
std::string const & a_productID ) const
virtual

Returns the multi-group boundaries for the requested label and product.

Parameters
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_productID[in] ID for the requested product.
Returns
List of multi-group boundaries.

Implements GIDI::Protare.

Definition at line 827 of file GIDI_protare.cc.

827 {
828
829 Styles::HeatedMultiGroup const *heatedMultiGroupStyle1 = m_styles.get<Styles::HeatedMultiGroup>( a_temperatureInfo.heatedMultiGroup( ) );
830
831 return( heatedMultiGroupStyle1->groupBoundaries( a_productID ) );
832}

◆ hasFission()

bool GIDI::ProtareSingle::hasFission ( ) const
virtual

Returns true if at least one reaction contains a fission channel.

Returns
true if at least one reaction contains a fission channel.

Implements GIDI::Protare.

Definition at line 858 of file GIDI_protare.cc.

858 {
859
860 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
861 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
862
863 if( !reaction1->active( ) ) continue;
864 if( reaction1->hasFission( ) ) return( true );
865 }
866 return( false );
867}

Referenced by delayedNeutronProducts().

◆ incompleteParticles()

void GIDI::ProtareSingle::incompleteParticles ( Transporting::Settings const & a_settings,
std::set< std::string > & a_incompleteParticles ) const
virtual

Calls the incompleteParticles method for each active reaction in the reactions and orphanProducts nodes.

Parameters
a_settings[in] Specifies the requested label.
a_incompleteParticles[out] The list of particles whose completeParticle method returns false.

Implements GIDI::Protare.

Definition at line 1513 of file GIDI_protare.cc.

1513 {
1514
1515 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1516 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
1517
1518 if( !reaction1->active( ) ) continue;
1519 reaction1->incompleteParticles( a_settings, a_incompleteParticles );
1520 }
1521
1522 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
1523 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
1524
1525 if( !reaction1->active( ) ) continue;
1526 reaction1->incompleteParticles( a_settings, a_incompleteParticles );
1527 }
1528}

◆ incompleteReaction() [1/2]

Reaction * GIDI::ProtareSingle::incompleteReaction ( std::size_t a_index)
inline

Returns the a_index - 1 reaction.

Definition at line 4827 of file GIDI.hpp.

◆ incompleteReaction() [2/2]

Reaction const * GIDI::ProtareSingle::incompleteReaction ( std::size_t a_index) const
inline

Returns the a_index - 1 reaction.

Definition at line 4828 of file GIDI.hpp.

◆ incompleteReactions() [1/2]

Suite & GIDI::ProtareSingle::incompleteReactions ( )
inline

Returns a reference to the m_incompleteReactions member.

Definition at line 4752 of file GIDI.hpp.

◆ incompleteReactions() [2/2]

Suite const & GIDI::ProtareSingle::incompleteReactions ( ) const
inline

Returns a const reference to the m_incompleteReactions member.

Definition at line 4753 of file GIDI.hpp.

◆ incrementNumberOfLazyParsingHelperForms()

void GIDI::ProtareSingle::incrementNumberOfLazyParsingHelperForms ( )
inline

Definition at line 4737 of file GIDI.hpp.

4737{ ++m_numberOfLazyParsingHelperForms; }

◆ incrementNumberOfLazyParsingHelperFormsReplaced()

void GIDI::ProtareSingle::incrementNumberOfLazyParsingHelperFormsReplaced ( )
inline

Increments the m_numberOfLazyParsingHelperForms member of this by 1.

Definition at line 4739 of file GIDI.hpp.

4739{ ++m_numberOfLazyParsingHelperFormsReplaced; }

◆ interaction()

std::string const & GIDI::ProtareSingle::interaction ( LUPI_maybeUnused std::size_t a_index = 0) const
inline

Returns the value of the m_interaction member.

Definition at line 4787 of file GIDI.hpp.

◆ internalPoPs() [1/2]

PoPI::Database & GIDI::ProtareSingle::internalPoPs ( )
inline

Returns a reference to the m_internalPoPs member.

Definition at line 4806 of file GIDI.hpp.

◆ internalPoPs() [2/2]

PoPI::Database const & GIDI::ProtareSingle::internalPoPs ( ) const
inline

Returns a const reference to the m_internalPoPs member.

Definition at line 4805 of file GIDI.hpp.

◆ intid()

int GIDI::ProtareSingle::intid ( std::string const & a_id) const
virtual

Returns the intid for the requested particle or -1 if the particle is not in this PoPs database.

Parameters
a_id[in] The GNDS PoPs id for particle whose intd is requested.
Returns
C++ int for the requested particle or -1 if particle is not in PoPs.

Implements GIDI::Protare.

Definition at line 611 of file GIDI_protare.cc.

611 {
612
613 return( m_internalPoPs.intid( a_id ) );
614}

◆ isDelayedFissionNeutronComplete()

bool GIDI::ProtareSingle::isDelayedFissionNeutronComplete ( ) const
virtual

Returns false* if protare has delayed fission neutrons for an active reaction and they are not complete; otherwise, returns **true.

Returns
bool

Implements GIDI::Protare.

Definition at line 875 of file GIDI_protare.cc.

875 {
876
877 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
878 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
879
880 if( !reaction1->active( ) ) continue;
881 if( reaction1->hasFission( ) ) {
882 OutputChannel const *outputChannel = reaction1->outputChannel( );
883 if( !outputChannel->isDelayedFissionNeutronComplete( ) ) return( false );
884 }
885 }
886
887 return( true );
888}

◆ isPhotoAtomic()

bool GIDI::ProtareSingle::isPhotoAtomic ( ) const
inline

Returns the value of the m_isPhotoAtomic member.

Definition at line 4746 of file GIDI.hpp.

Referenced by GIDI::Product::Product(), and MCGIDI::ProtareSingle::ProtareSingle().

◆ isTNSL_ProtareSingle()

bool GIDI::ProtareSingle::isTNSL_ProtareSingle ( ) const
inlinevirtual

Returns true if the instance is a ProtareSingle instance with only TNSL data and false otherwise.

Reimplemented from GIDI::Protare.

Definition at line 4745 of file GIDI.hpp.

Referenced by GIDI::Product::Product().

◆ libraries()

std::vector< std::string > GIDI::ProtareSingle::libraries ( LUPI_maybeUnused std::size_t a_index = 0) const
inline

Returns the libraries that this resided in.

Definition at line 4785 of file GIDI.hpp.

◆ maximumLegendreOrder()

int GIDI::ProtareSingle::maximumLegendreOrder ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
std::string const & a_productID ) const
virtual

Determines the maximum Legredre order present in the multi-group transfer matrix for a give product for a give label.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_productID[in] The id of the requested product.
Returns
The maximum Legredre order. If no transfer matrix data are present for the requested product, -1 is returned.

Implements GIDI::Protare.

Definition at line 653 of file GIDI_protare.cc.

654 {
655
656 int _maximumLegendreOrder = -1;
657 ExcludeReactionsSet excludeReactionsSet;
658
659 if( useMultiGroupSummedData( a_settings, excludeReactionsSet ) ) {
660 _maximumLegendreOrder = m_multiGroupSummedReaction->maximumLegendreOrder( a_smr, a_settings, a_temperatureInfo, a_productID );
661 if( useMultiGroupSummedDelayedNeutronsData( a_settings ) ) {
662 _maximumLegendreOrder = std::max( _maximumLegendreOrder, m_multiGroupSummedDelayedNeutrons->maximumLegendreOrder(
663 a_smr, a_settings, a_temperatureInfo, a_productID ) );
664 } }
665 else {
666 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
667 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
668
669 if( !reaction1->active( ) ) continue;
670 int r_maximumLegendreOrder = reaction1->maximumLegendreOrder( a_smr, a_settings, a_temperatureInfo, a_productID );
671 if( r_maximumLegendreOrder > _maximumLegendreOrder ) _maximumLegendreOrder = r_maximumLegendreOrder;
672 }
673
674 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
675 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
676
677 if( !reaction1->active( ) ) continue;
678 int r_maximumLegendreOrder = reaction1->maximumLegendreOrder( a_smr, a_settings, a_temperatureInfo, a_productID );
679 if( r_maximumLegendreOrder > _maximumLegendreOrder ) _maximumLegendreOrder = r_maximumLegendreOrder;
680 }
681 }
682
683 return( _maximumLegendreOrder );
684}
std::set< std::size_t > ExcludeReactionsSet
Definition GIDI.hpp:47

◆ muCutoffForCoulombPlusNuclearElastic()

stringAndDoublePairs GIDI::ProtareSingle::muCutoffForCoulombPlusNuclearElastic ( ) const
virtual
Returns
A list of label, mu cutoff pairs.

Implements GIDI::Protare.

Definition at line 1465 of file GIDI_protare.cc.

1465 {
1466
1467 stringAndDoublePairs muCutoffs;
1468
1469 for( std::size_t i1 = 0; i1 < m_styles.size( ); ++i1 ) {
1470 Styles::Base const *style1 = m_styles.get<Styles::Base>( i1 );
1471
1472 if( style1->moniker( ) == GIDI_CoulombPlusNuclearElasticMuCutoffStyleChars ) {
1473 Styles::CoulombPlusNuclearElasticMuCutoff const *style2 = static_cast<Styles::CoulombPlusNuclearElasticMuCutoff const *>( style1 );
1474
1475 stringAndDoublePair labelMu( style2->label( ), style2->muCutoff( ) );
1476
1477 muCutoffs.push_back( std::move( labelMu ) );
1478 }
1479 }
1480
1481 return( muCutoffs );
1482}
#define GIDI_CoulombPlusNuclearElasticMuCutoffStyleChars
Definition GIDI.hpp:243
std::pair< std::string, double > stringAndDoublePair
Definition GIDI.hpp:485
std::vector< stringAndDoublePair > stringAndDoublePairs
Definition GIDI.hpp:486

◆ multiGroupAvailableEnergy()

Vector GIDI::ProtareSingle::multiGroupAvailableEnergy ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total available energy for the requested label. This is a cross section weighted available energy summed over all reactions.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the available energy.
Returns
The requested multi-group available energy as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1208 of file GIDI_protare.cc.

1209 {
1210
1211 Vector vector( 0 );
1212
1213 if( useMultiGroupSummedData( a_settings, a_reactionsToExclude ) ) {
1214 vector = m_multiGroupSummedReaction->multiGroupAvailableEnergy( a_smr, a_settings, a_temperatureInfo ); }
1215 else {
1216 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1217 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
1218
1219 if( reaction1 != nullptr ) vector += reaction1->multiGroupAvailableEnergy( a_smr, a_settings, a_temperatureInfo );
1220 }
1221 }
1222
1223 return( vector );
1224}
Reaction const * reactionToMultiGroup(Transporting::MG const &a_settings, std::size_t a_index, ExcludeReactionsSet const &a_reactionsToExclude) const

Referenced by multiGroupDepositionEnergy().

◆ multiGroupAvailableMomentum()

Vector GIDI::ProtareSingle::multiGroupAvailableMomentum ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total available momentum for the requested label. This is a cross section weighted available momentum summed over all reactions.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the available momentum.
Returns
The requested multi-group available momentum as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1343 of file GIDI_protare.cc.

1344 {
1345
1346 Vector vector( 0 );
1347
1348 if( useMultiGroupSummedData( a_settings, a_reactionsToExclude ) ) {
1349 vector = m_multiGroupSummedReaction->multiGroupAvailableMomentum( a_smr, a_settings, a_temperatureInfo ); }
1350 else {
1351 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1352 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
1353
1354 if( reaction1 != nullptr ) vector += reaction1->multiGroupAvailableMomentum( a_smr, a_settings, a_temperatureInfo );
1355 }
1356 }
1357
1358 return( vector );
1359}

Referenced by multiGroupDepositionMomentum().

◆ multiGroupAverageEnergy()

Vector GIDI::ProtareSingle::multiGroupAverageEnergy ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
std::string const & a_productID,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total average energy for the requested label for the requested product. This is a cross section weighted average energy summed over all reactions.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_productID[in] Particle id for the requested product.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the average energy.
Returns
The requested multi-group average energy as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1239 of file GIDI_protare.cc.

1241 {
1242
1243 Vector vector( 0 );
1244
1245 if( useMultiGroupSummedData( a_settings, a_reactionsToExclude ) ) {
1246 vector = m_multiGroupSummedReaction->multiGroupAverageEnergy( a_smr, a_settings, a_temperatureInfo, a_productID );
1247 if( useMultiGroupSummedDelayedNeutronsData( a_settings ) ) {
1248 vector += m_multiGroupSummedDelayedNeutrons->multiGroupAverageEnergy( a_smr, a_settings, a_temperatureInfo, a_productID );
1249 } }
1250 else {
1251 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1252 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
1253
1254 if( reaction1 != nullptr ) vector += reaction1->multiGroupAverageEnergy( a_smr, a_settings, a_temperatureInfo, a_productID );
1255 }
1256
1257 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
1258 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
1259
1260 if( !reaction1->active( ) ) continue;
1261 vector += reaction1->multiGroupAverageEnergy( a_smr, a_settings, a_temperatureInfo, a_productID );
1262 }
1263 }
1264
1265 return( vector );
1266}

Referenced by multiGroupDepositionEnergy().

◆ multiGroupAverageMomentum()

Vector GIDI::ProtareSingle::multiGroupAverageMomentum ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
std::string const & a_productID,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total average momentum for the requested label for the requested product. This is a cross section weighted average momentum summed over all reactions.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_productID[in] Particle id for the requested product.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the average momentum.
Returns
The requested multi-group average momentum as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1374 of file GIDI_protare.cc.

1376 {
1377
1378 Vector vector( 0 );
1379
1380 if( useMultiGroupSummedData( a_settings, a_reactionsToExclude ) ) {
1381 vector = m_multiGroupSummedReaction->multiGroupAverageMomentum( a_smr, a_settings, a_temperatureInfo, a_productID );
1382 if( useMultiGroupSummedDelayedNeutronsData( a_settings ) ) {
1383 vector += m_multiGroupSummedDelayedNeutrons->multiGroupAverageMomentum( a_smr, a_settings, a_temperatureInfo, a_productID );
1384 } }
1385 else {
1386 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1387 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
1388
1389 if( reaction1 != nullptr ) vector += reaction1->multiGroupAverageMomentum( a_smr, a_settings, a_temperatureInfo, a_productID );
1390 }
1391
1392 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
1393 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
1394
1395 if( !reaction1->active( ) ) continue;
1396 vector += reaction1->multiGroupAverageMomentum( a_smr, a_settings, a_temperatureInfo, a_productID );
1397 }
1398 }
1399
1400 return( vector );
1401}

Referenced by multiGroupDepositionMomentum().

◆ multiGroupCrossSection()

Vector GIDI::ProtareSingle::multiGroupCrossSection ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {},
std::string const & a_label = "" ) const
virtual

Returns the multi-group, total cross section for the requested label. This is summed over all reactions.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the cross section.
a_label[in] If not an empty string, this is used as the label for the form to return and the a_temperatureInfo labels are ignored.
Returns
The requested multi-group cross section as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 940 of file GIDI_protare.cc.

941 {
942
943 Vector vector;
944
945 if( useMultiGroupSummedData( a_settings, a_reactionsToExclude ) ) {
946 vector = m_multiGroupSummedReaction->multiGroupCrossSection( a_smr, a_settings, a_temperatureInfo, a_label ); }
947 else {
948 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
949 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
950
951 if( reaction1 != nullptr ) vector += reaction1->multiGroupCrossSection( a_smr, a_settings, a_temperatureInfo, a_label );
952 }
953 }
954
955 return( vector );
956}

◆ multiGroupDepositionEnergy()

Vector GIDI::ProtareSingle::multiGroupDepositionEnergy ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
Transporting::Particles const & a_particles,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total deposition energy for the requested label. This is a cross section weighted deposition energy summed over all reactions. The deposition energy is calculated by subtracting the average energy from each transportable particle from the available energy. The list of transportable particles is specified via the list of particle specified in the a_settings argument.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label and the products that are transported.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_particles[in] The list of particles to be transported.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the deposition energy.
Returns
The requested multi-group deposition energy as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1282 of file GIDI_protare.cc.

1284 {
1285
1286 bool atLeastOneReactionHasAllParticlesTracked = false;
1287 Vector vector;
1288
1289 if( a_settings.zeroDepositionIfAllProductsTracked( ) ) {
1290 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1291 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
1292
1293 if( reaction1 != nullptr ) {
1294 if( reaction1->areAllProductsTracked( a_particles ) ) {
1295 atLeastOneReactionHasAllParticlesTracked = true;
1296 break;
1297 }
1298 }
1299 }
1300 }
1301
1302 if( atLeastOneReactionHasAllParticlesTracked ) {
1303 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1304 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
1305
1306 if( reaction1 != nullptr ) vector += reaction1->multiGroupDepositionEnergy( a_smr, a_settings, a_temperatureInfo, a_particles );
1307 }
1308 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
1309 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
1310
1311 if( !reaction1->active( ) ) continue;
1312 vector += reaction1->multiGroupDepositionEnergy( a_smr, a_settings, a_temperatureInfo, a_particles );
1313 } }
1314 else {
1315 std::map<std::string, Transporting::Particle> const &products( a_particles.particles( ) );
1316 vector = multiGroupAvailableEnergy( a_smr, a_settings, a_temperatureInfo, a_reactionsToExclude );
1317 Vector availableEnergy( vector );
1318
1319 for( std::map<std::string, Transporting::Particle>::const_iterator iter = products.begin( ); iter != products.end( ); ++iter ) {
1320 vector -= multiGroupAverageEnergy( a_smr, a_settings, a_temperatureInfo, iter->first, a_reactionsToExclude );
1321 }
1322
1323 for( std::size_t index = 0; index < availableEnergy.size( ); ++index ) { // Check for values that should probably be 0.0.
1324 if( std::fabs( vector[index] ) < std::fabs( 1e-14 * availableEnergy[index] ) ) vector[index] = 0.0;
1325 }
1326 }
1327
1328 return( vector );
1329}
Vector multiGroupAvailableEnergy(LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupAverageEnergy(LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, std::string const &a_productID, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const

◆ multiGroupDepositionMomentum()

Vector GIDI::ProtareSingle::multiGroupDepositionMomentum ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
Transporting::Particles const & a_particles,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total deposition momentum for the requested label. This is a cross section weighted deposition momentum summed over all reactions. The deposition momentum is calculated by subtracting the average momentum from each transportable particle from the available momentum. The list of transportable particles is specified via the list of particle specified in the a_settings argument.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_particles[in] The list of particles to be transported.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the deposition momentum.
Returns
The requested multi-group deposition momentum as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1417 of file GIDI_protare.cc.

1419 {
1420
1421 std::map<std::string, Transporting::Particle> const &products( a_particles.particles( ) );
1422 Vector vector = multiGroupAvailableMomentum( a_smr, a_settings, a_temperatureInfo, a_reactionsToExclude );
1423
1424 for( std::map<std::string, Transporting::Particle>::const_iterator iter = products.begin( ); iter != products.end( ); ++iter ) {
1425 vector -= multiGroupAverageMomentum( a_smr, a_settings, a_temperatureInfo, iter->first, a_reactionsToExclude );
1426 }
1427
1428 return( vector );
1429}
Vector multiGroupAvailableMomentum(LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector multiGroupAverageMomentum(LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, std::string const &a_productID, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const

◆ multiGroupFissionGammaMultiplicity()

Vector GIDI::ProtareSingle::multiGroupFissionGammaMultiplicity ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total fission gamma multiplicity for the requested label. This is a cross section weighted multiplicity.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the multiplicity.
Returns
The requested multi-group fission neutron multiplicity as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1036 of file GIDI_protare.cc.

1037 {
1038
1039 Vector vector( 0 );
1040
1041 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1042 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
1043
1044 if( !reaction1->active( ) ) continue;
1045 if( reaction1->hasFission( ) ) vector += reaction1->multiGroupMultiplicity( a_smr, a_settings, a_temperatureInfo, PoPI::IDs::photon );
1046 }
1047
1048 return( vector );
1049}
static std::string const photon
Definition PoPI.hpp:162

◆ multiGroupFissionMatrix()

Matrix GIDI::ProtareSingle::multiGroupFissionMatrix ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
Transporting::Particles const & a_particles,
std::size_t a_order,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Like ProtareSingle::multiGroupProductMatrix, but only returns the fission neutron, transfer matrix.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label and if delayed neutrons should be included.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_particles[in] The list of particles to be transported.
a_order[in] Requested product matrix, Legendre order.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the fission matrix.
Returns
The requested multi-group neutron fission matrix as a GIDI::Matrix.

Implements GIDI::Protare.

Definition at line 1143 of file GIDI_protare.cc.

1145 {
1146
1147 Matrix matrix( 0, 0 );
1148
1149 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1150 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
1151
1152 if( !reaction1->active( ) ) continue;
1153 if( reaction1->hasFission( ) ) matrix += reaction1->multiGroupFissionMatrix( a_smr, a_settings, a_temperatureInfo, a_particles, a_order );
1154 }
1155
1156 return( matrix );
1157}

◆ multiGroupFissionNeutronMultiplicity()

Vector GIDI::ProtareSingle::multiGroupFissionNeutronMultiplicity ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total fission neutron multiplicity for the requested label. This is a cross section weighted multiplicity.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the cross multiplicity.
Returns
The requested multi-group fission neutron multiplicity as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1010 of file GIDI_protare.cc.

1011 {
1012
1013 Vector vector( 0 );
1014
1015 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1016 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
1017
1018 if( !reaction1->active( ) ) continue;
1019 if( reaction1->hasFission( ) ) vector += reaction1->multiGroupMultiplicity( a_smr, a_settings, a_temperatureInfo, PoPI::IDs::neutron );
1020 }
1021
1022 return( vector );
1023}
static std::string const neutron
Definition PoPI.hpp:164

◆ multiGroupGain()

Vector GIDI::ProtareSingle::multiGroupGain ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
std::string const & a_productID,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, gain for the requested particle and label. This is a cross section weighted gain summed over all reactions.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_productID[in] The particle PoPs' id for the whose gain is to be calculated.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the gain.
Returns
The requested multi-group gain as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1443 of file GIDI_protare.cc.

1445 {
1446
1447 Vector vector( 0 );
1448 std::string const projectile_ID = projectile( ).ID( );
1449
1450 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1451 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
1452
1453 if( reaction1 != nullptr ) vector += reaction1->multiGroupGain( a_smr, a_settings, a_temperatureInfo, a_productID, projectile_ID );
1454 }
1455
1456 return( vector );
1457}
std::string const & ID() const
Definition GIDI.hpp:760
ParticleInfo const & projectile() const
Definition GIDI.hpp:4541

◆ multiGroupInverseSpeed()

Vector GIDI::ProtareSingle::multiGroupInverseSpeed ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo ) const
virtual

Returns the inverse speeds for the requested label. The label must be for a heated multi-group style.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
Returns
List of inverse speeds.

Implements GIDI::Protare.

Definition at line 844 of file GIDI_protare.cc.

845 {
846
847 Styles::HeatedMultiGroup const *heatedMultiGroupStyle1 = m_styles.get<Styles::HeatedMultiGroup>( a_temperatureInfo.heatedMultiGroup( ) );
848
849 return( heatedMultiGroupStyle1->inverseSpeedData( ) );
850}

◆ multiGroupMultiplicity()

Vector GIDI::ProtareSingle::multiGroupMultiplicity ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
std::string const & a_productID,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total multiplicity for the requested label for the requested product. This is a cross section weighted multiplicity.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_productID[in] Id for the requested product.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the multiplicity.
Returns
The requested multi-group multiplicity as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 970 of file GIDI_protare.cc.

972 {
973
974 Vector vector;
975
976 if( useMultiGroupSummedData( a_settings, a_reactionsToExclude ) ) {
977 vector = m_multiGroupSummedReaction->multiGroupMultiplicity( a_smr, a_settings, a_temperatureInfo, a_productID );
978 if( useMultiGroupSummedDelayedNeutronsData( a_settings ) ) {
979 vector += m_multiGroupSummedDelayedNeutrons->multiGroupMultiplicity( a_smr, a_settings, a_temperatureInfo, a_productID );
980 } }
981 else {
982 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
983 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
984
985 if( reaction1 != nullptr ) vector += reaction1->multiGroupMultiplicity( a_smr, a_settings, a_temperatureInfo, a_productID );
986 }
987
988 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
989 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
990
991 if( !reaction1->active( ) ) continue;
992 vector += reaction1->multiGroupMultiplicity( a_smr, a_settings, a_temperatureInfo, a_productID );
993 }
994 }
995
996 return( vector );
997}

◆ multiGroupProductMatrix()

Matrix GIDI::ProtareSingle::multiGroupProductMatrix ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
Transporting::Particles const & a_particles,
std::string const & a_productID,
std::size_t a_order,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total product matrix for the requested label for the requested product id for the requested Legendre order. If no data are found, an empty GIDI::Matrix is returned.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label and if delayed neutrons should be included.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_particles[in] The list of particles to be transported.
a_productID[in] PoPs id for the requested product.
a_order[in] Requested product matrix, Legendre order.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the product matrix.
Returns
The requested multi-group product matrix as a GIDI::Matrix.

Implements GIDI::Protare.

Definition at line 1101 of file GIDI_protare.cc.

1103 {
1104
1105 Matrix matrix( 0, 0 );
1106
1107 if( useMultiGroupSummedData( a_settings, a_reactionsToExclude ) ) {
1108 matrix = m_multiGroupSummedReaction->multiGroupProductMatrix( a_smr, a_settings, a_temperatureInfo, a_particles, a_productID, a_order );
1109 if( useMultiGroupSummedDelayedNeutronsData( a_settings ) ) {
1110 matrix += m_multiGroupSummedDelayedNeutrons->multiGroupProductMatrix( a_smr, a_settings, a_temperatureInfo, a_particles, a_productID, a_order );
1111 } }
1112 else {
1113 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1114 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
1115
1116 if( reaction1 != nullptr ) matrix += reaction1->multiGroupProductMatrix( a_smr, a_settings, a_temperatureInfo, a_particles, a_productID, a_order );
1117 }
1118
1119 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
1120 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
1121
1122 if( !reaction1->active( ) ) continue;
1123 matrix += reaction1->multiGroupProductMatrix( a_smr, a_settings, a_temperatureInfo, a_particles, a_productID, a_order );
1124 }
1125 }
1126
1127 return( matrix );
1128}

Referenced by multiGroupTransportCorrection().

◆ multiGroupQ()

Vector GIDI::ProtareSingle::multiGroupQ ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
bool a_final,
bool a_effectivePhotoAtomic = true,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group, total Q for the requested label. This is a cross section weighted Q summed over all reactions

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_final[in] If false, only the Q for the primary reactions are return, otherwise, the Q for the final reactions.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the Q.
Returns
The requested multi-group Q as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1064 of file GIDI_protare.cc.

1066 {
1067
1068 Vector vector( 0 );
1069
1070 if( useMultiGroupSummedData( a_settings, a_reactionsToExclude ) ) {
1071 vector = m_multiGroupSummedReaction->multiGroupQ( a_smr, a_settings, a_temperatureInfo, a_final );
1072 if( useMultiGroupSummedDelayedNeutronsData( a_settings ) ) {
1073 vector += m_multiGroupSummedDelayedNeutrons->multiGroupQ( a_smr, a_settings, a_temperatureInfo, a_final );
1074 } }
1075 else {
1076 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
1077 Reaction const *reaction1 = reactionToMultiGroup( a_settings, i1, a_reactionsToExclude );
1078
1079 if( reaction1 != nullptr ) vector += reaction1->multiGroupQ( a_smr, a_settings, a_temperatureInfo, a_final );
1080 }
1081 }
1082
1083 return( vector );
1084}

◆ multiGroupSummedDelayedNeutrons()

OutputChannel const * GIDI::ProtareSingle::multiGroupSummedDelayedNeutrons ( ) const
inline

Returns the m_multiGroupSummedReaction member which is a pointer.

Definition at line 4769 of file GIDI.hpp.

◆ multiGroupSummedReaction()

Reaction const * GIDI::ProtareSingle::multiGroupSummedReaction ( ) const
inline

Returns the m_multiGroupSummedReaction member which is a pointer.

Definition at line 4768 of file GIDI.hpp.

◆ multiGroupTransportCorrection()

Vector GIDI::ProtareSingle::multiGroupTransportCorrection ( LUPI::StatusMessageReporting & a_smr,
Transporting::MG const & a_settings,
Styles::TemperatureInfo const & a_temperatureInfo,
Transporting::Particles const & a_particles,
std::size_t a_order,
TransportCorrectionType a_transportCorrectionType,
double a_temperature,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the multi-group transport correction for the requested label. The transport correction is calculated from the transfer matrix for the projectile id for the Legendre order of a_order + 1.

Parameters
a_smr[Out] If errors are not to be thrown, then the error is reported via this instance.
a_settings[in] Specifies the requested label.
a_temperatureInfo[in] Specifies the temperature and labels use to lookup the requested data.
a_particles[in] The list of particles to be transported.
a_order[in] Maximum Legendre order for transport. The returned transport correction is for the next higher Legender order.
a_transportCorrectionType[in] Requested transport correction type.
a_temperature[in] The temperature of the flux to use when collapsing. Pass to the GIDI::collapse method.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the transport correction.
Returns
The requested multi-group transport correction as a GIDI::Vector.

Implements GIDI::Protare.

Definition at line 1175 of file GIDI_protare.cc.

1177 {
1178
1179 if( a_transportCorrectionType == TransportCorrectionType::None ) return( Vector( 0 ) );
1180
1181 Matrix matrix( multiGroupProductMatrix( a_smr, a_settings, a_temperatureInfo, a_particles, projectile( ).ID( ), a_order + 1 ) );
1182 Matrix matrixCollapsed = collapse( matrix, a_settings, a_particles, a_temperature, projectile( ).ID( ) );
1183 std::size_t size = matrixCollapsed.size( );
1184 std::vector<double> transportCorrection1( size, 0 );
1185
1186 if( a_transportCorrectionType == TransportCorrectionType::None ) {
1187 }
1188 else if( a_transportCorrectionType == TransportCorrectionType::Pendlebury ) {
1189 for( std::size_t index = 0; index < size; ++index ) transportCorrection1[index] = matrixCollapsed[index][index]; }
1190 else {
1191 throw Exception( "Unsupported transport correction: only None and Pendlebury (i.e., Pendlebury/Underhill) are currently supported." );
1192 }
1193 return( Vector( transportCorrection1 ) );
1194}
Matrix multiGroupProductMatrix(LUPI::StatusMessageReporting &a_smr, Transporting::MG const &a_settings, Styles::TemperatureInfo const &a_temperatureInfo, Transporting::Particles const &a_particles, std::string const &a_productID, std::size_t a_order, ExcludeReactionsSet const &a_reactionsToExclude=ExcludeReactionsSet {}) const
Vector collapse(Vector const &a_vector, Transporting::Settings const &a_settings, Transporting::Particles const &a_particles, double a_temperature)

◆ nuclearPlusCoulombInterferenceOnlyReaction()

Reaction const * GIDI::ProtareSingle::nuclearPlusCoulombInterferenceOnlyReaction ( ) const
inline

Returns the m_nuclearPlusCoulombInterferenceOnlyReaction member which is a pointer.

Definition at line 4763 of file GIDI.hpp.

◆ nuclideGammaBranchStateInfos()

PoPI::NuclideGammaBranchStateInfos const & GIDI::ProtareSingle::nuclideGammaBranchStateInfos ( ) const
inline

Returns the value of the m_nuclideGammaBranchStateInfos member.

Definition at line 4731 of file GIDI.hpp.

Referenced by GIDI::Functions::Branching1d::Branching1d().

◆ numberOfInactiveReactions()

std::size_t GIDI::ProtareSingle::numberOfInactiveReactions ( ) const

The method returns the number of reactions of this that have been deactivated.

Returns
The number of deactivated reaction of this.

Definition at line 780 of file GIDI_protare.cc.

780 {
781
782 std::size_t count = 0;
783
784 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
785 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
786
787 if( !reaction1->active( ) ) ++count;
788 }
789
790 return( count );
791}

◆ numberOfIncompleteReactions()

std::size_t GIDI::ProtareSingle::numberOfIncompleteReactions ( ) const
inline

Returns the number of incomplete reactions in the Protare.

Definition at line 4826 of file GIDI.hpp.

◆ numberOfLazyParsingHelperForms()

int GIDI::ProtareSingle::numberOfLazyParsingHelperForms ( ) const
inlinevirtual

Returns the value of the m_numberOfLazyParsingHelperForms member.

Implements GIDI::Protare.

Definition at line 4789 of file GIDI.hpp.

◆ numberOfLazyParsingHelperFormsReplaced()

int GIDI::ProtareSingle::numberOfLazyParsingHelperFormsReplaced ( ) const
inlinevirtual

Returns the value of the m_numberOfLazyParsingHelperFormsReplaced member.

Implements GIDI::Protare.

Definition at line 4791 of file GIDI.hpp.

◆ numberOfOrphanProducts()

std::size_t GIDI::ProtareSingle::numberOfOrphanProducts ( ) const
inlinevirtual

Returns the number of orphan product reactions in the Protare.

Implements GIDI::Protare.

Definition at line 4822 of file GIDI.hpp.

◆ numberOfProtares()

std::size_t GIDI::ProtareSingle::numberOfProtares ( ) const
inlinevirtual

Returns 1.

Implements GIDI::Protare.

Definition at line 4777 of file GIDI.hpp.

◆ numberOfReactions()

std::size_t GIDI::ProtareSingle::numberOfReactions ( ) const
inlinevirtual

Returns the number of reactions in the Protare.

Implements GIDI::Protare.

Definition at line 4815 of file GIDI.hpp.

◆ onlyRutherfordScatteringPresent()

bool GIDI::ProtareSingle::onlyRutherfordScatteringPresent ( ) const
inline

Returns the value of m_onlyRutherfordScatteringPresent.

Definition at line 4761 of file GIDI.hpp.

◆ orphanProduct() [1/2]

Reaction * GIDI::ProtareSingle::orphanProduct ( std::size_t a_index)
inlinevirtual

Returns the a_index - 1 orphan product reaction.

Implements GIDI::Protare.

Definition at line 4823 of file GIDI.hpp.

◆ orphanProduct() [2/2]

Reaction const * GIDI::ProtareSingle::orphanProduct ( std::size_t a_index) const
inlinevirtual

Returns the a_index - 1 orphan product reaction.

Implements GIDI::Protare.

Definition at line 4824 of file GIDI.hpp.

◆ orphanProducts() [1/2]

Suite & GIDI::ProtareSingle::orphanProducts ( )
inline

Returns a reference to the m_orphanProducts member.

Definition at line 4750 of file GIDI.hpp.

◆ orphanProducts() [2/2]

Suite const & GIDI::ProtareSingle::orphanProducts ( ) const
inline

Returns a const reference to the m_orphanProducts member.

Definition at line 4751 of file GIDI.hpp.

◆ parseEvaluatedTargetInfo()

void GIDI::ProtareSingle::parseEvaluatedTargetInfo ( HAPI::Node const & a_node)

This method parses the targetInfo node in the evaluated style node into the m_targetInfo member of this.

Parameters
a_node[in] The protare (i.e., reactionSuite) node to be parsed and used to construct a Protare.

Definition at line 1594 of file GIDI_protare.cc.

1594 {
1595
1596 m_targetInfo.parseEvaluatedTargetInfo( a_node );
1597}

Referenced by GIDI::Styles::Evaluated::Evaluated().

◆ photoAtomicIncoherentDoppler()

Suite const & GIDI::ProtareSingle::photoAtomicIncoherentDoppler ( ) const
inline

Definition at line 4771 of file GIDI.hpp.

4771{ return( m_photoAtomicIncoherentDoppler ); }

◆ productIDs()

void GIDI::ProtareSingle::productIDs ( std::set< std::string > & a_ids,
Transporting::Particles const & a_particles,
bool a_transportablesOnly ) const
virtual

Fills in a std::set with a unique list of all product indices produced by reactions and orphanProducts. If a_transportablesOnly is true, only transportable product incides are return.

Parameters
a_ids[out] The unique list of product indices.
a_particles[in] The list of particles to be transported.
a_transportablesOnly[in] If true, only transportable product indices are added in the list.

Implements GIDI::Protare.

Definition at line 625 of file GIDI_protare.cc.

625 {
626
627 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
628 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 );
629
630 if( !reaction1->active( ) ) continue;
631 reaction1->productIDs( a_ids, a_particles, a_transportablesOnly );
632 }
633
634 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
635 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
636
637 if( !reaction1->active( ) ) continue;
638 reaction1->productIDs( a_ids, a_particles, a_transportablesOnly );
639 }
640}

◆ projectileEnergyMax()

double GIDI::ProtareSingle::projectileEnergyMax ( ) const
inline

Definition at line 4744 of file GIDI.hpp.

4744{ return( m_projectileEnergyMax ); }

◆ projectileEnergyMin()

double GIDI::ProtareSingle::projectileEnergyMin ( ) const
inline

Definition at line 4743 of file GIDI.hpp.

4743{ return( m_projectileEnergyMin ); }

◆ projectileFrame()

Frame GIDI::ProtareSingle::projectileFrame ( LUPI_maybeUnused std::size_t a_index = 0) const
inline

Returns the value of the m_projectileFrame member.

Definition at line 4788 of file GIDI.hpp.

Referenced by toXMLList().

◆ protare() [1/2]

ProtareSingle * GIDI::ProtareSingle::protare ( std::size_t a_index)
virtual

Returns the pointer representing the protare (i.e., this) if a_index is 0 and nullptr otherwise.

Parameters
a_index[in] Must always be 0.
Returns
Returns the pointer representing this.

Implements GIDI::Protare.

Definition at line 583 of file GIDI_protare.cc.

583 {
584
585 if( a_index != 0 ) return( nullptr );
586 return( this );
587}

Referenced by MCGIDI::ProtareSingle::ProtareSingle().

◆ protare() [2/2]

ProtareSingle const * GIDI::ProtareSingle::protare ( std::size_t a_index) const
virtual

Returns the pointer representing the protare (i.e., this) if a_index is 0 and nullptr otherwise.

Parameters
a_index[in] Must always be 0.
Returns
Returns the pointer representing this.

Implements GIDI::Protare.

Definition at line 597 of file GIDI_protare.cc.

597 {
598
599 if( a_index != 0 ) return( nullptr );
600 return( this );
601}

◆ protareType()

ProtareType GIDI::ProtareSingle::protareType ( ) const
inlinevirtual

Returns the type of the protare.

Implements GIDI::Protare.

Definition at line 4776 of file GIDI.hpp.

◆ reaction() [1/3]

Reaction * GIDI::ProtareSingle::reaction ( std::size_t a_index)
inlinevirtual

Returns the a_index - 1 reaction.

Implements GIDI::Protare.

Definition at line 4816 of file GIDI.hpp.

◆ reaction() [2/3]

Reaction const * GIDI::ProtareSingle::reaction ( std::size_t a_index) const
inlinevirtual

Returns the a_index - 1 reaction.

Implements GIDI::Protare.

Definition at line 4817 of file GIDI.hpp.

◆ reaction() [3/3]

Reaction const * GIDI::ProtareSingle::reaction ( std::size_t a_index,
Transporting::MG const & a_settings,
ExcludeReactionsSet const & a_reactionsToExclude = ExcludeReactionsSet {} ) const
virtual

Returns the (a_index+1)th reaction or nullptr. If the indexed reaction is deactivated or exlucded, a nullptr is returned.

Parameters
a_index[in] The index of the requested reaction.
a_settings[in] Specifies the requested label.
a_reactionsToExclude[in] A list of reaction indices that are to be ignored when calculating the cross section.
Returns
The (a_index+1)th reaction or nullptr.

Implements GIDI::Protare.

Definition at line 768 of file GIDI_protare.cc.

769 {
770
771 return( reactionToMultiGroup( a_settings, a_index, a_reactionsToExclude ) );
772}

◆ reactions() [1/2]

Suite & GIDI::ProtareSingle::reactions ( )
inline

Returns a reference to the m_reactions member.

Definition at line 4748 of file GIDI.hpp.

◆ reactions() [2/2]

Suite const & GIDI::ProtareSingle::reactions ( ) const
inline

Returns a const reference to the m_reactions member.

Definition at line 4749 of file GIDI.hpp.

◆ reactionToMultiGroup()

Reaction const * GIDI::ProtareSingle::reactionToMultiGroup ( Transporting::MG const & a_settings,
std::size_t a_index,
ExcludeReactionsSet const & a_reactionsToExclude ) const

Returns a_reaction except when Rutherford scattering present and only nuclear + Coulomb interference is wanted. Otherwise returns m_nuclearPlusCoulombInterferenceOnlyReaction which may be a nullptr.

Parameters
a_settings[in] Specifies user options.
a_reaction[in] Reaction pointer to return if check passes.
Returns
Const reaction pointer or m_nuclearPlusCoulombInterferenceOnlyReaction.

Definition at line 564 of file GIDI_protare.cc.

565 {
566
567 Reaction const *reaction1 = m_reactions.get<Reaction>( a_index );
568
569 if( !reaction1->active( ) ) return( nullptr );
570 if( a_reactionsToExclude.find( a_index ) != a_reactionsToExclude.end( ) ) return( nullptr );
571
572 return( checkIf_nuclearPlusCoulombInterferenceWanted( a_settings, reaction1 ) );
573}
Reaction const * checkIf_nuclearPlusCoulombInterferenceWanted(Transporting::MG const &a_settings, Reaction const *a_reaction) const

Referenced by multiGroupAvailableEnergy(), multiGroupAvailableMomentum(), multiGroupAverageEnergy(), multiGroupAverageMomentum(), multiGroupCrossSection(), multiGroupDepositionEnergy(), multiGroupGain(), multiGroupMultiplicity(), multiGroupProductMatrix(), multiGroupQ(), and reaction().

◆ realFileName()

std::string const & GIDI::ProtareSingle::realFileName ( LUPI_maybeUnused std::size_t a_index = 0) const
inline

Returns the value of the m_realFileName member.

Definition at line 4783 of file GIDI.hpp.

◆ RutherfordScatteringPresent()

bool GIDI::ProtareSingle::RutherfordScatteringPresent ( ) const
inline

Returns the value of m_RutherfordScatteringPresent.

Definition at line 4759 of file GIDI.hpp.

◆ saveAs()

void GIDI::ProtareSingle::saveAs ( std::string const & a_fileName) const

Write this to a file in GNDS/XML format.

Parameters
a_fileName[in] Name of file to save XML lines to.

Definition at line 1536 of file GIDI_protare.cc.

1536 {
1537
1538 GUPI::WriteInfo writeInfo;
1539
1540 toXMLList( writeInfo, "" );
1541
1542 std::ofstream fileio;
1543 fileio.open( a_fileName.c_str( ) );
1544 for( std::list<std::string>::iterator iter = writeInfo.m_lines.begin( ); iter != writeInfo.m_lines.end( ); ++iter ) {
1545 fileio << *iter << std::endl;
1546 }
1547 fileio.close( );
1548}
void toXMLList(GUPI::WriteInfo &a_writeInfo, std::string const &a_indent="") const
std::list< std::string > m_lines
Definition GUPI.hpp:45

◆ setDataManager()

void GIDI::ProtareSingle::setDataManager ( HAPI::DataManager * a_dataManager)
inline

Sets the member m_dataManager to a_dataManager.

Definition at line 4735 of file GIDI.hpp.

◆ style() [1/2]

Styles::Base & GIDI::ProtareSingle::style ( std::string const & a_label)
inlinevirtual

Returns the style with label a_label.

Implements GIDI::Protare.

Definition at line 4800 of file GIDI.hpp.

◆ style() [2/2]

Styles::Base const & GIDI::ProtareSingle::style ( std::string const & a_label) const
inline

Returns the const style with label a_label.

Definition at line 4801 of file GIDI.hpp.

◆ styles() [1/2]

◆ styles() [2/2]

Styles::Suite const & GIDI::ProtareSingle::styles ( ) const
inlinevirtual

Returns a const reference to the m_styles member.

Implements GIDI::Protare.

Definition at line 4803 of file GIDI.hpp.

◆ sums() [1/2]

Sums::Sums & GIDI::ProtareSingle::sums ( )
inline

Returns a reference to the m_sums member.

Definition at line 4755 of file GIDI.hpp.

◆ sums() [2/2]

Sums::Sums const & GIDI::ProtareSingle::sums ( ) const
inline

Returns a reference to the m_sums member.

Definition at line 4756 of file GIDI.hpp.

◆ targetInfo()

TargetInfo::TargetInfo const & GIDI::ProtareSingle::targetInfo ( ) const
inline

Increments the m_numberOfLazyParsingHelperFormsReplaced member of this by 1.

Returns a const reference to the m_targetInfo member.

Definition at line 4742 of file GIDI.hpp.

◆ temperatures()

Styles::TemperatureInfos GIDI::ProtareSingle::temperatures ( ) const
virtual

Returns a list of all process temperature data. For each temeprature, the labels for its

  • heated cross section data,
  • gridded cross section data,
  • multi-group data, and
  • multi-group upscatter data.

are returned. If no data are present for a give data type (e.g., gridded cross section, multi-group upscatter), its label is an empty std::string.

Returns
The list of temperatures and their labels via an Styles::TemperatureInfos instance. The Styles::TemperatureInfos class has many (if not all) the method of a std::vector.

Implements GIDI::Protare.

Definition at line 700 of file GIDI_protare.cc.

700 {
701
702 std::size_t size( m_styles.size( ) );
703 Styles::TemperatureInfos temperature_infos;
704
705 for( std::size_t i1 = 0; i1 < size; ++i1 ) {
706 Styles::Base const *style1 = m_styles.get<Styles::Base>( i1 );
707
708 if( style1->moniker( ) == GIDI_heatedStyleChars ) {
709 PhysicalQuantity const &temperature = style1->temperature( );
710 std::string heated_cross_section( style1->label( ) );
711 std::string gridded_cross_section( "" );
712 std::string URR_probability_tables( "" );
713 std::string heated_multi_group( "" );
714 std::string Sn_elastic_upscatter( "" );
715
716 for( std::size_t i2 = 0; i2 < size; ++i2 ) {
717 Styles::Base const *style2 = m_styles.get<Styles::Base>( i2 );
718
719 if( style2->moniker( ) == GIDI_multiGroupStyleChars ) continue;
720 if( style2->temperature( ).value( ) != temperature.value( ) ) continue;
721
722 if( style2->moniker( ) == GIDI_griddedCrossSectionStyleChars ) {
723 gridded_cross_section = style2->label( ); }
724 else if( style2->moniker( ) == GIDI_URR_probabilityTablesStyleChars ) {
725 URR_probability_tables = style2->label( ); }
726 else if( style2->moniker( ) == GIDI_SnElasticUpScatterStyleChars ) {
727 Sn_elastic_upscatter = style2->label( ); }
728 else if( style2->moniker( ) == GIDI_heatedMultiGroupStyleChars ) {
729 heated_multi_group = style2->label( );
730 }
731 }
732 temperature_infos.push_back( Styles::TemperatureInfo( temperature, heated_cross_section, gridded_cross_section, URR_probability_tables,
733 heated_multi_group, Sn_elastic_upscatter ) );
734 }
735 }
736
737 std::sort( temperature_infos.begin( ), temperature_infos.end( ), sortTemperatures );
738
739 return( temperature_infos );
740}
#define GIDI_multiGroupStyleChars
Definition GIDI.hpp:246
#define GIDI_griddedCrossSectionStyleChars
Definition GIDI.hpp:251
#define GIDI_URR_probabilityTablesStyleChars
Definition GIDI.hpp:252
#define GIDI_heatedMultiGroupStyleChars
Definition GIDI.hpp:253
#define GIDI_heatedStyleChars
Definition GIDI.hpp:250
#define GIDI_SnElasticUpScatterStyleChars
Definition GIDI.hpp:254
std::vector< Styles::TemperatureInfo > TemperatureInfos
Definition GIDI.hpp:3440

◆ thresholdFactor()

double GIDI::ProtareSingle::thresholdFactor ( ) const
inlinevirtual

Returns the value of the m_thresholdFactor member.

Implements GIDI::Protare.

Definition at line 4793 of file GIDI.hpp.

◆ toXMLList()

void GIDI::ProtareSingle::toXMLList ( GUPI::WriteInfo & a_writeInfo,
std::string const & a_indent = "" ) const
virtual

Fills the argument a_writeInfo with the XML lines that represent this. Recursively enters each sub-node.

Parameters
a_writeInfo[in/out] Instance containing incremental indentation and other information and stores the appended lines.
a_indent[in] The amount to indent this node.

Reimplemented from GUPI::Ancestry.

Definition at line 1557 of file GIDI_protare.cc.

1557 {
1558
1559 std::string indent2 = a_writeInfo.incrementalIndent( a_indent );
1560 std::string header = LUPI_XML_verionEncoding;
1561 std::string attributes;
1562
1563 a_writeInfo.push_back( header );
1564
1565 attributes = a_writeInfo.addAttribute( GIDI_projectileChars, projectile( ).ID( ) );
1566 attributes += a_writeInfo.addAttribute( GIDI_targetChars, GNDS_target( ).ID( ) );
1567 attributes += a_writeInfo.addAttribute( GIDI_evaluationChars, evaluation( ) );
1568 attributes += a_writeInfo.addAttribute( GIDI_formatChars, m_formatVersion.format( ) );
1569 attributes += a_writeInfo.addAttribute( GIDI_projectileFrameChars, frameToString( projectileFrame( ) ) );
1570
1571 a_writeInfo.addNodeStarter( a_indent, moniker( ), attributes );
1572
1573 m_externalFiles.toXMLList( a_writeInfo, indent2 );
1574 m_styles.toXMLList( a_writeInfo, indent2 );
1575
1576 std::vector<std::string> pops_XMLList;
1577 m_internalPoPs.toXMLList( pops_XMLList, indent2 );
1578 for( std::vector<std::string>::iterator iter = pops_XMLList.begin( ); iter != pops_XMLList.end( ); ++iter ) a_writeInfo.push_back( *iter );
1579
1580 m_reactions.toXMLList( a_writeInfo, indent2 );
1581 m_orphanProducts.toXMLList( a_writeInfo, indent2 );
1582 m_sums.toXMLList( a_writeInfo, indent2 );
1583 m_fissionComponents.toXMLList( a_writeInfo, indent2 );
1584
1585 a_writeInfo.addNodeEnder( moniker( ) );
1586}
#define GIDI_projectileChars
Definition GIDI.hpp:413
#define GIDI_evaluationChars
Definition GIDI.hpp:415
#define GIDI_targetChars
Definition GIDI.hpp:414
#define GIDI_projectileFrameChars
Definition GIDI.hpp:419
#define GIDI_formatChars
Definition GIDI.hpp:167
#define LUPI_XML_verionEncoding
Definition LUPI.hpp:28
std::string const & evaluation(LUPI_maybeUnused std::size_t a_index=0) const
Definition GIDI.hpp:4786
Frame projectileFrame(LUPI_maybeUnused std::size_t a_index=0) const
Definition GIDI.hpp:4788
ParticleInfo const & GNDS_target() const
Definition GIDI.hpp:4547
std::string const & moniker() const
Definition GUPI.hpp:102
void push_back(std::string const &a_line)
Definition GUPI.hpp:53
void addNodeEnder(std::string const &a_moniker)
Definition GUPI.hpp:59
std::string incrementalIndent(std::string const &indent)
Definition GUPI.hpp:52
void addNodeStarter(std::string const &indent, std::string const &a_moniker, std::string const &a_attributes="")
Definition GUPI.hpp:55
std::string addAttribute(std::string const &a_name, std::string const &a_value) const
Definition GUPI.hpp:60
std::string frameToString(Frame a_frame)
Definition GIDI_misc.cc:375

Referenced by saveAs().

◆ updateReactionIndices()

void GIDI::ProtareSingle::updateReactionIndices ( std::size_t a_offset) const
virtual

Re-indexs the reactions in the reactions, orphanProducts and fissionComponents suites.

Implements GIDI::Protare.

Definition at line 798 of file GIDI_protare.cc.

798 {
799
800 for( std::size_t i1 = 0; i1 < m_reactions.size( ); ++i1 ) {
801 Reaction const *reaction1 = m_reactions.get<Reaction>( i1 + a_offset );
802
803 reaction1->setReactionIndex( i1 );
804 }
805 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
806 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
807
808 reaction1->setReactionIndex( i1 );
809 }
810 for( std::size_t i1 = 0; i1 < m_orphanProducts.size( ); ++i1 ) {
811 Reaction const *reaction1 = m_orphanProducts.get<Reaction>( i1 );
812
813 reaction1->setReactionIndex( i1 );
814 }
815}

Referenced by MCGIDI::ProtareSingle::ProtareSingle().


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