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