Geant4 Cross Reference |
1 // 2 // ******************************************************************** 3 // * License and Disclaimer * 4 // * * 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. * 10 // * * 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitation of liability. * 17 // * * 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************************************** 25 // 26 // neutron_hp -- source file 27 // J.P. Wellisch, Nov-1996 28 // A prototype of the low energy neutron transport model. 29 // 30 // 070523 bug fix for G4FPE_DEBUG on by A. Howard ( and T. Koi) 31 // 32 // P. Arce, June-2014 Conversion neutron_hp to particle_hp 33 // V. Ivanchenko July-2023 converted back 34 // 35 #include "G4NeutronHPCapture.hh" 36 37 #include "G4IonTable.hh" 38 #include "G4NeutronHPCaptureFS.hh" 39 #include "G4ParticleHPDeExGammas.hh" 40 #include "G4ParticleHPManager.hh" 41 #include "G4ParticleTable.hh" 42 #include "G4ParticleHPThermalBoost.hh" 43 #include "G4SystemOfUnits.hh" 44 #include "G4Threading.hh" 45 46 G4NeutronHPCapture::G4NeutronHPCapture() : G4HadronicInteraction("NeutronHPCapture") 47 { 48 SetMinEnergy(0.0); 49 SetMaxEnergy(20. * MeV); 50 } 51 52 G4NeutronHPCapture::~G4NeutronHPCapture() 53 { 54 if (!G4Threading::IsWorkerThread()) { 55 if (theCapture != nullptr) { 56 for (auto& ite : *theCapture) { 57 delete ite; 58 } 59 theCapture->clear(); 60 } 61 } 62 } 63 64 G4HadFinalState* 65 G4NeutronHPCapture::ApplyYourself(const G4HadProjectile& aTrack, 66 G4Nucleus& aNucleus) 67 { 68 G4ParticleHPManager::GetInstance()->OpenReactionWhiteBoard(); 69 const G4Material* theMaterial = aTrack.GetMaterial(); 70 auto n = (G4int)theMaterial->GetNumberOfElements(); 71 std::size_t index = theMaterial->GetElement(0)->GetIndex(); 72 if (n != 1) { 73 auto xSec = new G4double[n]; 74 G4double sum = 0; 75 G4int i; 76 const G4double* NumAtomsPerVolume = theMaterial->GetVecNbOfAtomsPerVolume(); 77 G4double rWeight; 78 G4ParticleHPThermalBoost aThermalE; 79 for (i = 0; i < n; ++i) { 80 index = theMaterial->GetElement(i)->GetIndex(); 81 rWeight = NumAtomsPerVolume[i]; 82 xSec[i] = ((*theCapture)[index])->GetXsec( 83 aThermalE.GetThermalEnergy(aTrack, theMaterial->GetElement(i), 84 theMaterial->GetTemperature())); 85 xSec[i] *= rWeight; 86 sum += xSec[i]; 87 } 88 G4double random = G4UniformRand(); 89 G4double running = 0; 90 for (i = 0; i < n; ++i) { 91 running += xSec[i]; 92 index = theMaterial->GetElement(i)->GetIndex(); 93 // if(random<=running/sum) break; 94 if (sum == 0 || random <= running / sum) break; 95 } 96 if (i == n) i = std::max(0, n - 1); 97 delete[] xSec; 98 } 99 100 G4HadFinalState* result = ((*theCapture)[index])->ApplyYourself(aTrack); 101 102 // Overwrite target parameters 103 aNucleus.SetParameters(G4ParticleHPManager::GetInstance()->GetReactionWhiteBoard()->GetTargA(), 104 G4ParticleHPManager::GetInstance()->GetReactionWhiteBoard()->GetTargZ()); 105 const G4Element* target_element = (*G4Element::GetElementTable())[index]; 106 const G4Isotope* target_isotope = nullptr; 107 auto iele = (G4int)target_element->GetNumberOfIsotopes(); 108 for (G4int j = 0; j != iele; ++j) { 109 target_isotope = target_element->GetIsotope(j); 110 if (target_isotope->GetN() 111 == G4ParticleHPManager::GetInstance()->GetReactionWhiteBoard()->GetTargA()) 112 break; 113 } 114 // G4cout << "Target Material of this reaction is " << theMaterial->GetName() << G4endl; 115 // G4cout << "Target Element of this reaction is " << target_element->GetName() << G4endl; 116 // G4cout << "Target Isotope of this reaction is " << target_isotope->GetName() << G4endl; 117 aNucleus.SetIsotope(target_isotope); 118 119 G4ParticleHPManager::GetInstance()->CloseReactionWhiteBoard(); 120 return result; 121 } 122 123 const std::pair<G4double, G4double> 124 G4NeutronHPCapture::GetFatalEnergyCheckLevels() const 125 { 126 // max energy non-conservation is mass of heavy nucleus 127 return std::pair<G4double, G4double>(10.0 * perCent, 350.0 * CLHEP::GeV); 128 } 129 130 G4int G4NeutronHPCapture::GetVerboseLevel() const 131 { 132 return G4ParticleHPManager::GetInstance()->GetVerboseLevel(); 133 } 134 135 void G4NeutronHPCapture::SetVerboseLevel(G4int newValue) 136 { 137 G4ParticleHPManager::GetInstance()->SetVerboseLevel(newValue); 138 } 139 140 void G4NeutronHPCapture::BuildPhysicsTable(const G4ParticleDefinition&) 141 { 142 G4ParticleHPManager* hpmanager = G4ParticleHPManager::GetInstance(); 143 144 theCapture = hpmanager->GetCaptureFinalStates(); 145 146 if (G4Threading::IsMasterThread()) { 147 if (theCapture == nullptr) 148 theCapture = new std::vector<G4ParticleHPChannel*>; 149 150 if (numEle == (G4int)G4Element::GetNumberOfElements()) return; 151 152 if (theCapture->size() == G4Element::GetNumberOfElements()) { 153 numEle = (G4int)G4Element::GetNumberOfElements(); 154 return; 155 } 156 157 if (G4FindDataDir("G4NEUTRONHPDATA") == nullptr) 158 throw G4HadronicException( 159 __FILE__, __LINE__, 160 "Please setenv G4NEUTRONHPDATA to point to the neutron cross-section files."); 161 dirName = G4FindDataDir("G4NEUTRONHPDATA"); 162 G4String tString = "/Capture"; 163 dirName = dirName + tString; 164 165 auto theFS = new G4NeutronHPCaptureFS; 166 for (G4int i = numEle; i < (G4int)G4Element::GetNumberOfElements(); ++i) { 167 theCapture->push_back(new G4ParticleHPChannel); 168 ((*theCapture)[i])->Init((*(G4Element::GetElementTable()))[i], dirName); 169 ((*theCapture)[i])->Register(theFS); 170 } 171 delete theFS; 172 hpmanager->RegisterCaptureFinalStates(theCapture); 173 } 174 numEle = (G4int)G4Element::GetNumberOfElements(); 175 } 176 177 void G4NeutronHPCapture::ModelDescription(std::ostream& outFile) const 178 { 179 outFile << "High Precision model based on Evaluated Nuclear Data Files (ENDF)" 180 << " for radiative capture reaction of neutrons below 20 MeV"; 181 } 182