76 G4int IterationsLimit = 100000;
87 if (theMeanMult <= MaxAverageMultiplicity) {
90 theChannel = theMicrocanonicalEnsemble->
ChooseAandZ(theFragment);
105 theChannel = theMacrocanonicalEnsemble->
ChooseAandZ(theFragment);
109 if (!ChannelOk)
delete theChannel;
112 }
while (!ChannelOk);
117 theResult->push_back(
new G4Fragment(theFragment));
118 delete theMicrocanonicalEnsemble;
119 if (theMacrocanonicalEnsemble != 0)
delete theMacrocanonicalEnsemble;
142 }
while (Iterations++ < IterationsLimit );
147 if (Iterations >= IterationsLimit)
148 throw G4HadronicException(__FILE__, __LINE__,
"G4StatMF::BreakItUp: Was not possible to solve for temperature of breaking channel");
164 G4FragmentVector::iterator j;
165 for (j = theResult->begin(); j != theResult->end(); j++)
166 FragmentsEnergy += (*j)->GetMomentum().
e();
167 SavedScaleFactor = ScaleFactor;
168 ScaleFactor = InitialMomentum.
e()/FragmentsEnergy;
170 for (j = theResult->begin(); j != theResult->end(); j++) {
171 ScaledMomentum = ScaleFactor * (*j)->GetMomentum().vect();
172 G4double Mass = (*j)->GetMomentum().m();
174 NewMomentum.
setVect(ScaledMomentum);
175 NewMomentum.
setE(std::sqrt(ScaledMomentum.
mag2()+Mass*Mass));
176 (*j)->SetMomentum(NewMomentum);
179 }
while (ScaleFactor > 1.0+1.e-5 && std::abs(ScaleFactor-SavedScaleFactor)/ScaleFactor > 1.e-10);
183 G4FragmentVector::iterator i;
184 for (i = theResult->begin(); i != theResult->end(); i++) {
187 (*i)->SetMomentum(FourMom);
188 (*i)->SetCreatorModelID(
_secID);
192 delete theMicrocanonicalEnsemble;
193 if (theMacrocanonicalEnsemble != 0)
delete theMacrocanonicalEnsemble;
220 }
else if (Da < 0.0) {
223 if (T < 0.001*
MeV)
return false;
227 Db = (U - TotalEnergy)/U;
237 Db = (U - TotalEnergy)/U;
246 for (
G4int j = 0; j < 1000; j++) {
248 if (std::abs(Ta-Tc) <=
eps) {
273 Temperature = (Ta+T)*0.5;
static const G4double eps
std::vector< G4Fragment * > G4FragmentVector
static constexpr double MeV
Hep3Vector boostVector() const
HepLorentzVector & boost(double, double, double)
void setVect(const Hep3Vector &)
G4double GetExcitationEnergy() const
const G4LorentzVector & GetMomentum() const
static G4double GetMassExcess(const G4int A, const G4int Z)
static G4int GetModelID(const G4int modelIndex)
G4bool CheckFragments(void)
G4double GetFragmentsEnergy(G4double T) const
size_t GetMultiplicity(void)
G4StatMFChannel * ChooseAandZ(const G4Fragment &theFragment)
G4StatMFChannel * ChooseAandZ(const G4Fragment &theFragment)
static G4double GetMaxAverageMultiplicity(G4int A)
static G4double GetCoulomb()
G4FragmentVector * BreakItUp(const G4Fragment &theNucleus)
G4VStatMFEnsemble * _theEnsemble
G4double CalcEnergy(G4int A, G4int Z, const G4StatMFChannel *aChannel, G4double T)
G4bool FindTemperatureOfBreakingChannel(const G4Fragment &theFragment, const G4StatMFChannel *aChannel, G4double &Temperature)
G4double GetMeanTemperature(void) const
G4double GetMeanMultiplicity(void) const
T max(const T t1, const T t2)
brief Return the largest of the two arguments