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

#include <G4NeutronHPThermalScattering.hh>

+ Inheritance diagram for G4NeutronHPThermalScattering:

Public Member Functions

 G4NeutronHPThermalScattering ()
 
 ~G4NeutronHPThermalScattering ()
 
G4HadFinalState * ApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &aTargetNucleus)
 
virtual const std::pair< G4double, G4double > GetFatalEnergyCheckLevels () const
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4HadFinalState * ApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &, G4Nucleus &)
 
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)
 
const G4HadronicInteraction * GetMyPointer () const
 
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
 
G4bool operator== (const G4HadronicInteraction &right) const
 
G4bool operator!= (const G4HadronicInteraction &right) const
 
virtual const std::pair< G4double, G4double > GetFatalEnergyCheckLevels () const
 
virtual std::pair< G4double, G4double > GetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void ModelDescription (std::ostream &outFile) const
 

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 64 of file G4NeutronHPThermalScattering.hh.

Constructor & Destructor Documentation

◆ G4NeutronHPThermalScattering()

G4NeutronHPThermalScattering::G4NeutronHPThermalScattering ( )

Definition at line 49 of file G4NeutronHPThermalScattering.cc.

50 :G4HadronicInteraction("NeutronHPThermalScattering")
51{
52
53 theHPElastic = new G4NeutronHPElastic();
54
55 SetMinEnergy( 0.*eV );
56 SetMaxEnergy( 4*eV );
57 theXSection = new G4NeutronHPThermalScatteringData();
58 theXSection->BuildPhysicsTable( *(G4Neutron::Neutron()) );
59
60 buildPhysicsTable();
61
62}
G4HadronicInteraction::SetMinEnergy
void SetMinEnergy(G4double anEnergy)
Definition: G4HadronicInteraction.hh:90
G4HadronicInteraction::SetMaxEnergy
void SetMaxEnergy(const G4double anEnergy)
Definition: G4HadronicInteraction.hh:103
G4NeutronHPThermalScatteringData::BuildPhysicsTable
void BuildPhysicsTable(const G4ParticleDefinition &)
Definition: G4NeutronHPThermalScatteringData.cc:204
G4Neutron::Neutron
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104

◆ ~G4NeutronHPThermalScattering()

G4NeutronHPThermalScattering::~G4NeutronHPThermalScattering ( )

Definition at line 66 of file G4NeutronHPThermalScattering.cc.

67{
68
69 for ( std::map < G4int , std::map < G4double , std::vector < E_isoAng* >* >* >::iterator it = incoherentFSs.begin() ; it != incoherentFSs.end() ; it++ )
70 {
71 std::map < G4double , std::vector < E_isoAng* >* >::iterator itt;
72 for ( itt = it->second->begin() ; itt != it->second->end() ; itt++ )
73 {
74 std::vector< E_isoAng* >::iterator ittt;
75 for ( ittt = itt->second->begin(); ittt != itt->second->end() ; ittt++ )
76 {
77 delete *ittt;
78 }
79 delete itt->second;
80 }
81 delete it->second;
82 }
83
84 for ( std::map < G4int , std::map < G4double , std::vector < std::pair< G4double , G4double >* >* >* >::iterator it = coherentFSs.begin() ; it != coherentFSs.end() ; it++ )
85 {
86 std::map < G4double , std::vector < std::pair< G4double , G4double >* >* >::iterator itt;
87 for ( itt = it->second->begin() ; itt != it->second->end() ; itt++ )
88 {
89 std::vector < std::pair< G4double , G4double >* >::iterator ittt;
90 for ( ittt = itt->second->begin(); ittt != itt->second->end() ; ittt++ )
91 {
92 delete *ittt;
93 }
94 delete itt->second;
95 }
96 delete it->second;
97 }
98
99 for ( std::map < G4int , std::map < G4double , std::vector < E_P_E_isoAng* >* >* >::iterator it = inelasticFSs.begin() ; it != inelasticFSs.end() ; it++ )
100 {
101 std::map < G4double , std::vector < E_P_E_isoAng* >* >::iterator itt;
102 for ( itt = it->second->begin() ; itt != it->second->end() ; itt++ )
103 {
104 std::vector < E_P_E_isoAng* >::iterator ittt;
105 for ( ittt = itt->second->begin(); ittt != itt->second->end() ; ittt++ )
106 {
107 std::vector < E_isoAng* >::iterator it4;
108 for ( it4 = (*ittt)->vE_isoAngle.begin() ; it4 != (*ittt)->vE_isoAngle.end() ; it4++ )
109 {
110 delete *it4;
111 }
112 delete *ittt;
113 }
114 delete itt->second;
115 }
116 delete it->second;
117 }
118
119 delete theHPElastic;
120 delete theXSection;
121}
G4double
double G4double
Definition: G4Types.hh:64
G4int
int G4int
Definition: G4Types.hh:66

Member Function Documentation

◆ ApplyYourself()

G4HadFinalState * G4NeutronHPThermalScattering::ApplyYourself ( const G4HadProjectile &  aTrack,
G4Nucleus &  aTargetNucleus 
)
virtual

Implements G4HadronicInteraction.

Definition at line 302 of file G4NeutronHPThermalScattering.cc.

303{
304
305// Select Element > Reaction >
306
307 const G4Material * theMaterial = aTrack.GetMaterial();
308 G4double aTemp = theMaterial->GetTemperature();
309 G4int n = theMaterial->GetNumberOfElements();
310 //static const G4ElementTable* theElementTable = G4Element::GetElementTable();
311
312 G4bool findThermalElement = false;
313 G4int ielement;
314 const G4Element* theElement = NULL;
315 for ( G4int i = 0; i < n ; i++ )
316 {
317 theElement = theMaterial->GetElement(i);
318 //Select target element
319 if ( aNucleus.GetZ_asInt() == (G4int)(theElement->GetZ() + 0.5 ) )
320 {
321 //Check Applicability of Thermal Scattering
322 if ( getTS_ID( NULL , theElement ) != -1 )
323 {
324 ielement = getTS_ID( NULL , theElement );
325 findThermalElement = true;
326 break;
327 }
328 else if ( getTS_ID( theMaterial , theElement ) != -1 )
329 {
330 ielement = getTS_ID( theMaterial , theElement );
331 findThermalElement = true;
332 break;
333 }
334 }
335 }
336
337 if ( findThermalElement == true )
338 {
339
340// Select Reaction (Inelastic, coherent, incoherent)
341
342 G4ParticleDefinition* pd = const_cast< G4ParticleDefinition* >( aTrack.GetDefinition() );
343 G4DynamicParticle* dp = new G4DynamicParticle ( pd , aTrack.Get4Momentum() );
344 G4double total = theXSection->GetCrossSection( dp , theElement , theMaterial );
345 G4double inelastic = theXSection->GetInelasticCrossSection( dp , theElement , theMaterial );
346
347
348 G4double random = G4UniformRand();
349 if ( random <= inelastic/total )
350 {
351 // Inelastic
352
353 // T_L and T_H
354 std::map < G4double , std::vector< E_P_E_isoAng* >* >::iterator it;
355 std::vector<G4double> v_temp;
356 v_temp.clear();
357 for ( it = inelasticFSs.find( ielement )->second->begin() ; it != inelasticFSs.find( ielement )->second->end() ; it++ )
358 {
359 v_temp.push_back( it->first );
360 }
361
362// T_L T_H
363 std::pair < G4double , G4double > tempLH = find_LH ( aTemp , &v_temp );
364//
365// For T_L aNEP_EPM_TL and T_H aNEP_EPM_TH
366//
367 std::vector< E_P_E_isoAng* >* vNEP_EPM_TL = 0;
368 std::vector< E_P_E_isoAng* >* vNEP_EPM_TH = 0;
369
370 if ( tempLH.first != 0.0 && tempLH.second != 0.0 )
371 {
372 vNEP_EPM_TL = inelasticFSs.find( ielement )->second->find ( tempLH.first/kelvin )->second;
373 vNEP_EPM_TH = inelasticFSs.find( ielement )->second->find ( tempLH.second/kelvin )->second;
374 }
375 else if ( tempLH.first == 0.0 )
376 {
377 std::map < G4double , std::vector< E_P_E_isoAng* >* >::iterator itm;
378 itm = inelasticFSs.find( ielement )->second->begin();
379 vNEP_EPM_TL = itm->second;
380 itm++;
381 vNEP_EPM_TH = itm->second;
382 }
383 else if ( tempLH.second == 0.0 )
384 {
385 std::map < G4double , std::vector< E_P_E_isoAng* >* >::iterator itm;
386 itm = inelasticFSs.find( ielement )->second->end();
387 itm--;
388 vNEP_EPM_TH = itm->second;
389 itm--;
390 vNEP_EPM_TL = itm->second;
391 }
392
393//
394
395 G4double rand_for_sE = G4UniformRand();
396
397 std::pair< G4double , E_isoAng > TL = create_sE_and_EPM_from_pE_and_vE_P_E_isoAng ( rand_for_sE , aTrack.GetKineticEnergy() , vNEP_EPM_TL );
398 std::pair< G4double , E_isoAng > TH = create_sE_and_EPM_from_pE_and_vE_P_E_isoAng ( rand_for_sE , aTrack.GetKineticEnergy() , vNEP_EPM_TH );
399
400 G4double sE;
401 sE = get_linear_interpolated ( aTemp , std::pair < G4double , G4double > ( tempLH.first , TL.first ) , std::pair < G4double , G4double > ( tempLH.second , TH.first ) );
402 E_isoAng anE_isoAng;
403 if ( TL.second.n == TH.second.n )
404 {
405 anE_isoAng.energy = sE;
406 anE_isoAng.n = TL.second.n;
407 for ( G4int i=0 ; i < anE_isoAng.n ; i++ )
408 {
409 G4double angle;
410 angle = get_linear_interpolated ( aTemp , std::pair< G4double , G4double > ( tempLH.first , TL.second.isoAngle[ i ] ) , std::pair< G4double , G4double > ( tempLH.second , TH.second.isoAngle[ i ] ) );
411 anE_isoAng.isoAngle.push_back( angle );
412 }
413 }
414 else
415 {
416 std::cout << "Do not Suuport yet." << std::endl;
417 }
418
419 //set
420 theParticleChange.SetEnergyChange( sE );
421 G4double mu = getMu( &anE_isoAng );
422 theParticleChange.SetMomentumChange( 0.0 , std::sqrt ( 1 - mu*mu ) , mu );
423
424 }
425 //else if ( random <= ( inelastic + theXSection->GetCoherentCrossSection( dp , (*theElementTable)[ ielement ] , aTemp ) ) / total )
426 else if ( random <= ( inelastic + theXSection->GetCoherentCrossSection( dp , theElement , theMaterial ) ) / total )
427 {
428 // Coherent Elastic
429
430 G4double E = aTrack.GetKineticEnergy();
431
432 // T_L and T_H
433 std::map < G4double , std::vector< std::pair< G4double , G4double >* >* >::iterator it;
434 std::vector<G4double> v_temp;
435 v_temp.clear();
436 for ( it = coherentFSs.find( ielement )->second->begin() ; it != coherentFSs.find( ielement )->second->end() ; it++ )
437 {
438 v_temp.push_back( it->first );
439 }
440
441// T_L T_H
442 std::pair < G4double , G4double > tempLH = find_LH ( aTemp , &v_temp );
443//
444//
445// For T_L anEPM_TL and T_H anEPM_TH
446//
447 std::vector< std::pair< G4double , G4double >* >* pvE_p_TL = 0;
448 std::vector< std::pair< G4double , G4double >* >* pvE_p_TH = 0;
449
450 if ( tempLH.first != 0.0 && tempLH.second != 0.0 )
451 {
452 pvE_p_TL = coherentFSs.find( ielement )->second->find ( tempLH.first/kelvin )->second;
453 pvE_p_TH = coherentFSs.find( ielement )->second->find ( tempLH.first/kelvin )->second;
454 }
455 else if ( tempLH.first == 0.0 )
456 {
457 pvE_p_TL = coherentFSs.find( ielement )->second->find ( v_temp[ 0 ] )->second;
458 pvE_p_TH = coherentFSs.find( ielement )->second->find ( v_temp[ 1 ] )->second;
459 }
460 else if ( tempLH.second == 0.0 )
461 {
462 pvE_p_TL = coherentFSs.find( ielement )->second->find ( v_temp.back() )->second;
463 std::vector< G4double >::iterator itv;
464 itv = v_temp.end();
465 itv--;
466 itv--;
467 pvE_p_TL = coherentFSs.find( ielement )->second->find ( *itv )->second;
468 }
469
470
471 std::vector< G4double > vE_T;
472 std::vector< G4double > vp_T;
473
474 G4int n1 = pvE_p_TL->size();
475 //G4int n2 = pvE_p_TH->size();
476
477 for ( G4int i=1 ; i < n1 ; i++ )
478 {
479 if ( (*pvE_p_TL)[i]->first != (*pvE_p_TH)[i]->first ) abort();
480 vE_T.push_back ( (*pvE_p_TL)[i]->first );
481 vp_T.push_back ( get_linear_interpolated ( aTemp , std::pair< G4double , G4double > ( tempLH.first , (*pvE_p_TL)[i]->second ) , std::pair< G4double , G4double > ( tempLH.second , (*pvE_p_TL)[i]->second ) ) );
482 }
483
484 G4int j = 0;
485 for ( G4int i = 1 ; i < n ; i++ )
486 {
487 if ( E/eV < vE_T[ i ] )
488 {
489 j = i-1;
490 break;
491 }
492 }
493
494 G4double rand_for_mu = G4UniformRand();
495
496 G4int k = 0;
497 for ( G4int i = 1 ; i < j ; i++ )
498 {
499 G4double Pi = vp_T[ i ] / vp_T[ j ];
500 if ( rand_for_mu < Pi )
501 {
502 k = i-1;
503 break;
504 }
505 }
506
507 //G4double Ei = vE_T[ j ];
508 G4double Ei = vE_T[ k ];
509
510 G4double mu = 1 - 2 * Ei / (E/eV) ;
511 //111102
512 if ( mu < -1.0 ) mu = -1.0;
513
514 theParticleChange.SetEnergyChange( E );
515 theParticleChange.SetMomentumChange( 0.0 , std::sqrt ( 1 - mu*mu ) , mu );
516
517
518 }
519 else
520 {
521 // InCoherent Elastic
522
523 // T_L and T_H
524 std::map < G4double , std::vector < E_isoAng* >* >::iterator it;
525 std::vector<G4double> v_temp;
526 v_temp.clear();
527 for ( it = incoherentFSs.find( ielement )->second->begin() ; it != incoherentFSs.find( ielement )->second->end() ; it++ )
528 {
529 v_temp.push_back( it->first );
530 }
531
532// T_L T_H
533 std::pair < G4double , G4double > tempLH = find_LH ( aTemp , &v_temp );
534
535//
536// For T_L anEPM_TL and T_H anEPM_TH
537//
538
539 E_isoAng anEPM_TL_E;
540 E_isoAng anEPM_TH_E;
541
542 if ( tempLH.first != 0.0 && tempLH.second != 0.0 )
543 {
544 anEPM_TL_E = create_E_isoAng_from_energy ( aTrack.GetKineticEnergy() , incoherentFSs.find( ielement )->second->find ( tempLH.first/kelvin )->second );
545 anEPM_TH_E = create_E_isoAng_from_energy ( aTrack.GetKineticEnergy() , incoherentFSs.find( ielement )->second->find ( tempLH.second/kelvin )->second );
546 }
547 else if ( tempLH.first == 0.0 )
548 {
549 anEPM_TL_E = create_E_isoAng_from_energy ( aTrack.GetKineticEnergy() , incoherentFSs.find( ielement )->second->find ( v_temp[ 0 ] )->second );
550 anEPM_TH_E = create_E_isoAng_from_energy ( aTrack.GetKineticEnergy() , incoherentFSs.find( ielement )->second->find ( v_temp[ 1 ] )->second );
551 }
552 else if ( tempLH.second == 0.0 )
553 {
554 anEPM_TH_E = create_E_isoAng_from_energy ( aTrack.GetKineticEnergy() , incoherentFSs.find( ielement )->second->find ( v_temp.back() )->second );
555 std::vector< G4double >::iterator itv;
556 itv = v_temp.end();
557 itv--;
558 itv--;
559 anEPM_TL_E = create_E_isoAng_from_energy ( aTrack.GetKineticEnergy() , incoherentFSs.find( ielement )->second->find ( *itv )->second );
560 }
561
562 // E_isoAng for aTemp and aTrack.GetKineticEnergy()
563 E_isoAng anEPM_T_E;
564
565 if ( anEPM_TL_E.n == anEPM_TH_E.n )
566 {
567 anEPM_T_E.n = anEPM_TL_E.n;
568 for ( G4int i=0 ; i < anEPM_TL_E.n ; i++ )
569 {
570 G4double angle;
571 angle = get_linear_interpolated ( aTemp , std::pair< G4double , G4double > ( tempLH.first , anEPM_TL_E.isoAngle[ i ] ) , std::pair< G4double , G4double > ( tempLH.second , anEPM_TH_E.isoAngle[ i ] ) );
572 anEPM_T_E.isoAngle.push_back( angle );
573 }
574 }
575 else
576 {
577 std::cout << "Do not Suuport yet." << std::endl;
578 }
579
580 // Decide mu
581 G4double mu = getMu ( &anEPM_T_E );
582
583 // Set Final State
584 theParticleChange.SetEnergyChange( aTrack.GetKineticEnergy() ); // No energy change in Elastic
585 theParticleChange.SetMomentumChange( 0.0 , std::sqrt ( 1 - mu*mu ) , mu );
586
587 }
588 delete dp;
589
590 return &theParticleChange;
591
592 }
593 else
594 {
595 // Not thermal element
596 // Neutron HP will handle
597 return theHPElastic -> ApplyYourself( aTrack, aNucleus );
598 }
599
600}
G4bool
bool G4bool
Definition: G4Types.hh:67
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:53
G4Element::GetZ
G4double GetZ() const
Definition: G4Element.hh:131
G4HadFinalState::SetEnergyChange
void SetEnergyChange(G4double anEnergy)
Definition: G4HadFinalState.cc:42
G4HadFinalState::SetMomentumChange
void SetMomentumChange(const G4ThreeVector &aV)
Definition: G4HadFinalState.cc:54
G4HadProjectile::GetMaterial
const G4Material * GetMaterial() const
Definition: G4HadProjectile.hh:76
G4HadProjectile::GetDefinition
const G4ParticleDefinition * GetDefinition() const
Definition: G4HadProjectile.hh:81
G4HadProjectile::GetKineticEnergy
G4double GetKineticEnergy() const
Definition: G4HadProjectile.hh:106
G4HadProjectile::Get4Momentum
const G4LorentzVector & Get4Momentum() const
Definition: G4HadProjectile.hh:86
G4Material::GetTemperature
G4double GetTemperature() const
Definition: G4Material.hh:181
G4Material::GetElement
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:201
G4Material::GetNumberOfElements
size_t GetNumberOfElements() const
Definition: G4Material.hh:185
G4NeutronHPThermalScatteringData::GetInelasticCrossSection
G4double GetInelasticCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
Definition: G4NeutronHPThermalScatteringData.cc:441
G4NeutronHPThermalScatteringData::GetCrossSection
G4double GetCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
Definition: G4NeutronHPThermalScatteringData.cc:419
G4NeutronHPThermalScatteringData::GetCoherentCrossSection
G4double GetCoherentCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
Definition: G4NeutronHPThermalScatteringData.cc:450
G4NeutronHPThermalScattering::ApplyYourself
G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &aTargetNucleus)
Definition: G4NeutronHPThermalScattering.cc:302
E_isoAng
Definition: G4NeutronHPThermalScattering.hh:49
E_isoAng::isoAngle
std::vector< G4double > isoAngle
Definition: G4NeutronHPThermalScattering.hh:52
E_isoAng::energy
G4double energy
Definition: G4NeutronHPThermalScattering.hh:50
E_isoAng::n
G4int n
Definition: G4NeutronHPThermalScattering.hh:51

Referenced by ApplyYourself().

◆ GetFatalEnergyCheckLevels()

const std::pair< G4double, G4double > G4NeutronHPThermalScattering::GetFatalEnergyCheckLevels ( ) const
virtual

Reimplemented from G4HadronicInteraction.

Definition at line 999 of file G4NeutronHPThermalScattering.cc.

1000{
1001 //return std::pair<G4double, G4double>(10*perCent,10*GeV);
1002 return std::pair<G4double, G4double>(10*perCent,DBL_MAX);
1003}
DBL_MAX
#define DBL_MAX
Definition: templates.hh:83
The documentation for this class was generated from the following files: