Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
Public Member Functions | List of all members
G4ParticleHPThermalScattering Class Reference

#include <G4ParticleHPThermalScattering.hh>

+ Inheritance diagram for G4ParticleHPThermalScattering:

Public Member Functions

 G4ParticleHPThermalScattering ()
 
 ~G4ParticleHPThermalScattering () override
 
G4HadFinalState * ApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &aTargetNucleus) override
 
const std::pair< G4double, G4double > GetFatalEnergyCheckLevels () const override
 
void AddUserThermalScatteringFile (G4String, G4String)
 
void BuildPhysicsTable (const G4ParticleDefinition &) override
 
void ModelDescription (std::ostream &outFile) const override
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4int GetVerboseLevel () const
 
void SetVerboseLevel (G4int value)
 
const G4String & GetModelName () const
 
void DeActivateFor (const G4Material *aMaterial)
 
void ActivateFor (const G4Material *aMaterial)
 
void DeActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Element *anElement)
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4bool IsBlocked (const G4Element *anElement) const
 
void SetRecoilEnergyThreshold (G4double val)
 
G4double GetRecoilEnergyThreshold () const
 
virtual std::pair< G4double, G4double > GetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void InitialiseModel ()
 
 G4HadronicInteraction (const G4HadronicInteraction &right)=delete
 
const G4HadronicInteraction & operator= (const G4HadronicInteraction &right)=delete
 
G4bool operator== (const G4HadronicInteraction &right) const =delete
 
G4bool operator!= (const G4HadronicInteraction &right) const =delete
 

Additional Inherited Members

- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 84 of file G4ParticleHPThermalScattering.hh.

Constructor & Destructor Documentation

◆ G4ParticleHPThermalScattering()

G4ParticleHPThermalScattering::G4ParticleHPThermalScattering ( )

Definition at line 56 of file G4ParticleHPThermalScattering.cc.

57 : G4HadronicInteraction("NeutronHPThermalScattering")
58{
59 theHPElastic = new G4ParticleHPElastic();
60
61 SetMinEnergy(0. * eV);
62 SetMaxEnergy(4 * eV);
63 theXSection = new G4ParticleHPThermalScatteringData();
64
65 nMaterial = 0;
66 nElement = 0;
67}
G4HadronicInteraction::SetMinEnergy
void SetMinEnergy(G4double anEnergy)
Definition G4HadronicInteraction.hh:89
G4HadronicInteraction::G4HadronicInteraction
G4HadronicInteraction(const G4String &modelName="HadronicModel")
Definition G4HadronicInteraction.cc:41
G4HadronicInteraction::SetMaxEnergy
void SetMaxEnergy(const G4double anEnergy)
Definition G4HadronicInteraction.hh:102

◆ ~G4ParticleHPThermalScattering()

G4ParticleHPThermalScattering::~G4ParticleHPThermalScattering ( )
override

Definition at line 69 of file G4ParticleHPThermalScattering.cc.

70{
71 delete theHPElastic;
72}

Member Function Documentation

◆ AddUserThermalScatteringFile()

void G4ParticleHPThermalScattering::AddUserThermalScatteringFile ( G4String nameG4Element,
G4String filename )

Definition at line 1175 of file G4ParticleHPThermalScattering.cc.

1177{
1178 names.AddThermalElement(nameG4Element, filename);
1179 theXSection->AddUserThermalScatteringFile(nameG4Element, filename);
1180 buildPhysicsTable();
1181}

◆ ApplyYourself()

G4HadFinalState * G4ParticleHPThermalScattering::ApplyYourself ( const G4HadProjectile & aTrack,
G4Nucleus & aTargetNucleus )
overridevirtual

Reimplemented from G4HadronicInteraction.

Definition at line 304 of file G4ParticleHPThermalScattering.cc.

306{
307 // Select Element > Reaction >
308
309 const G4Material* theMaterial = aTrack.GetMaterial();
310 G4double aTemp = theMaterial->GetTemperature();
311 auto n = (G4int)theMaterial->GetNumberOfElements();
312
313 G4bool findThermalElement = false;
314 G4int ielement;
315 const G4Element* theElement = nullptr;
316 for (G4int i = 0; i < n; ++i) {
317 theElement = theMaterial->GetElement(i);
318 // Select target element
319 if (aNucleus.GetZ_asInt() == (G4int)(theElement->GetZ() + 0.5)) {
320 // Check Applicability of Thermal Scattering
321 if (getTS_ID(nullptr, theElement) != -1) {
322 ielement = getTS_ID(nullptr, theElement);
323 findThermalElement = true;
324 break;
325 }
326 if (getTS_ID(theMaterial, theElement) != -1) {
327 ielement = getTS_ID(theMaterial, theElement);
328 findThermalElement = true;
329 break;
330 }
331 }
332 }
333
334 if (findThermalElement) {
335 // Select Reaction (Inelastic, coherent, incoherent)
336 const G4ParticleDefinition* pd = aTrack.GetDefinition();
337 auto dp = new G4DynamicParticle(pd, aTrack.Get4Momentum());
338 G4double total = theXSection->GetCrossSection(dp, theElement, theMaterial);
339 G4double inelastic = theXSection->GetInelasticCrossSection(dp, theElement, theMaterial);
340
341 G4double random = G4UniformRand();
342 if (random <= inelastic / total) {
343 // Inelastic
344
345 std::vector<G4double> v_temp;
346 v_temp.clear();
347 for (auto it = inelasticFSs->find(ielement)->second->cbegin();
348 it != inelasticFSs->find(ielement)->second->cend(); ++it)
349 {
350 v_temp.push_back(it->first);
351 }
352
353 std::pair<G4double, G4double> tempLH = find_LH(aTemp, &v_temp);
354 //
355 // For T_L aNEP_EPM_TL and T_H aNEP_EPM_TH
356 //
357 std::vector<E_P_E_isoAng*>* vNEP_EPM_TL = nullptr;
358 std::vector<E_P_E_isoAng*>* vNEP_EPM_TH = nullptr;
359
360 if (tempLH.first != 0.0 && tempLH.second != 0.0) {
361 vNEP_EPM_TL = inelasticFSs->find(ielement)->second->find(tempLH.first / kelvin)->second;
362 vNEP_EPM_TH = inelasticFSs->find(ielement)->second->find(tempLH.second / kelvin)->second;
363 }
364 else if (tempLH.first == 0.0) {
365 auto itm = inelasticFSs->find(ielement)->second->cbegin();
366 vNEP_EPM_TL = itm->second;
367 ++itm;
368 vNEP_EPM_TH = itm->second;
369 tempLH.first = tempLH.second;
370 tempLH.second = itm->first;
371 }
372 else if (tempLH.second == 0.0) {
373 auto itm = inelasticFSs->find(ielement)->second->cend();
374 --itm;
375 vNEP_EPM_TH = itm->second;
376 --itm;
377 vNEP_EPM_TL = itm->second;
378 tempLH.second = tempLH.first;
379 tempLH.first = itm->first;
380 }
381
382 G4double sE = 0., mu = 1.0;
383
384 // New Geant4 method - Stochastic temperature interpolation of the final state
385 // (continuous temperature interpolation was used previously)
386 std::pair<G4double, G4double> secondaryParam;
387 G4double rand_temp = G4UniformRand();
388 if (rand_temp < (aTemp - tempLH.first) / (tempLH.second - tempLH.first))
389 secondaryParam = sample_inelastic_E_mu(aTrack.GetKineticEnergy(), vNEP_EPM_TH);
390 else
391 secondaryParam = sample_inelastic_E_mu(aTrack.GetKineticEnergy(), vNEP_EPM_TL);
392
393 sE = secondaryParam.first;
394 mu = secondaryParam.second;
395
396 // set
397 theParticleChange.SetEnergyChange(sE);
398 G4double phi = CLHEP::twopi * G4UniformRand();
399 G4double sint = std::sqrt(1 - mu * mu);
400 theParticleChange.SetMomentumChange(sint * std::cos(phi), sint * std::sin(phi), mu);
401 }
402 else if (random
403 <= (inelastic + theXSection->GetCoherentCrossSection(dp, theElement, theMaterial))
404 / total)
405 {
406 // Coherent Elastic
407
408 G4double E = aTrack.GetKineticEnergy();
409
410 // T_L and T_H
411 std::vector<G4double> v_temp;
412 v_temp.clear();
413 for (auto it = coherentFSs->find(ielement)->second->cbegin();
414 it != coherentFSs->find(ielement)->second->cend(); ++it)
415 {
416 v_temp.push_back(it->first);
417 }
418
419 // T_L T_H
420 std::pair<G4double, G4double> tempLH = find_LH(aTemp, &v_temp);
421 //
422 //
423 // For T_L anEPM_TL and T_H anEPM_TH
424 //
425 std::vector<std::pair<G4double, G4double>*>* pvE_p_TL = nullptr;
426 std::vector<std::pair<G4double, G4double>*>* pvE_p_TH = nullptr;
427
428 if (tempLH.first != 0.0 && tempLH.second != 0.0) {
429 pvE_p_TL = coherentFSs->find(ielement)->second->find(tempLH.first / kelvin)->second;
430 pvE_p_TH = coherentFSs->find(ielement)->second->find(tempLH.first / kelvin)->second;
431 }
432 else if (tempLH.first == 0.0) {
433 pvE_p_TL = coherentFSs->find(ielement)->second->find(v_temp[0])->second;
434 pvE_p_TH = coherentFSs->find(ielement)->second->find(v_temp[1])->second;
435 tempLH.first = tempLH.second;
436 tempLH.second = v_temp[1];
437 }
438 else if (tempLH.second == 0.0) {
439 pvE_p_TH = coherentFSs->find(ielement)->second->find(v_temp.back())->second;
440 auto itv = v_temp.cend();
441 --itv;
442 --itv;
443 pvE_p_TL = coherentFSs->find(ielement)->second->find(*itv)->second;
444 tempLH.second = tempLH.first;
445 tempLH.first = *itv;
446 }
447 else {
448 // tempLH.first == 0.0 && tempLH.second
449 throw G4HadronicException(
450 __FILE__, __LINE__,
451 "A problem is found in Thermal Scattering Data! Unexpected temperature values in data");
452 }
453
454 std::vector<G4double> vE_T;
455 std::vector<G4double> vp_T;
456
457 auto n1 = (G4int)pvE_p_TL->size();
458
459 // New Geant4 method - Stochastic interpolation of the final state
460 std::vector<std::pair<G4double, G4double>*>* pvE_p_T_sampled;
461 G4double rand_temp = G4UniformRand();
462 if (rand_temp < (aTemp - tempLH.first) / (tempLH.second - tempLH.first))
463 pvE_p_T_sampled = pvE_p_TH;
464 else
465 pvE_p_T_sampled = pvE_p_TL;
466
467 // 171005 fix bug, contribution from H.N. TRAN@CEA
468 for (G4int i = 0; i < n1; ++i) {
469 vE_T.push_back((*pvE_p_T_sampled)[i]->first);
470 vp_T.push_back((*pvE_p_T_sampled)[i]->second);
471 }
472
473 G4int j = 0;
474 for (G4int i = 1; i < n1; ++i) {
475 if (E / eV < vE_T[i]) {
476 j = i - 1;
477 break;
478 }
479 }
480
481 G4double rand_for_mu = G4UniformRand();
482
483 G4int k = 0;
484 for (G4int i = 0; i <= j; ++i) {
485 G4double Pi = vp_T[i] / vp_T[j];
486 if (rand_for_mu < Pi) {
487 k = i;
488 break;
489 }
490 }
491
492 G4double Ei = vE_T[k];
493
494 G4double mu = 1 - 2 * Ei / (E / eV);
495
496 if (mu < -1.0) mu = -1.0;
497
498 theParticleChange.SetEnergyChange(E);
499 G4double phi = CLHEP::twopi * G4UniformRand();
500 G4double sint = std::sqrt(1 - mu * mu);
501 theParticleChange.SetMomentumChange(sint * std::cos(phi), sint * std::sin(phi), mu);
502 }
503 else {
504 // InCoherent Elastic
505
506 // T_L and T_H
507 std::vector<G4double> v_temp;
508 v_temp.clear();
509 for (auto it = incoherentFSs->find(ielement)->second->cbegin();
510 it != incoherentFSs->find(ielement)->second->cend(); ++it)
511 {
512 v_temp.push_back(it->first);
513 }
514
515 // T_L T_H
516 std::pair<G4double, G4double> tempLH = find_LH(aTemp, &v_temp);
517
518 //
519 // For T_L anEPM_TL and T_H anEPM_TH
520 //
521
522 E_isoAng anEPM_TL_E;
523 E_isoAng anEPM_TH_E;
524
525 if (tempLH.first != 0.0 && tempLH.second != 0.0) {
526 // Interpolate TL and TH
527 anEPM_TL_E = create_E_isoAng_from_energy(
528 aTrack.GetKineticEnergy(),
529 incoherentFSs->find(ielement)->second->find(tempLH.first / kelvin)->second);
530 anEPM_TH_E = create_E_isoAng_from_energy(
531 aTrack.GetKineticEnergy(),
532 incoherentFSs->find(ielement)->second->find(tempLH.second / kelvin)->second);
533 }
534 else if (tempLH.first == 0.0) {
535 // Extrapolate T0 and T1
536 anEPM_TL_E = create_E_isoAng_from_energy(
537 aTrack.GetKineticEnergy(),
538 incoherentFSs->find(ielement)->second->find(v_temp[0])->second);
539 anEPM_TH_E = create_E_isoAng_from_energy(
540 aTrack.GetKineticEnergy(),
541 incoherentFSs->find(ielement)->second->find(v_temp[1])->second);
542 tempLH.first = tempLH.second;
543 tempLH.second = v_temp[1];
544 }
545 else if (tempLH.second == 0.0) {
546 // Extrapolate Tmax-1 and Tmax
547 anEPM_TH_E = create_E_isoAng_from_energy(
548 aTrack.GetKineticEnergy(),
549 incoherentFSs->find(ielement)->second->find(v_temp.back())->second);
550 auto itv = v_temp.cend();
551 --itv;
552 --itv;
553 anEPM_TL_E = create_E_isoAng_from_energy(
554 aTrack.GetKineticEnergy(), incoherentFSs->find(ielement)->second->find(*itv)->second);
555 tempLH.second = tempLH.first;
556 tempLH.first = *itv;
557 }
558
559 // E_isoAng for aTemp and aTrack.GetKineticEnergy()
560 G4double mu = 1.0;
561
562 // New Geant4 method - Stochastic interpolation of the final state
563 E_isoAng anEPM_T_E_sampled;
564 G4double rand_temp = G4UniformRand();
565 if (rand_temp < (aTemp - tempLH.first) / (tempLH.second - tempLH.first))
566 anEPM_T_E_sampled = anEPM_TH_E;
567 else
568 anEPM_T_E_sampled = anEPM_TL_E;
569
570 mu = getMu(&anEPM_T_E_sampled);
571
572 // Set Final State
573 theParticleChange.SetEnergyChange(aTrack.GetKineticEnergy()); // No energy change in Elastic
574 G4double phi = CLHEP::twopi * G4UniformRand();
575 G4double sint = std::sqrt(1 - mu * mu);
576 theParticleChange.SetMomentumChange(sint * std::cos(phi), sint * std::sin(phi), mu);
577 }
578 delete dp;
579
580 return &theParticleChange;
581 }
582 // Not thermal element
583 // Neutron HP will handle
584 return theHPElastic->ApplyYourself(aTrack, aNucleus,
585 true); // L. Thulliez 2021/05/04 (CEA-Saclay)
586}
G4double
double G4double
Definition G4Types.hh:83
G4bool
bool G4bool
Definition G4Types.hh:86
G4int
int G4int
Definition G4Types.hh:85
G4UniformRand
#define G4UniformRand()
Definition Randomize.hh:52
G4Element::GetZ
G4double GetZ() const
Definition G4Element.hh:119
G4HadFinalState::SetEnergyChange
void SetEnergyChange(G4double anEnergy)
Definition G4HadFinalState.cc:39
G4HadFinalState::SetMomentumChange
void SetMomentumChange(const G4ThreeVector &aV)
Definition G4HadFinalState.hh:56
G4HadProjectile::GetMaterial
const G4Material * GetMaterial() const
Definition G4HadProjectile.hh:81
G4HadProjectile::GetDefinition
const G4ParticleDefinition * GetDefinition() const
Definition G4HadProjectile.hh:86
G4HadProjectile::GetKineticEnergy
G4double GetKineticEnergy() const
Definition G4HadProjectile.hh:116
G4HadProjectile::Get4Momentum
const G4LorentzVector & Get4Momentum() const
Definition G4HadProjectile.hh:91
G4Material::GetTemperature
G4double GetTemperature() const
Definition G4Material.hh:176
G4Material::GetElement
const G4Element * GetElement(G4int iel) const
Definition G4Material.hh:192
G4Material::GetNumberOfElements
std::size_t GetNumberOfElements() const
Definition G4Material.hh:180
G4ParticleHPElastic::ApplyYourself
G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &aTargetNucleus) override
Definition G4ParticleHPElastic.cc:63
G4ParticleHPThermalScatteringData::GetCrossSection
G4double GetCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
Definition G4ParticleHPThermalScatteringData.cc:341
G4ParticleHPThermalScatteringData::GetCoherentCrossSection
G4double GetCoherentCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
Definition G4ParticleHPThermalScatteringData.cc:373
G4ParticleHPThermalScatteringData::GetInelasticCrossSection
G4double GetInelasticCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
Definition G4ParticleHPThermalScatteringData.cc:362
G4INCL::CrossSections::total
G4double total(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:61
E_isoAng
Definition G4ParticleHPThermalScattering.hh:53

◆ BuildPhysicsTable()

void G4ParticleHPThermalScattering::BuildPhysicsTable ( const G4ParticleDefinition & particle)
overridevirtual

Reimplemented from G4HadronicInteraction.

Definition at line 121 of file G4ParticleHPThermalScattering.cc.

122{
123 buildPhysicsTable();
124 theHPElastic->BuildPhysicsTable(particle);
125}
G4ParticleHPElastic::BuildPhysicsTable
void BuildPhysicsTable(const G4ParticleDefinition &) override
Definition G4ParticleHPElastic.cc:158

◆ GetFatalEnergyCheckLevels()

const std::pair< G4double, G4double > G4ParticleHPThermalScattering::GetFatalEnergyCheckLevels ( ) const
overridevirtual

Reimplemented from G4HadronicInteraction.

Definition at line 1169 of file G4ParticleHPThermalScattering.cc.

1170{
1171 // max energy non-conservation is mass of heavy nucleus
1172 return std::pair<G4double, G4double>(10.0 * perCent, 350.0 * CLHEP::GeV);
1173}

◆ ModelDescription()

void G4ParticleHPThermalScattering::ModelDescription ( std::ostream & outFile) const
overridevirtual

Reimplemented from G4HadronicInteraction.

Definition at line 1197 of file G4ParticleHPThermalScattering.cc.

1198{
1199 outFile << "High Precision model based on thermal scattering data in\n"
1200 << "evaluated nuclear data libraries for neutrons below 5eV\n"
1201 << "on specific materials\n";
1202}
The documentation for this class was generated from the following files: