Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4NeutronHPElasticData Class Reference

#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 ()
 
- Public Member Functions inherited from G4VCrossSectionDataSet
 G4VCrossSectionDataSet (const G4String &nam="")
 
virtual ~G4VCrossSectionDataSet ()
 
virtual G4bool IsElementApplicable (const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
 
virtual G4bool IsIsoApplicable (const G4DynamicParticle *, G4int Z, G4int A, const G4Element *elm=0, const G4Material *mat=0)
 
G4double GetCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
G4double ComputeCrossSection (const G4DynamicParticle *, const G4Element *, const G4Material *mat=0)
 
virtual G4double GetElementCrossSection (const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
 
virtual G4double GetIsoCrossSection (const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
 
virtual G4IsotopeSelectIsotope (const G4Element *, G4double kinEnergy)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &)
 
virtual void CrossSectionDescription (std::ostream &) const
 
void SetVerboseLevel (G4int value)
 
G4double GetMinKinEnergy () const
 
void SetMinKinEnergy (G4double value)
 
G4double GetMaxKinEnergy () const
 
void SetMaxKinEnergy (G4double value)
 
const G4StringGetName () const
 

Additional Inherited Members

- Protected Member Functions inherited from G4VCrossSectionDataSet
void SetName (const G4String &)
 
- Protected Attributes inherited from G4VCrossSectionDataSet
G4int verboseLevel
 

Detailed Description

Definition at line 50 of file G4NeutronHPElasticData.hh.

Constructor & Destructor Documentation

◆ G4NeutronHPElasticData()

G4NeutronHPElasticData::G4NeutronHPElasticData ( )

Definition at line 44 of file G4NeutronHPElasticData.cc.

45:G4VCrossSectionDataSet("NeutronHPElasticXS")
46{
47 SetMinKinEnergy( 0*MeV );
48 SetMaxKinEnergy( 20*MeV );
49
50 ke_cache = 0.0;
51 xs_cache = 0.0;
52 element_cache = NULL;
53 material_cache = NULL;
54
55 theCrossSections = 0;
56 onFlightDB = true;
58}
void BuildPhysicsTable(const G4ParticleDefinition &)
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
void SetMaxKinEnergy(G4double value)
void SetMinKinEnergy(G4double value)

◆ ~G4NeutronHPElasticData()

G4NeutronHPElasticData::~G4NeutronHPElasticData ( )

Definition at line 60 of file G4NeutronHPElasticData.cc.

61{
62 if ( theCrossSections != 0 ) theCrossSections->clearAndDestroy();
63 delete theCrossSections;
64}
void clearAndDestroy()

Member Function Documentation

◆ BuildPhysicsTable()

void G4NeutronHPElasticData::BuildPhysicsTable ( const G4ParticleDefinition aP)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 106 of file G4NeutronHPElasticData.cc.

107{
108
109 if(&aP!=G4Neutron::Neutron())
110 throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!");
111
112//080428
113 if ( getenv( "G4NEUTRONHP_NEGLECT_DOPPLER" ) )
114 {
115 G4cout << "Find environment variable of \"G4NEUTRONHP_NEGLECT_DOPPLER\"." << G4endl;
116 G4cout << "On the fly Doppler broadening will be neglect in the cross section calculation of elastic scattering of neutrons (<20MeV)." << G4endl;
117 onFlightDB = false;
118 }
119
120 size_t numberOfElements = G4Element::GetNumberOfElements();
121// TKDB
122 //if ( theCrossSections == 0 ) theCrossSections = new G4PhysicsTable( numberOfElements );
123 if ( theCrossSections == NULL )
124 theCrossSections = new G4PhysicsTable( numberOfElements );
125 else
126 theCrossSections->clearAndDestroy();
127
128 // make a PhysicsVector for each element
129
130 static const G4ElementTable *theElementTable = G4Element::GetElementTable();
131 for( size_t i=0; i<numberOfElements; ++i )
132 {
134 Instance()->MakePhysicsVector((*theElementTable)[i], this);
135 theCrossSections->push_back(physVec);
136 }
137}
std::vector< G4Element * > G4ElementTable
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
static size_t GetNumberOfElements()
Definition: G4Element.cc:406
static const G4ElementTable * GetElementTable()
Definition: G4Element.cc:399
G4PhysicsVector * MakePhysicsVector(G4Element *thE, G4NeutronHPFissionData *theP)
static G4NeutronHPData * Instance()
void push_back(G4PhysicsVector *)

Referenced by G4NeutronHPElasticData().

◆ DumpPhysicsTable()

void G4NeutronHPElasticData::DumpPhysicsTable ( const G4ParticleDefinition aP)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 139 of file G4NeutronHPElasticData.cc.

140{
141 if(&aP!=G4Neutron::Neutron())
142 throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!");
143
144//
145// Dump element based cross section
146// range 10e-5 eV to 20 MeV
147// 10 point per decade
148// in barn
149//
150
151 G4cout << G4endl;
152 G4cout << G4endl;
153 G4cout << "Elastic Cross Section of Neutron HP"<< G4endl;
154 G4cout << "(Pointwise cross-section at 0 Kelvin.)" << G4endl;
155 G4cout << G4endl;
156 G4cout << "Name of Element" << G4endl;
157 G4cout << "Energy[eV] XS[barn]" << G4endl;
158 G4cout << G4endl;
159
160 size_t numberOfElements = G4Element::GetNumberOfElements();
161 static const G4ElementTable *theElementTable = G4Element::GetElementTable();
162
163 for ( size_t i = 0 ; i < numberOfElements ; ++i )
164 {
165
166 G4cout << (*theElementTable)[i]->GetName() << G4endl;
167
168 G4int ie = 0;
169
170 for ( ie = 0 ; ie < 130 ; ie++ )
171 {
172 G4double eKinetic = 1.0e-5 * std::pow ( 10.0 , ie/10.0 ) *eV;
173 G4bool outOfRange = false;
174
175 if ( eKinetic < 20*MeV )
176 {
177 G4cout << eKinetic/eV << " " << (*((*theCrossSections)(i))).GetValue(eKinetic, outOfRange)/barn << G4endl;
178 }
179
180 }
181
182 G4cout << G4endl;
183 }
184
185
186// G4cout << "G4NeutronHPElasticData::DumpPhysicsTable still to be implemented"<<G4endl;
187}
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
bool G4bool
Definition: G4Types.hh:67

◆ EnableOnFlightDopplerBroadening()

void G4NeutronHPElasticData::EnableOnFlightDopplerBroadening ( )
inline

Definition at line 83 of file G4NeutronHPElasticData.hh.

83{ onFlightDB = true; };

◆ GetCrossSection()

G4double G4NeutronHPElasticData::GetCrossSection ( const G4DynamicParticle aP,
const G4Element anE,
G4double  aT 
)

Definition at line 194 of file G4NeutronHPElasticData.cc.

196{
197 G4double result = 0;
198 G4bool outOfRange;
199 G4int index = anE->GetIndex();
200
201 // prepare neutron
202 G4double eKinetic = aP->GetKineticEnergy();
203
204 // T. K.
205// if ( getenv( "G4NEUTRONHP_NEGLECT_DOPPLER" ) )
206//080428
207 if ( !onFlightDB )
208 {
209 G4double factor = 1.0;
210 if ( eKinetic < aT * k_Boltzmann )
211 {
212 // below 0.1 eV neutrons
213 // Have to do some, but now just igonre.
214 // Will take care after performance check.
215 // factor = factor * targetV;
216 }
217 return ( (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) )* factor;
218 }
219
220 G4ReactionProduct theNeutron( aP->GetDefinition() );
221 theNeutron.SetMomentum( aP->GetMomentum() );
222 theNeutron.SetKineticEnergy( eKinetic );
223
224 // prepare thermal nucleus
225 G4Nucleus aNuc;
226 G4double eps = 0.0001;
227 G4double theA = anE->GetN();
228 G4double theZ = anE->GetZ();
229 G4double eleMass;
230
231
232 eleMass = ( G4NucleiProperties::GetNuclearMass( static_cast<G4int>(theA+eps) , static_cast<G4int>(theZ+eps) )
234
235 G4ReactionProduct boosted;
236 G4double aXsection;
237
238 // MC integration loop
239 G4int counter = 0;
240 G4double buffer = 0;
241 G4int size = G4int(std::max(10., aT/60*kelvin));
242 G4ThreeVector neutronVelocity = 1./G4Neutron::Neutron()->GetPDGMass()*theNeutron.GetMomentum();
243 G4double neutronVMag = neutronVelocity.mag();
244
245 while(counter == 0 || std::abs(buffer-result/std::max(1,counter)) > 0.03*buffer)
246 {
247 if(counter) buffer = result/counter;
248 while (counter<size)
249 {
250 counter ++;
251 G4ReactionProduct aThermalNuc = aNuc.GetThermalNucleus(eleMass, aT);
252 boosted.Lorentz(theNeutron, aThermalNuc);
253 G4double theEkin = boosted.GetKineticEnergy();
254 aXsection = (*((*theCrossSections)(index))).GetValue(theEkin, outOfRange);
255 // velocity correction.
256 G4ThreeVector targetVelocity = 1./aThermalNuc.GetMass()*aThermalNuc.GetMomentum();
257 aXsection *= (targetVelocity-neutronVelocity).mag()/neutronVMag;
258 result += aXsection;
259 }
260 size += size;
261 }
262 result /= counter;
263/*
264 // Checking impact of G4NEUTRONHP_NEGLECT_DOPPLER
265 G4cout << " result " << result << " "
266 << (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) << " "
267 << (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) /result << G4endl;
268*/
269 return result;
270}
double mag() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
G4ThreeVector GetMomentum() const
G4double GetZ() const
Definition: G4Element.hh:131
size_t GetIndex() const
Definition: G4Element.hh:182
G4double GetN() const
Definition: G4Element.hh:134
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4ReactionProduct GetThermalNucleus(G4double aMass, G4double temp=-1) const
Definition: G4Nucleus.cc:130
G4double GetKineticEnergy() const
G4ThreeVector GetMomentum() const
void Lorentz(const G4ReactionProduct &p1, const G4ReactionProduct &p2)
G4double GetMass() const
#define buffer
Definition: xmlparse.cc:611

Referenced by GetIsoCrossSection().

◆ GetIsoCrossSection()

G4double G4NeutronHPElasticData::GetIsoCrossSection ( const G4DynamicParticle dp,
G4int  ,
G4int  ,
const G4Isotope ,
const G4Element element,
const G4Material material 
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 80 of file G4NeutronHPElasticData.cc.

85{
86 if ( dp->GetKineticEnergy() == ke_cache && element == element_cache && material == material_cache ) return xs_cache;
87
88 ke_cache = dp->GetKineticEnergy();
89 element_cache = element;
90 material_cache = material;
91 G4double xs = GetCrossSection( dp , element , material->GetTemperature() );
92 xs_cache = xs;
93 return xs;
94}
G4double GetTemperature() const
Definition: G4Material.hh:181
G4double GetCrossSection(const G4DynamicParticle *, const G4Element *, G4double aT)

◆ IgnoreOnFlightDopplerBroadening()

void G4NeutronHPElasticData::IgnoreOnFlightDopplerBroadening ( )
inline

Definition at line 82 of file G4NeutronHPElasticData.hh.

82{ onFlightDB = false; };

◆ IsIsoApplicable()

G4bool G4NeutronHPElasticData::IsIsoApplicable ( const G4DynamicParticle dp,
G4int  ,
G4int  ,
const G4Element ,
const G4Material  
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 66 of file G4NeutronHPElasticData.cc.

70{
71
72 G4double eKin = dp->GetKineticEnergy();
73 if ( eKin > GetMaxKinEnergy()
74 || eKin < GetMinKinEnergy()
75 || dp->GetDefinition() != G4Neutron::Neutron() ) return false;
76
77 return true;
78}

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