Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ParticleHPElasticData Class Reference

#include <G4ParticleHPElasticData.hh>

Inheritance diagram for G4ParticleHPElasticData:
Collaboration diagram for G4ParticleHPElasticData:

Public Member Functions

 G4ParticleHPElasticData ()
 
 ~G4ParticleHPElasticData ()
 
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 ()
 
void SetVerboseLevel (G4int)
 
G4int GetVerboseLevel () const
 
virtual void CrossSectionDescription (std::ostream &) const
 
- Public Member Functions inherited from G4VCrossSectionDataSet
 G4VCrossSectionDataSet (const G4String &nam="")
 
virtual ~G4VCrossSectionDataSet ()
 
virtual G4bool IsElementApplicable (const G4DynamicParticle *, G4int Z, 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 G4IsotopeSelectIsotope (const G4Element *, G4double kinEnergy)
 
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 51 of file G4ParticleHPElasticData.hh.

Constructor & Destructor Documentation

G4ParticleHPElasticData::G4ParticleHPElasticData ( )

Definition at line 49 of file G4ParticleHPElasticData.cc.

50 :G4VCrossSectionDataSet("NeutronHPElasticXS")
51 {
52  SetMinKinEnergy( 0*MeV );
53  SetMaxKinEnergy( 20*MeV );
54 
55  theCrossSections = 0;
56  onFlightDB = true;
57  instanceOfWorker = false;
59  instanceOfWorker = true;
60  }
61  element_cache = NULL;
62  material_cache = NULL;
63  ke_cache = 0.0;
64  xs_cache = 0.0;
65 // BuildPhysicsTable( *G4Neutron::Neutron() );
66 }
G4VCrossSectionDataSet(const G4String &nam="")
void SetMinKinEnergy(G4double value)
G4bool IsWorkerThread()
Definition: G4Threading.cc:145
void SetMaxKinEnergy(G4double value)
static constexpr double MeV
Definition: G4SIunits.hh:214

Here is the call graph for this function:

G4ParticleHPElasticData::~G4ParticleHPElasticData ( )

Definition at line 68 of file G4ParticleHPElasticData.cc.

69 {
70  if ( theCrossSections != NULL && instanceOfWorker != true ) {
71  theCrossSections->clearAndDestroy();
72  delete theCrossSections;
73  theCrossSections = NULL;
74  }
75 }
void clearAndDestroy()

Here is the call graph for this function:

Member Function Documentation

void G4ParticleHPElasticData::BuildPhysicsTable ( const G4ParticleDefinition aP)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 117 of file G4ParticleHPElasticData.cc.

118 {
119 
120  if(&aP!=G4Neutron::Neutron())
121  throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!");
122 
123 //080428
124  if ( G4ParticleHPManager::GetInstance()->GetNeglectDoppler() )
125  {
126  G4cout << "Find a flag of \"G4NEUTRONHP_NEGLECT_DOPPLER\"." << G4endl;
127  G4cout << "On the fly Doppler broadening will be neglect in the cross section calculation of elastic scattering of neutrons (<20MeV)." << G4endl;
128  onFlightDB = false;
129  }
130 
131  if ( G4Threading::IsWorkerThread() ) {
133  return;
134  }
135 
136  size_t numberOfElements = G4Element::GetNumberOfElements();
137 // TKDB
138  //if ( theCrossSections == 0 ) theCrossSections = new G4PhysicsTable( numberOfElements );
139  if ( theCrossSections == NULL )
140  theCrossSections = new G4PhysicsTable( numberOfElements );
141  else
142  theCrossSections->clearAndDestroy();
143 
144  // make a PhysicsVector for each element
145 
146  static G4ThreadLocal G4ElementTable *theElementTable = 0 ; if (!theElementTable) theElementTable= G4Element::GetElementTable();
147  for( size_t i=0; i<numberOfElements; ++i )
148  {
150  Instance(G4Neutron::Neutron())->MakePhysicsVector((*theElementTable)[i], this);
151  theCrossSections->push_back(physVec);
152  }
153 
155 }
static G4ParticleHPManager * GetInstance()
void push_back(G4PhysicsVector *)
#define G4ThreadLocal
Definition: tls.hh:89
void RegisterElasticCrossSections(G4PhysicsTable *val)
G4GLOB_DLL std::ostream G4cout
static size_t GetNumberOfElements()
Definition: G4Element.cc:405
G4PhysicsTable * GetElasticCrossSections()
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
G4bool IsWorkerThread()
Definition: G4Threading.cc:145
G4PhysicsVector * MakePhysicsVector(G4Element *thE, G4ParticleHPFissionData *theP)
static G4ParticleHPData * Instance(G4ParticleDefinition *projectile)
#define G4endl
Definition: G4ios.hh:61
std::vector< G4Element * > G4ElementTable
static G4ElementTable * GetElementTable()
Definition: G4Element.cc:398
void clearAndDestroy()

Here is the call graph for this function:

void G4ParticleHPElasticData::CrossSectionDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 299 of file G4ParticleHPElasticData.cc.

300 {
301  outFile << "High Precision cross data based on Evaluated Nuclear Data Files (ENDF) for elastic reaction of neutrons below 20MeV\n" ;
302 }
void G4ParticleHPElasticData::DumpPhysicsTable ( const G4ParticleDefinition aP)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 157 of file G4ParticleHPElasticData.cc.

158 {
159  if(&aP!=G4Neutron::Neutron())
160  throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!");
161 
162 //
163 // Dump element based cross section
164 // range 10e-5 eV to 20 MeV
165 // 10 point per decade
166 // in barn
167 //
168 
169  G4cout << G4endl;
170  G4cout << G4endl;
171  G4cout << "Elastic Cross Section of Neutron HP"<< G4endl;
172  G4cout << "(Pointwise cross-section at 0 Kelvin.)" << G4endl;
173  G4cout << G4endl;
174  G4cout << "Name of Element" << G4endl;
175  G4cout << "Energy[eV] XS[barn]" << G4endl;
176  G4cout << G4endl;
177 
178  size_t numberOfElements = G4Element::GetNumberOfElements();
179  static G4ThreadLocal G4ElementTable *theElementTable = 0 ; if (!theElementTable) theElementTable= G4Element::GetElementTable();
180 
181  for ( size_t i = 0 ; i < numberOfElements ; ++i )
182  {
183 
184  G4cout << (*theElementTable)[i]->GetName() << G4endl;
185 
186  G4int ie = 0;
187 
188  for ( ie = 0 ; ie < 130 ; ie++ )
189  {
190  G4double eKinetic = 1.0e-5 * G4Pow::GetInstance()->powA ( 10.0 , ie/10.0 ) *eV;
191  G4bool outOfRange = false;
192 
193  if ( eKinetic < 20*MeV )
194  {
195  G4cout << eKinetic/eV << " " << (*((*theCrossSections)(i))).GetValue(eKinetic, outOfRange)/barn << G4endl;
196  }
197 
198  }
199 
200  G4cout << G4endl;
201  }
202 
203 
204 // G4cout << "G4ParticleHPElasticData::DumpPhysicsTable still to be implemented"<<G4endl;
205 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
#define G4ThreadLocal
Definition: tls.hh:89
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
static size_t GetNumberOfElements()
Definition: G4Element.cc:405
bool G4bool
Definition: G4Types.hh:79
static constexpr double eV
Definition: G4SIunits.hh:215
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
std::vector< G4Element * > G4ElementTable
double G4double
Definition: G4Types.hh:76
static constexpr double barn
Definition: G4SIunits.hh:105
static G4ElementTable * GetElementTable()
Definition: G4Element.cc:398

Here is the call graph for this function:

void G4ParticleHPElasticData::EnableOnFlightDopplerBroadening ( )
inline

Definition at line 84 of file G4ParticleHPElasticData.hh.

84 { onFlightDB = true; };
G4double G4ParticleHPElasticData::GetCrossSection ( const G4DynamicParticle aP,
const G4Element anE,
G4double  aT 
)

Definition at line 213 of file G4ParticleHPElasticData.cc.

214 {
215  G4double result = 0;
216  G4bool outOfRange;
217  G4int index = anE->GetIndex();
218 
219  // prepare neutron
220  G4double eKinetic = aP->GetKineticEnergy();
221 
222  if ( !onFlightDB )
223  {
224  //NEGLECT_DOPPLER_B.
225  G4double factor = 1.0;
226  if ( eKinetic < aT * k_Boltzmann )
227  {
228  // below 0.1 eV neutrons
229  // Have to do some, but now just igonre.
230  // Will take care after performance check.
231  // factor = factor * targetV;
232  }
233  return ( (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) )* factor;
234  }
235 
236  G4ReactionProduct theNeutron( aP->GetDefinition() );
237  theNeutron.SetMomentum( aP->GetMomentum() );
238  theNeutron.SetKineticEnergy( eKinetic );
239 
240  // prepare thermal nucleus
241  G4Nucleus aNuc;
242  G4double eps = 0.0001;
243  G4double theA = anE->GetN();
244  G4double theZ = anE->GetZ();
245  G4double eleMass;
246 
247 
248  eleMass = ( G4NucleiProperties::GetNuclearMass( static_cast<G4int>(theA+eps) , static_cast<G4int>(theZ+eps) )
250 
251  G4ReactionProduct boosted;
252  G4double aXsection;
253 
254  // MC integration loop
255  G4int counter = 0;
256  G4double buffer = 0;
257  G4int size = G4int(std::max(10., aT/60*kelvin));
258  G4ThreeVector neutronVelocity = 1./G4Neutron::Neutron()->GetPDGMass()*theNeutron.GetMomentum();
259  G4double neutronVMag = neutronVelocity.mag();
260 
261  while(counter == 0 || std::abs(buffer-result/std::max(1,counter)) > 0.03*buffer) // Loop checking, 11.05.2015, T. Koi
262  {
263  if(counter) buffer = result/counter;
264  while (counter<size) // Loop checking, 11.05.2015, T. Koi
265  {
266  counter ++;
267  G4ReactionProduct aThermalNuc = aNuc.GetThermalNucleus(eleMass, aT);
268  boosted.Lorentz(theNeutron, aThermalNuc);
269  G4double theEkin = boosted.GetKineticEnergy();
270  aXsection = (*((*theCrossSections)(index))).GetValue(theEkin, outOfRange);
271  // velocity correction.
272  G4ThreeVector targetVelocity = 1./aThermalNuc.GetMass()*aThermalNuc.GetMomentum();
273  aXsection *= (targetVelocity-neutronVelocity).mag()/neutronVMag;
274  result += aXsection;
275  }
276  size += size;
277  }
278  result /= counter;
279 /*
280  // Checking impact of G4NEUTRONHP_NEGLECT_DOPPLER
281  G4cout << " result " << result << " "
282  << (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) << " "
283  << (*((*theCrossSections)(index))).GetValue(eKinetic, outOfRange) /result << G4endl;
284 */
285  return result;
286 }
G4double G4ParticleHPJENDLHEData::G4double result
static G4double GetNuclearMass(const G4double A, const G4double Z)
void Lorentz(const G4ReactionProduct &p1, const G4ReactionProduct &p2)
G4double GetKineticEnergy() const
G4double GetN() const
Definition: G4Element.hh:135
void SetMomentum(const G4double x, const G4double y, const G4double z)
G4double GetZ() const
Definition: G4Element.hh:131
#define buffer
Definition: xmlparse.cc:628
static const G4double eps
G4ParticleDefinition * GetDefinition() const
G4ReactionProduct GetThermalNucleus(G4double aMass, G4double temp=-1) const
Definition: G4Nucleus.cc:143
int G4int
Definition: G4Types.hh:78
float k_Boltzmann
Definition: hepunit.py:299
bool G4bool
Definition: G4Types.hh:79
size_t GetIndex() const
Definition: G4Element.hh:182
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static constexpr double kelvin
Definition: G4SIunits.hh:281
G4double GetKineticEnergy() const
G4double GetPDGMass() const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4ThreeVector GetMomentum() const
double G4double
Definition: G4Types.hh:76
double mag() const
G4double GetMass() const
G4ThreeVector GetMomentum() const

Here is the call graph for this function:

Here is the caller graph for this function:

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

Reimplemented from G4VCrossSectionDataSet.

Definition at line 91 of file G4ParticleHPElasticData.cc.

96 {
97  if ( dp->GetKineticEnergy() == ke_cache && element == element_cache && material == material_cache ) return xs_cache;
98 
99  ke_cache = dp->GetKineticEnergy();
100  element_cache = element;
101  material_cache = material;
102  G4double xs = GetCrossSection( dp , element , material->GetTemperature() );
103  xs_cache = xs;
104  return xs;
105 }
G4double GetKineticEnergy() const
G4double GetCrossSection(const G4DynamicParticle *, const G4Element *, G4double aT)
string material
Definition: eplot.py:19
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4int G4ParticleHPElasticData::GetVerboseLevel ( ) const
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 289 of file G4ParticleHPElasticData.cc.

290 {
292 }
static G4ParticleHPManager * GetInstance()

Here is the call graph for this function:

void G4ParticleHPElasticData::IgnoreOnFlightDopplerBroadening ( )
inline

Definition at line 83 of file G4ParticleHPElasticData.hh.

83 { onFlightDB = false; };
G4bool G4ParticleHPElasticData::IsIsoApplicable ( const G4DynamicParticle dp,
G4int  ,
G4int  ,
const G4Element ,
const G4Material  
)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 77 of file G4ParticleHPElasticData.cc.

81 {
82 
83  G4double eKin = dp->GetKineticEnergy();
84  if ( eKin > GetMaxKinEnergy()
85  || eKin < GetMinKinEnergy()
86  || dp->GetDefinition() != G4Neutron::Neutron() ) return false;
87 
88  return true;
89 }
G4double GetKineticEnergy() const
G4ParticleDefinition * GetDefinition() const
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

void G4ParticleHPElasticData::SetVerboseLevel ( G4int  newValue)
virtual

Reimplemented from G4VCrossSectionDataSet.

Definition at line 295 of file G4ParticleHPElasticData.cc.

296 {
298 }
static G4ParticleHPManager * GetInstance()

Here is the call graph for this function:


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