#include <G4Fissioner.hh>
Inheritance diagram for G4Fissioner:
Public Member Functions | |
G4Fissioner () | |
virtual | ~G4Fissioner () |
void | collide (G4InuclParticle *bullet, G4InuclParticle *target, G4CollisionOutput &output) |
Definition at line 47 of file G4Fissioner.hh.
G4Fissioner::G4Fissioner | ( | ) |
virtual G4Fissioner::~G4Fissioner | ( | ) | [inline, virtual] |
void G4Fissioner::collide | ( | G4InuclParticle * | bullet, | |
G4InuclParticle * | target, | |||
G4CollisionOutput & | output | |||
) | [virtual] |
Implements G4VCascadeCollider.
Definition at line 60 of file G4Fissioner.cc.
References G4FissionStore::addConfig(), G4CollisionOutput::addOutgoingNuclei(), G4InuclSpecialFunctions::bindingEnergy(), G4InuclSpecialFunctions::bindingEnergyAsymptotic(), C1, G4FissionStore::clear(), G4InuclParticle::Fissioner, G4InuclSpecialFunctions::G4cbrt(), G4cerr, G4cout, G4endl, G4lrint(), G4FissionStore::generateConfiguration(), G4InuclSpecialFunctions::generateWithRandomAngles(), G4InuclNuclei::getA(), G4InuclNuclei::getExitationEnergy(), G4InuclParticle::getMass(), G4InuclNuclei::getNucleiMass(), G4InuclNuclei::getZ(), G4InuclSpecialFunctions::inuclRndm(), G4InuclSpecialFunctions::nucleiLevelDensity(), G4InuclSpecialFunctions::randomGauss(), G4FissionStore::setVerboseLevel(), G4FissionStore::size(), G4CascadeColliderBase::validateOutput(), and G4VCascadeCollider::verboseLevel.
Referenced by G4EquilibriumEvaporator::collide().
00062 { 00063 if (verboseLevel) { 00064 G4cout << " >>> G4Fissioner::collide" << G4endl; 00065 } 00066 00067 // const G4int itry_max = 1000; 00068 00069 G4InuclNuclei* nuclei_target = dynamic_cast<G4InuclNuclei*>(target); 00070 if (!nuclei_target) { 00071 G4cerr << " >>> G4Fissioner -> target is not nuclei " << G4endl; 00072 return; 00073 } 00074 00075 if (verboseLevel > 1) 00076 G4cout << " Fissioner input\n" << *nuclei_target << G4endl; 00077 00078 // Initialize buffer for fission possibilities 00079 fissionStore.setVerboseLevel(verboseLevel); 00080 fissionStore.clear(); 00081 00082 G4int A = nuclei_target->getA(); 00083 G4int Z = nuclei_target->getZ(); 00084 00085 G4double EEXS = nuclei_target->getExitationEnergy(); 00086 G4double mass_in = nuclei_target->getMass(); 00087 G4double e_in = mass_in; // Mass includes excitation 00088 G4double PARA = 0.055 * G4cbrt(A*A) * (G4cbrt(A-Z) + G4cbrt(Z)); 00089 G4double TEM = std::sqrt(EEXS / PARA); 00090 G4double TETA = 0.494 * G4cbrt(A) * TEM; 00091 00092 TETA = TETA / std::sinh(TETA); 00093 00094 if (A < 246) PARA += (nucleiLevelDensity(A) - PARA) * TETA; 00095 00096 G4int A1 = A/2 + 1; 00097 G4int Z1; 00098 G4int A2 = A - A1; 00099 00100 G4double ALMA = -1000.0; 00101 G4double DM1 = bindingEnergy(A,Z); 00102 G4double EVV = EEXS - DM1; 00103 G4double DM2 = bindingEnergyAsymptotic(A, Z); 00104 G4double DTEM = (A < 220 ? 0.5 : 1.15); 00105 00106 TEM += DTEM; 00107 00108 G4double AL1[2] = { -0.15, -0.15 }; 00109 G4double BET1[2] = { 0.05, 0.05 }; 00110 00111 G4double R12 = G4cbrt(A1) + G4cbrt(A2); 00112 00113 for (G4int i = 0; i < 50 && A1 > 30; i++) { 00114 A1--; 00115 A2 = A - A1; 00116 G4double X3 = 1.0 / G4cbrt(A1); 00117 G4double X4 = 1.0 / G4cbrt(A2); 00118 Z1 = G4lrint(getZopt(A1, A2, Z, X3, X4, R12) - 1.); 00119 G4double EDEF1[2]; 00120 G4int Z2 = Z - Z1; 00121 G4double VPOT, VCOUL; 00122 00123 potentialMinimization(VPOT, EDEF1, VCOUL, A1, A2, Z1, Z2, AL1, BET1, R12); 00124 00125 G4double DM3 = bindingEnergy(A1,Z1); 00126 G4double DM4 = bindingEnergyAsymptotic(A1, Z1); 00127 G4double DM5 = bindingEnergy(A2,Z2); 00128 G4double DM6 = bindingEnergyAsymptotic(A2, Z2); 00129 G4double DMT1 = DM4 + DM6 - DM2; 00130 G4double DMT = DM3 + DM5 - DM1; 00131 G4double EZL = EEXS + DMT - VPOT; 00132 00133 if(EZL > 0.0) { // generate fluctuations 00134 // faster, using randomGauss 00135 G4double C1 = std::sqrt(getC2(A1, A2, X3, X4, R12) / TEM); 00136 G4double DZ = randomGauss(C1); 00137 00138 DZ = DZ > 0.0 ? DZ + 0.5 : -std::fabs(DZ - 0.5); 00139 Z1 += G4int(DZ); 00140 Z2 -= G4int(DZ); 00141 00142 G4double DEfin = randomGauss(TEM); 00143 G4double EZ = (DMT1 + (DMT - DMT1) * TETA - VPOT + DEfin) / TEM; 00144 00145 if (EZ >= ALMA) ALMA = EZ; 00146 G4double EK = VCOUL + DEfin + 0.5 * TEM; 00147 G4double EV = EVV + bindingEnergy(A1,Z1) + bindingEnergy(A2,Z2) - EK; 00148 00149 if (EV > 0.0) fissionStore.addConfig(A1, Z1, EZ, EK, EV); 00150 }; 00151 }; 00152 00153 G4int store_size = fissionStore.size(); 00154 if (store_size == 0) return; // No fission products 00155 00156 G4FissionConfiguration config = 00157 fissionStore.generateConfiguration(ALMA, inuclRndm()); 00158 00159 A1 = G4int(config.afirst); 00160 A2 = A - A1; 00161 Z1 = G4int(config.zfirst); 00162 00163 G4int Z2 = Z - Z1; 00164 00165 G4double mass1 = G4InuclNuclei::getNucleiMass(A1,Z1); 00166 G4double mass2 = G4InuclNuclei::getNucleiMass(A2,Z2); 00167 G4double EK = config.ekin; 00168 G4double pmod = std::sqrt(0.001 * EK * mass1 * mass2 / mass_in); 00169 00170 G4LorentzVector mom1 = generateWithRandomAngles(pmod, mass1); 00171 G4LorentzVector mom2; mom2.setVectM(-mom1.vect(), mass2); 00172 00173 G4double e_out = mom1.e() + mom2.e(); 00174 G4double EV = 1000.0 * (e_in - e_out) / A; 00175 if (EV <= 0.0) return; // No fission energy 00176 00177 G4double EEXS1 = EV*A1; 00178 G4double EEXS2 = EV*A2; 00179 00180 G4InuclNuclei nuclei1(mom1, A1, Z1, EEXS1, G4InuclParticle::Fissioner); 00181 G4InuclNuclei nuclei2(mom2, A2, Z2, EEXS2, G4InuclParticle::Fissioner); 00182 00183 // Pass only last two nuclear fragments 00184 static std::vector<G4InuclNuclei> frags(2); // Always the same size! 00185 frags[0] = nuclei1; 00186 frags[1] = nuclei2; 00187 validateOutput(0, target, frags); // Check energy conservation 00188 00189 output.addOutgoingNuclei(frags); 00190 }