#include <G4NeutronHPorLElasticData.hh>
Inheritance diagram for G4NeutronHPorLElasticData:
Public Member Functions | |
G4NeutronHPorLElasticData () | |
G4NeutronHPorLElasticData (G4NeutronHPChannel *, std::set< G4String > *) | |
~G4NeutronHPorLElasticData () | |
G4bool | IsIsoApplicable (const G4DynamicParticle *, G4int, G4int, const G4Element *, const G4Material *) |
G4double | GetIsoCrossSection (const G4DynamicParticle *, G4int, G4int, const G4Isotope *, const G4Element *, const G4Material *) |
G4double | GetCrossSection (const G4DynamicParticle *, const G4Element *, G4double aT) |
void | BuildPhysicsTable (const G4ParticleDefinition &) |
void | DumpPhysicsTable (const G4ParticleDefinition &) |
Definition at line 55 of file G4NeutronHPorLElasticData.hh.
G4NeutronHPorLElasticData::G4NeutronHPorLElasticData | ( | ) |
Definition at line 42 of file G4NeutronHPorLElasticData.cc.
References G4VCrossSectionDataSet::SetMaxKinEnergy(), and G4VCrossSectionDataSet::SetMinKinEnergy().
00043 { 00044 SetMinKinEnergy( 0*MeV ); 00045 SetMaxKinEnergy( 20*MeV ); 00046 00047 ke_cache = 0.0; 00048 xs_cache = 0.0; 00049 element_cache = NULL; 00050 material_cache = NULL; 00051 // BuildPhysicsTable(*G4Neutron::Neutron()); 00052 }
G4NeutronHPorLElasticData::G4NeutronHPorLElasticData | ( | G4NeutronHPChannel * | , | |
std::set< G4String > * | ||||
) |
Definition at line 91 of file G4NeutronHPorLElasticData.cc.
References G4VCrossSectionDataSet::SetMaxKinEnergy(), and G4VCrossSectionDataSet::SetMinKinEnergy().
00092 :G4VCrossSectionDataSet("NeutronHPorLElasticXS") 00093 { 00094 theElasticChannel = pChannel; 00095 unavailable_elements = pSet; 00096 00097 SetMinKinEnergy( 0*MeV ); 00098 SetMaxKinEnergy( 20*MeV ); 00099 00100 ke_cache = 0.0; 00101 xs_cache = 0.0; 00102 element_cache = NULL; 00103 material_cache = NULL; 00104 }
G4NeutronHPorLElasticData::~G4NeutronHPorLElasticData | ( | ) |
void G4NeutronHPorLElasticData::BuildPhysicsTable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 118 of file G4NeutronHPorLElasticData.cc.
References G4Neutron::Neutron().
00119 { 00120 if( &aP!=G4Neutron::Neutron() ) 00121 throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!"); 00122 }
void G4NeutronHPorLElasticData::DumpPhysicsTable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 126 of file G4NeutronHPorLElasticData.cc.
References G4Neutron::Neutron().
00127 { 00128 if(&aP!=G4Neutron::Neutron()) 00129 throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!"); 00130 // G4cout << "G4NeutronHPorLElasticData::DumpPhysicsTable still to be implemented"<<G4endl; 00131 }
G4double G4NeutronHPorLElasticData::GetCrossSection | ( | const G4DynamicParticle * | , | |
const G4Element * | , | |||
G4double | aT | |||
) |
Definition at line 141 of file G4NeutronHPorLElasticData.cc.
References buffer, G4DynamicParticle::GetDefinition(), G4Element::GetIndex(), G4ReactionProduct::GetKineticEnergy(), G4DynamicParticle::GetKineticEnergy(), G4ReactionProduct::GetMass(), G4ReactionProduct::GetMomentum(), G4DynamicParticle::GetMomentum(), G4Element::GetN(), G4NucleiProperties::GetNuclearMass(), G4ParticleDefinition::GetPDGMass(), G4Nucleus::GetThermalNucleus(), G4NeutronHPChannel::GetXsec(), G4Element::GetZ(), G4ReactionProduct::Lorentz(), and G4Neutron::Neutron().
Referenced by GetIsoCrossSection().
00142 { 00143 00144 //G4cout << "Choice G4NeutronHPorLElasticData for element " << anE->GetName() << G4endl; 00145 G4double result = 0; 00146 // G4bool outOfRange; 00147 G4int index = anE->GetIndex(); 00148 00149 // prepare neutron 00150 G4double eKinetic = aP->GetKineticEnergy(); 00151 G4ReactionProduct theNeutron( aP->GetDefinition() ); 00152 theNeutron.SetMomentum( aP->GetMomentum() ); 00153 theNeutron.SetKineticEnergy( eKinetic ); 00154 00155 // prepare thermal nucleus 00156 G4Nucleus aNuc; 00157 G4double eps = 0.0001; 00158 G4double theA = anE->GetN(); 00159 G4double theZ = anE->GetZ(); 00160 G4double eleMass; 00161 eleMass = ( G4NucleiProperties::GetNuclearMass(static_cast<G4int>(theA+eps), static_cast<G4int>(theZ+eps)) 00162 ) / G4Neutron::Neutron()->GetPDGMass(); 00163 00164 G4ReactionProduct boosted; 00165 G4double aXsection; 00166 00167 // MC integration loop 00168 G4int counter = 0; 00169 G4double buffer = 0; 00170 G4int size = G4int(std::max(10., aT/60*kelvin)); 00171 G4ThreeVector neutronVelocity = 1./G4Neutron::Neutron()->GetPDGMass()*theNeutron.GetMomentum(); 00172 G4double neutronVMag = neutronVelocity.mag(); 00173 while(counter == 0 || std::abs(buffer-result/std::max(1,counter)) > 0.03*buffer) 00174 { 00175 if(counter) buffer = result/counter; 00176 while (counter<size) 00177 { 00178 counter ++; 00179 G4ReactionProduct aThermalNuc = aNuc.GetThermalNucleus(eleMass, aT); 00180 boosted.Lorentz(theNeutron, aThermalNuc); 00181 G4double theEkin = boosted.GetKineticEnergy(); 00182 //aXsection = (*((*theCrossSections)(index))).GetValue(theEkin, outOfRange); 00183 aXsection = theElasticChannel[index].GetXsec( theEkin ); 00184 // velocity correction. 00185 G4ThreeVector targetVelocity = 1./aThermalNuc.GetMass()*aThermalNuc.GetMomentum(); 00186 aXsection *= (targetVelocity-neutronVelocity).mag()/neutronVMag; 00187 result += aXsection; 00188 } 00189 size += size; 00190 } 00191 result /= counter; 00192 return result; 00193 }
G4double G4NeutronHPorLElasticData::GetIsoCrossSection | ( | const G4DynamicParticle * | , | |
G4int | , | |||
G4int | , | |||
const G4Isotope * | , | |||
const G4Element * | , | |||
const G4Material * | ||||
) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 73 of file G4NeutronHPorLElasticData.cc.
References GetCrossSection(), and G4DynamicParticle::GetKineticEnergy().
00078 { 00079 if ( dp->GetKineticEnergy() == ke_cache && element == element_cache && material == material_cache ) return xs_cache; 00080 00081 ke_cache = dp->GetKineticEnergy(); 00082 element_cache = element; 00083 material_cache = material; 00084 G4double xs = GetCrossSection( dp , element , material->GetTemperature() ); 00085 xs_cache = xs; 00086 return xs; 00087 //return GetCrossSection( dp , element , material->GetTemperature() ); 00088 }
G4bool G4NeutronHPorLElasticData::IsIsoApplicable | ( | const G4DynamicParticle * | , | |
G4int | , | |||
G4int | , | |||
const G4Element * | , | |||
const G4Material * | ||||
) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 59 of file G4NeutronHPorLElasticData.cc.
References G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetKineticEnergy(), G4VCrossSectionDataSet::GetMaxKinEnergy(), G4VCrossSectionDataSet::GetMinKinEnergy(), G4Element::GetName(), and G4Neutron::Neutron().
00063 { 00064 G4double eKin = dp->GetKineticEnergy(); 00065 if ( eKin > GetMaxKinEnergy() 00066 || eKin < GetMinKinEnergy() 00067 || dp->GetDefinition() != G4Neutron::Neutron() ) return false; 00068 if ( unavailable_elements->find( element->GetName() ) != unavailable_elements->end() ) return false; 00069 00070 return true; 00071 }