Geant4 Cross Reference |
1 // 1 2 // ******************************************* 3 // * License and Disclaimer 4 // * 5 // * The Geant4 software is copyright of th 6 // * the Geant4 Collaboration. It is provided 7 // * conditions of the Geant4 Software License 8 // * LICENSE and available at http://cern.ch/ 9 // * include a list of copyright holders. 10 // * 11 // * Neither the authors of this software syst 12 // * institutes,nor the agencies providing fin 13 // * work make any representation or warran 14 // * regarding this software system or assum 15 // * use. Please see the license in the file 16 // * for the full disclaimer and the limitatio 17 // * 18 // * This code implementation is the result 19 // * technical work of the GEANT4 collaboratio 20 // * By using, copying, modifying or distri 21 // * any work based on the software) you ag 22 // * use in resulting scientific publicati 23 // * acceptance of all terms of the Geant4 Sof 24 // ******************************************* 25 // 26 // 27 // 28 // Physics model class G4NeutronFissionVI 29 // Created: 03 October 2023 30 // Author V.Ivanchenko 31 // 32 33 #include "G4NeutronFissionVI.hh" 34 #include "G4HadronicInteractionRegistry.hh" 35 #include "G4SystemOfUnits.hh" 36 #include "G4ParticleDefinition.hh" 37 #include "G4Fragment.hh" 38 #include "G4FragmentVector.hh" 39 #include "G4NucleiProperties.hh" 40 #include "G4VEvaporation.hh" 41 #include "G4VEvaporationChannel.hh" 42 #include "G4VPreCompoundModel.hh" 43 #include "G4ExcitationHandler.hh" 44 #include "G4PhotonEvaporation.hh" 45 #include "G4DynamicParticle.hh" 46 #include "G4ParticleHPManager.hh" 47 #include "G4ParticleTable.hh" 48 #include "G4IonTable.hh" 49 #include "G4Electron.hh" 50 #include "G4Deuteron.hh" 51 #include "G4Triton.hh" 52 #include "G4He3.hh" 53 #include "G4Alpha.hh" 54 #include "Randomize.hh" 55 #include "G4RandomDirection.hh" 56 #include "G4HadronicParameters.hh" 57 #include "G4PhysicsModelCatalog.hh" 58 59 G4NeutronFissionVI::G4NeutronFissionVI() 60 : G4HadronicInteraction("nFissionVI"), 61 fManagerHP(G4ParticleHPManager::GetInstanc 62 minExcitation(0.1*CLHEP::keV), 63 emaxT(fManagerHP->GetMaxEnergyDoppler()), 64 lab4mom(0.,0.,0.,0.) 65 { 66 SetMinEnergy( 0.0*CLHEP::GeV ); 67 SetMaxEnergy( G4HadronicParameters::Instance 68 theTableOfIons = G4ParticleTable::GetParticl 69 } 70 71 G4NeutronFissionVI::~G4NeutronFissionVI() 72 { 73 if (fLocalHandler) { delete fHandler; } 74 } 75 76 void G4NeutronFissionVI::InitialiseModel() 77 { 78 if (fFission != nullptr && fHandler != nullp 79 G4HadronicInteraction* p = 80 G4HadronicInteractionRegistry::Instance()- 81 if (nullptr != p) { 82 fHandler = (static_cast<G4VPreCompoundMode 83 } 84 if (nullptr == fHandler) { 85 fHandler = new G4ExcitationHandler(); 86 fLocalHandler = true; 87 } 88 fHandler->Initialise(); 89 fFission = fHandler->GetEvaporation()->GetFi 90 91 G4DeexPrecoParameters* param = 92 G4NuclearLevelData::GetInstance()->GetPara 93 minExcitation = param->GetMinExcitation(); 94 secID = G4PhysicsModelCatalog::GetModelID("m 95 } 96 97 G4HadFinalState* G4NeutronFissionVI::ApplyYour 98 const G4HadProjectile& aTrack, G4Nucleus& 99 { 100 theParticleChange.Clear(); 101 theParticleChange.SetStatusChange(stopAndKil 102 G4double T = aTrack.GetMaterial()->GetTemper 103 104 G4int A = theNucleus.GetA_asInt(); 105 G4int Z = theNucleus.GetZ_asInt(); 106 107 G4double time = aTrack.GetGlobalTime(); 108 G4double ekin = aTrack.GetKineticEnergy(); 109 110 // Create initial state 111 G4double mass = G4NucleiProperties::GetNucle 112 113 // no Doppler broading 114 G4double factT = T/CLHEP::STP_Temperature; 115 116 if (ekin >= emaxT*factT || fManagerHP->GetNe 117 lab4mom.set(0.,0.,0.,mass); 118 119 } else { 120 G4double lambda = 1.0/(CLHEP::k_Boltzmann* 121 G4double erand = G4RandGamma::shoot(2.0, l 122 auto mom = G4RandomDirection()*std::sqrt(2 123 lab4mom.set(mom.x(), mom.y(), mom.z(), mas 124 } 125 126 lab4mom += aTrack.Get4Momentum(); 127 128 G4double M = lab4mom.mag(); 129 ++A; 130 mass = G4NucleiProperties::GetNuclearMass(A, 131 //G4cout << "Fission start: Z= " << Z << " A 132 // << " LabM= " << M << " Mcompound= " << 133 134 // protection against wrong kinematic 135 if (M < mass) { 136 G4double etot = std::max(mass, lab4mom.e() 137 G4double ptot = std::sqrt((etot - mass)*(e 138 G4ThreeVector v = lab4mom.vect().unit(); 139 lab4mom.set(v.x()*ptot,v.y()*ptot,v.z()*pt 140 } 141 142 G4Fragment* aFragment = new G4Fragment(A, Z, 143 144 if (verboseLevel > 1) { 145 G4cout << "G4NeutronFissionVI::ApplyYourse 146 << G4endl; 147 G4cout << aFragment << G4endl; 148 } 149 150 // 151 // Sample final state 152 // 153 fFission->GetEmissionProbability(aFragment); 154 G4Fragment* frag = fFission->EmittedFragment 155 G4ReactionProductVector* final = fHandler->B 156 if (nullptr != frag) { 157 G4ReactionProductVector* v = fHandler->Bre 158 for (auto & p : *v) { 159 final->push_back(p); 160 } 161 delete v; 162 delete frag; 163 } 164 165 if (verboseLevel > 1) { 166 G4cout << "G4NeutronFissionVI: " << fina 167 << " final particle secID= " << s 168 } 169 for (auto const & ptr : *final) { 170 G4double etot = ptr->GetTotalEnergy(); 171 const G4ParticleDefinition* theDef = ptr-> 172 ekin = std::max(0.0, etot - theDef->GetPDG 173 if (verboseLevel > 1) { 174 G4cout << theDef->GetParticleName() 175 << " Ekin(MeV)= " << ekin/MeV 176 << " p: " << ptr->GetMomentum() 177 << G4endl; 178 } 179 G4HadSecondary* news = new G4HadSecondary( 180 new G4DynamicParticle(theDef, ptr->GetMo 181 G4double timeF = std::max(ptr->GetFormatio 182 news->SetTime(time + timeF); 183 news->SetCreatorModelID(secID); 184 theParticleChange.AddSecondary(*news); 185 delete news; 186 delete ptr; 187 } 188 delete final; 189 delete aFragment; 190 191 //G4cout << "Fission done" << G4endl; 192 return &theParticleChange; 193 } 194 195