#include <G4NeutronHPInelasticData.hh>
Inheritance diagram for G4NeutronHPInelasticData:
Public Member Functions | |
G4NeutronHPInelasticData () | |
~G4NeutronHPInelasticData () | |
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 &) |
void | IgnoreOnFlightDopplerBroadening () |
void | EnableOnFlightDopplerBroadening () |
Definition at line 50 of file G4NeutronHPInelasticData.hh.
G4NeutronHPInelasticData::G4NeutronHPInelasticData | ( | ) |
Definition at line 44 of file G4NeutronHPInelasticData.cc.
References BuildPhysicsTable(), G4Neutron::Neutron(), G4VCrossSectionDataSet::SetMaxKinEnergy(), and G4VCrossSectionDataSet::SetMinKinEnergy().
00045 :G4VCrossSectionDataSet("NeutronHPInelasticXS") 00046 { 00047 00048 SetMinKinEnergy( 0*MeV ); 00049 SetMaxKinEnergy( 20*MeV ); 00050 00051 ke_cache = 0.0; 00052 xs_cache = 0.0; 00053 element_cache = NULL; 00054 material_cache = NULL; 00055 00056 onFlightDB = true; 00057 theCrossSections = 0; 00058 BuildPhysicsTable(*G4Neutron::Neutron()); 00059 }
G4NeutronHPInelasticData::~G4NeutronHPInelasticData | ( | ) |
Definition at line 61 of file G4NeutronHPInelasticData.cc.
References G4PhysicsTable::clearAndDestroy().
00062 { 00063 if ( theCrossSections != 0 ) theCrossSections->clearAndDestroy(); 00064 delete theCrossSections; 00065 }
void G4NeutronHPInelasticData::BuildPhysicsTable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 107 of file G4NeutronHPInelasticData.cc.
References G4PhysicsTable::clearAndDestroy(), G4cout, G4endl, G4Element::GetElementTable(), G4Element::GetNumberOfElements(), G4Neutron::Neutron(), and G4PhysicsTable::push_back().
Referenced by G4NeutronHPInelasticData().
00108 { 00109 if(&aP!=G4Neutron::Neutron()) 00110 throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!"); 00111 00112 //080428 00113 if ( getenv( "G4NEUTRONHP_NEGLECT_DOPPLER" ) ) 00114 { 00115 G4cout << "Find environment variable of \"G4NEUTRONHP_NEGLECT_DOPPLER\"." << G4endl; 00116 G4cout << "On the fly Doppler broadening will be neglect in the cross section calculation of inelastic scattering of neutrons (<20MeV)." << G4endl; 00117 onFlightDB = false; 00118 } 00119 00120 size_t numberOfElements = G4Element::GetNumberOfElements(); 00121 // theCrossSections = new G4PhysicsTable( numberOfElements ); 00122 // TKDB 00123 //if ( theCrossSections == 0 ) 00124 //{ theCrossSections = new G4PhysicsTable( numberOfElements ); } 00125 if ( theCrossSections == NULL ) 00126 theCrossSections = new G4PhysicsTable( numberOfElements ); 00127 else 00128 theCrossSections->clearAndDestroy(); 00129 00130 // make a PhysicsVector for each element 00131 00132 static const G4ElementTable *theElementTable = G4Element::GetElementTable(); 00133 for( size_t i=0; i<numberOfElements; ++i ) 00134 { 00135 G4PhysicsVector* physVec = G4NeutronHPData:: 00136 Instance()->MakePhysicsVector((*theElementTable)[i], this); 00137 theCrossSections->push_back(physVec); 00138 } 00139 }
void G4NeutronHPInelasticData::DumpPhysicsTable | ( | const G4ParticleDefinition & | ) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 141 of file G4NeutronHPInelasticData.cc.
References G4cout, G4endl, G4Element::GetElementTable(), G4VCrossSectionDataSet::GetName(), G4Element::GetNumberOfElements(), and G4Neutron::Neutron().
00142 { 00143 if(&aP!=G4Neutron::Neutron()) 00144 throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!"); 00145 00146 // 00147 // Dump element based cross section 00148 // range 10e-5 eV to 20 MeV 00149 // 10 point per decade 00150 // in barn 00151 // 00152 00153 G4cout << G4endl; 00154 G4cout << G4endl; 00155 G4cout << "Inelastic Cross Section of Neutron HP"<< G4endl; 00156 G4cout << "(Pointwise cross-section at 0 Kelvin.)" << G4endl; 00157 G4cout << G4endl; 00158 G4cout << "Name of Element" << G4endl; 00159 G4cout << "Energy[eV] XS[barn]" << G4endl; 00160 G4cout << G4endl; 00161 00162 size_t numberOfElements = G4Element::GetNumberOfElements(); 00163 static const G4ElementTable *theElementTable = G4Element::GetElementTable(); 00164 00165 for ( size_t i = 0 ; i < numberOfElements ; ++i ) 00166 { 00167 00168 G4cout << (*theElementTable)[i]->GetName() << G4endl; 00169 00170 G4int ie = 0; 00171 00172 for ( ie = 0 ; ie < 130 ; ie++ ) 00173 { 00174 G4double eKinetic = 1.0e-5 * std::pow ( 10.0 , ie/10.0 ) *eV; 00175 G4bool outOfRange = false; 00176 00177 if ( eKinetic < 20*MeV ) 00178 { 00179 G4cout << eKinetic/eV << " " << (*((*theCrossSections)(i))).GetValue(eKinetic, outOfRange)/barn << G4endl; 00180 } 00181 00182 } 00183 00184 G4cout << G4endl; 00185 } 00186 00187 //G4cout << "G4NeutronHPInelasticData::DumpPhysicsTable still to be implemented"<<G4endl; 00188 }
void G4NeutronHPInelasticData::EnableOnFlightDopplerBroadening | ( | ) | [inline] |
G4double G4NeutronHPInelasticData::GetCrossSection | ( | const G4DynamicParticle * | , | |
const G4Element * | , | |||
G4double | aT | |||
) |
Definition at line 193 of file G4NeutronHPInelasticData.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(), G4Element::GetZ(), G4ReactionProduct::Lorentz(), and G4Neutron::Neutron().
Referenced by GetIsoCrossSection().
00194 { 00195 G4double result = 0; 00196 G4bool outOfRange; 00197 G4int index = anE->GetIndex(); 00198 00199 // prepare neutron 00200 G4double eKinetic = aP->GetKineticEnergy(); 00201 00202 // T. K. 00203 //if ( getenv( "G4NEUTRONHP_NEGLECT_DOPPLER" ) ) 00204 //080428 00205 if ( !onFlightDB ) 00206 { 00207 G4double factor = 1.0; 00208 if ( eKinetic < aT * k_Boltzmann ) 00209 { 00210 // below 0.1 eV neutrons 00211 // Have to do some, but now just igonre. 00212 // Will take care after performance check. 00213 // factor = factor * targetV; 00214 } 00215 return ( (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) )* factor; 00216 } 00217 00218 G4ReactionProduct theNeutron( aP->GetDefinition() ); 00219 theNeutron.SetMomentum( aP->GetMomentum() ); 00220 theNeutron.SetKineticEnergy( eKinetic ); 00221 00222 // prepare thermal nucleus 00223 G4Nucleus aNuc; 00224 G4double eps = 0.0001; 00225 G4double theA = anE->GetN(); 00226 G4double theZ = anE->GetZ(); 00227 G4double eleMass; 00228 eleMass = ( G4NucleiProperties::GetNuclearMass(static_cast<G4int>(theA+eps), static_cast<G4int>(theZ+eps)) 00229 ) / G4Neutron::Neutron()->GetPDGMass(); 00230 00231 G4ReactionProduct boosted; 00232 G4double aXsection; 00233 00234 // MC integration loop 00235 G4int counter = 0; 00236 G4int failCount = 0; 00237 G4double buffer = 0; 00238 G4int size = G4int(std::max(10., aT/60*kelvin)); 00239 G4ThreeVector neutronVelocity = 1./G4Neutron::Neutron()->GetPDGMass()*theNeutron.GetMomentum(); 00240 G4double neutronVMag = neutronVelocity.mag(); 00241 00242 while(counter == 0 || std::abs(buffer-result/std::max(1,counter)) > 0.01*buffer) 00243 { 00244 if(counter) buffer = result/counter; 00245 while (counter<size) 00246 { 00247 counter ++; 00248 G4ReactionProduct aThermalNuc = aNuc.GetThermalNucleus(eleMass, aT); 00249 boosted.Lorentz(theNeutron, aThermalNuc); 00250 G4double theEkin = boosted.GetKineticEnergy(); 00251 aXsection = (*((*theCrossSections)(index))).GetValue(theEkin, outOfRange); 00252 if(aXsection <0) 00253 { 00254 if(failCount<1000) 00255 { 00256 failCount++; 00257 counter--; 00258 continue; 00259 } 00260 else 00261 { 00262 aXsection = 0; 00263 } 00264 } 00265 // velocity correction. 00266 G4ThreeVector targetVelocity = 1./aThermalNuc.GetMass()*aThermalNuc.GetMomentum(); 00267 aXsection *= (targetVelocity-neutronVelocity).mag()/neutronVMag; 00268 result += aXsection; 00269 } 00270 size += size; 00271 } 00272 result /= counter; 00273 /* 00274 // Checking impact of G4NEUTRONHP_NEGLECT_DOPPLER 00275 G4cout << " result " << result << " " 00276 << (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) << " " 00277 << (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) /result << G4endl; 00278 */ 00279 return result; 00280 }
G4double G4NeutronHPInelasticData::GetIsoCrossSection | ( | const G4DynamicParticle * | , | |
G4int | , | |||
G4int | , | |||
const G4Isotope * | , | |||
const G4Element * | , | |||
const G4Material * | ||||
) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 80 of file G4NeutronHPInelasticData.cc.
References GetCrossSection(), and G4DynamicParticle::GetKineticEnergy().
00085 { 00086 if ( dp->GetKineticEnergy() == ke_cache && element == element_cache && material == material_cache ) return xs_cache; 00087 00088 ke_cache = dp->GetKineticEnergy(); 00089 element_cache = element; 00090 material_cache = material; 00091 G4double xs = GetCrossSection( dp , element , material->GetTemperature() ); 00092 xs_cache = xs; 00093 return xs; 00094 //return GetCrossSection( dp , element , material->GetTemperature() ); 00095 }
void G4NeutronHPInelasticData::IgnoreOnFlightDopplerBroadening | ( | ) | [inline] |
G4bool G4NeutronHPInelasticData::IsIsoApplicable | ( | const G4DynamicParticle * | , | |
G4int | , | |||
G4int | , | |||
const G4Element * | , | |||
const G4Material * | ||||
) | [virtual] |
Reimplemented from G4VCrossSectionDataSet.
Definition at line 67 of file G4NeutronHPInelasticData.cc.
References G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetKineticEnergy(), G4VCrossSectionDataSet::GetMaxKinEnergy(), G4VCrossSectionDataSet::GetMinKinEnergy(), and G4Neutron::Neutron().
00071 { 00072 G4double eKin = dp->GetKineticEnergy(); 00073 if ( eKin > GetMaxKinEnergy() 00074 || eKin < GetMinKinEnergy() 00075 || dp->GetDefinition() != G4Neutron::Neutron() ) return false; 00076 00077 return true; 00078 }