#include <G4LFission.hh>
Inheritance diagram for G4LFission:
Public Member Functions | |
G4LFission (const G4String &name="G4LFission") | |
~G4LFission () | |
G4HadFinalState * | ApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus) |
virtual void | ModelDescription (std::ostream &outFile) const |
virtual const std::pair< G4double, G4double > | GetFatalEnergyCheckLevels () const |
Static Public Member Functions | |
static G4double | Atomas (const G4double A, const G4double Z) |
Definition at line 65 of file G4LFission.hh.
G4LFission::G4LFission | ( | const G4String & | name = "G4LFission" |
) |
Definition at line 50 of file G4LFission.cc.
References DBL_MAX, G4HadronicInteraction::SetMaxEnergy(), and G4HadronicInteraction::SetMinEnergy().
00051 : G4HadronicInteraction(name) 00052 { 00053 init(); 00054 SetMinEnergy(0.0*GeV); 00055 SetMaxEnergy(DBL_MAX); 00056 }
G4LFission::~G4LFission | ( | ) |
Definition at line 59 of file G4LFission.cc.
References G4HadFinalState::Clear(), and G4HadronicInteraction::theParticleChange.
00060 { 00061 theParticleChange.Clear(); 00062 }
G4HadFinalState * G4LFission::ApplyYourself | ( | const G4HadProjectile & | aTrack, | |
G4Nucleus & | targetNucleus | |||
) | [virtual] |
Implements G4HadronicInteraction.
Definition at line 94 of file G4LFission.cc.
References G4HadFinalState::AddSecondary(), Atomas(), G4HadFinalState::Clear(), G4cout, G4endl, G4UniformRand, G4Gamma::GammaDefinition(), G4HadProjectile::Get4Momentum(), G4Nucleus::GetA_asInt(), G4DynamicParticle::GetDefinition(), G4HadProjectile::GetDefinition(), G4HadProjectile::GetKineticEnergy(), G4HadSecondary::GetParticle(), G4ParticleDefinition::GetPDGCharge(), G4ParticleDefinition::GetPDGMass(), G4HadFinalState::GetSecondary(), G4DynamicParticle::GetTotalEnergy(), G4HadProjectile::GetTotalEnergy(), G4DynamicParticle::GetTotalMomentum(), G4HadProjectile::GetTotalMomentum(), G4Nucleus::GetZ_asInt(), G4Neutron::NeutronDefinition(), G4InuclParticleNames::nn, G4InuclParticleNames::pp, G4DynamicParticle::SetKineticEnergy(), G4DynamicParticle::SetMomentumDirection(), G4HadFinalState::SetStatusChange(), stopAndKill, G4HadronicInteraction::theParticleChange, and G4HadronicInteraction::verboseLevel.
Referenced by G4NeutronHPorLFissionModel::ApplyYourself().
00096 { 00097 theParticleChange.Clear(); 00098 const G4HadProjectile* aParticle = &aTrack; 00099 00100 G4double N = targetNucleus.GetA_asInt(); 00101 G4double Z = targetNucleus.GetZ_asInt(); 00102 theParticleChange.SetStatusChange(stopAndKill); 00103 00104 G4double P = aParticle->GetTotalMomentum()/MeV; 00105 G4double Px = aParticle->Get4Momentum().vect().x(); 00106 G4double Py = aParticle->Get4Momentum().vect().y(); 00107 G4double Pz = aParticle->Get4Momentum().vect().z(); 00108 G4double E = aParticle->GetTotalEnergy()/MeV; 00109 G4double E0 = aParticle->GetDefinition()->GetPDGMass()/MeV; 00110 G4double Q = aParticle->GetDefinition()->GetPDGCharge(); 00111 if (verboseLevel > 1) { 00112 G4cout << "G4LFission:ApplyYourself: incident particle:" << G4endl; 00113 G4cout << "P " << P << " MeV/c" << G4endl; 00114 G4cout << "Px " << Px << " MeV/c" << G4endl; 00115 G4cout << "Py " << Py << " MeV/c" << G4endl; 00116 G4cout << "Pz " << Pz << " MeV/c" << G4endl; 00117 G4cout << "E " << E << " MeV" << G4endl; 00118 G4cout << "mass " << E0 << " MeV" << G4endl; 00119 G4cout << "charge " << Q << G4endl; 00120 } 00121 // GHEISHA ADD operation to get total energy, mass, charge: 00122 if (verboseLevel > 1) { 00123 G4cout << "G4LFission:ApplyYourself: material:" << G4endl; 00124 G4cout << "A " << N << G4endl; 00125 G4cout << "Z " << Z << G4endl; 00126 G4cout << "atomic mass " << 00127 Atomas(N, Z) << "MeV" << G4endl; 00128 } 00129 E = E + Atomas(N, Z); 00130 G4double E02 = E*E - P*P; 00131 E0 = std::sqrt(std::abs(E02)); 00132 if (E02 < 0) E0 = -E0; 00133 Q = Q + Z; 00134 if (verboseLevel > 1) { 00135 G4cout << "G4LFission:ApplyYourself: total:" << G4endl; 00136 G4cout << "E " << E << " MeV" << G4endl; 00137 G4cout << "mass " << E0 << " MeV" << G4endl; 00138 G4cout << "charge " << Q << G4endl; 00139 } 00140 Px = -Px; 00141 Py = -Py; 00142 Pz = -Pz; 00143 00144 G4double e1 = aParticle->GetKineticEnergy()/MeV; 00145 if (e1 < 1.) e1 = 1.; 00146 00147 // Average number of neutrons 00148 G4double avern = 2.569 + 0.559*std::log(e1); 00149 G4bool photofission = 0; // For now 00150 // Take the following value if photofission is not included 00151 if (!photofission) avern = 2.569 + 0.900*std::log(e1); 00152 00153 // Average number of gammas 00154 G4double averg = 9.500 + 0.600*std::log(e1); 00155 00156 G4double ran = G4RandGauss::shoot(); 00157 // Number of neutrons 00158 G4int nn = static_cast<G4int>(avern + ran*1.23 + 0.5); 00159 ran = G4RandGauss::shoot(); 00160 // Number of gammas 00161 G4int ng = static_cast<G4int>(averg + ran*3. + 0.5); 00162 if (nn < 1) nn = 1; 00163 if (ng < 1) ng = 1; 00164 G4double exn = 0.; 00165 G4double exg = 0.; 00166 00167 // Make secondary neutrons and distribute kinetic energy 00168 G4DynamicParticle* aNeutron; 00169 G4int i; 00170 for (i = 1; i <= nn; i++) { 00171 ran = G4UniformRand(); 00172 G4int j; 00173 for (j = 1; j <= 10; j++) { 00174 if (ran < spneut[j-1]) goto label12; 00175 } 00176 j = 10; 00177 label12: 00178 ran = G4UniformRand(); 00179 G4double ekin = (j - 1)*1. + ran; 00180 exn = exn + ekin; 00181 aNeutron = new G4DynamicParticle(G4Neutron::NeutronDefinition(), 00182 G4ParticleMomentum(1.,0.,0.), 00183 ekin*MeV); 00184 theParticleChange.AddSecondary(aNeutron); 00185 } 00186 00187 // Make secondary gammas and distribute kinetic energy 00188 G4DynamicParticle* aGamma; 00189 for (i = 1; i <= ng; i++) { 00190 ran = G4UniformRand(); 00191 G4double ekin = -0.87*std::log(ran); 00192 exg = exg + ekin; 00193 aGamma = new G4DynamicParticle(G4Gamma::GammaDefinition(), 00194 G4ParticleMomentum(1.,0.,0.), 00195 ekin*MeV); 00196 theParticleChange.AddSecondary(aGamma); 00197 } 00198 00199 // Distribute momentum vectors and do Lorentz transformation 00200 00201 G4HadSecondary* theSecondary; 00202 00203 for (i = 1; i <= nn + ng; i++) { 00204 G4double ran1 = G4UniformRand(); 00205 G4double ran2 = G4UniformRand(); 00206 G4double cost = -1. + 2.*ran1; 00207 G4double sint = std::sqrt(std::abs(1. - cost*cost)); 00208 G4double phi = ran2*twopi; 00209 // G4cout << ran1 << " " << ran2 << G4endl; 00210 // G4cout << cost << " " << sint << " " << phi << G4endl; 00211 theSecondary = theParticleChange.GetSecondary(i - 1); 00212 G4double pp = theSecondary->GetParticle()->GetTotalMomentum()/MeV; 00213 G4double px = pp*sint*std::sin(phi); 00214 G4double py = pp*sint*std::cos(phi); 00215 G4double pz = pp*cost; 00216 // G4cout << pp << G4endl; 00217 // G4cout << px << " " << py << " " << pz << G4endl; 00218 G4double e = theSecondary->GetParticle()->GetTotalEnergy()/MeV; 00219 G4double e0 = theSecondary->GetParticle()->GetDefinition()->GetPDGMass()/MeV; 00220 00221 G4double a = px*Px + py*Py + pz*Pz; 00222 a = (a/(E + E0) - e)/E0; 00223 00224 px = px + a*Px; 00225 py = py + a*Py; 00226 pz = pz + a*Pz; 00227 G4double p2 = px*px + py*py + pz*pz; 00228 pp = std::sqrt(p2); 00229 e = std::sqrt(e0*e0 + p2); 00230 G4double ekin = e - theSecondary->GetParticle()->GetDefinition()->GetPDGMass()/MeV; 00231 theSecondary->GetParticle()->SetMomentumDirection(G4ParticleMomentum(px/pp, 00232 py/pp, 00233 pz/pp)); 00234 theSecondary->GetParticle()->SetKineticEnergy(ekin*MeV); 00235 } 00236 00237 return &theParticleChange; 00238 }
Definition at line 243 of file G4LFission.cc.
References G4Alpha::AlphaDefinition(), G4Deuteron::DeuteronDefinition(), G4Electron::ElectronDefinition(), G4ParticleDefinition::GetPDGMass(), G4Neutron::NeutronDefinition(), and G4Proton::ProtonDefinition().
Referenced by ApplyYourself().
00244 { 00245 G4double rmel = G4Electron::ElectronDefinition()->GetPDGMass()/MeV; 00246 G4double rmp = G4Proton::ProtonDefinition()->GetPDGMass()/MeV; 00247 G4double rmn = G4Neutron::NeutronDefinition()->GetPDGMass()/MeV; 00248 G4double rmd = G4Deuteron::DeuteronDefinition()->GetPDGMass()/MeV; 00249 G4double rma = G4Alpha::AlphaDefinition()->GetPDGMass()/MeV; 00250 00251 G4int ia = static_cast<G4int>(A + 0.5); 00252 if (ia < 1) return 0; 00253 G4int iz = static_cast<G4int>(Z + 0.5); 00254 if (iz < 0) return 0; 00255 if (iz > ia) return 0; 00256 00257 if (ia == 1) { 00258 if (iz == 0) return rmn; //neutron 00259 if (iz == 1) return rmp + rmel; //Hydrogen 00260 } 00261 else if (ia == 2 && iz == 1) { 00262 return rmd; //Deuteron 00263 } 00264 else if (ia == 4 && iz == 2) { 00265 return rma; //Alpha 00266 } 00267 00268 G4double mass = (A - Z)*rmn + Z*rmp + Z*rmel - 15.67*A 00269 + 17.23*std::pow(A, 2./3.) 00270 + 93.15*(A/2. - Z)*(A/2. - Z)/A 00271 + 0.6984523*Z*Z/std::pow(A, 1./3.); 00272 G4int ipp = (ia - iz)%2; 00273 G4int izz = iz%2; 00274 if (ipp == izz) mass = mass + (ipp + izz -1)*12.*std::pow(A, -0.5); 00275 00276 return mass; 00277 }
Reimplemented from G4HadronicInteraction.
Definition at line 279 of file G4LFission.cc.
00280 { 00281 // max energy non-conservation is mass of heavy nucleus 00282 return std::pair<G4double, G4double>(5*perCent,250*GeV); 00283 }
void G4LFission::ModelDescription | ( | std::ostream & | outFile | ) | const [virtual] |
Reimplemented from G4HadronicInteraction.
Definition at line 65 of file G4LFission.cc.
00066 { 00067 outFile << "G4LFission is one of the Low Energy Parameterized\n" 00068 << "(LEP) models used to implement neutron-induced fission of\n" 00069 << "nuclei. It is a re-engineered version of the GHEISHA code\n" 00070 << "of H. Fesefeldt which emits neutrons and gammas but no\n" 00071 << "nuclear fragments. The model is applicable to all incident\n" 00072 << "neutron energies.\n"; 00073 }