Attachment #30 for problem #1001

Lines 105-117 G4FragmentVector* G4VGammaDeexcitation:: Link Here

(-)geant4.9.2.b01/source/processes/hadronic/models/de_excitation/photon_evaporation/src/G4VGammaDeexcitation.cc (-1 / +10 lines)
105	Initialize();	106	Initialize();
106	G4FragmentVector* products = new G4FragmentVector;	107	G4FragmentVector* products = new G4FragmentVector;
107		108
108	while (CanDoTransition())	109	// --> acz 200808-27: Store the number of occured transitions, to limit the number of transitions to one per step
		110	// Reasons: we can only have one vaccant shell position and we want to give the user access at each transision
		111	G4int transitionNumber = 0;
		112
		113	while (transitionNumber == 0 && CanDoTransition())
		114	// <--
109	{	115	{
110	if (_verbose > 5) G4cout << "G4VGammaDeexcitation::DoChain - Looping" << G4endl;	116	if (_verbose > 5) G4cout << "G4VGammaDeexcitation::DoChain - Looping" << G4endl;
111		117
112	G4Fragment* gamma = GenerateGamma();	118	G4Fragment* gamma = GenerateGamma();
113	if (gamma != 0)	119	if (gamma != 0)
114	{	120	{
		121	// --> acz 200808-27
		122	transitionNumber++;
		123	// <--
		124
115	products->push_back(gamma);	125	products->push_back(gamma);
116	UpdateNucleus(gamma);	126	UpdateNucleus(gamma);
117	UpdateElectrons ();	127	UpdateElectrons ();

Lines 292-298 G4DecayProducts *G4NuclearDecayChannel:: Link Here

(-)geant4.9.2.b01/source/processes/hadronic/models/radioactive_decay/src/G4NuclearDecayChannel.cc (-3 / +22 lines)
292	// needed to hold the shell idex after ICM	292	// needed to hold the shell idex after ICM
293	G4int shellIndex = -1;	293	G4int shellIndex = -1;
294	//	294	//
295	if (daughterExcitation > 0.0)	295	// --> acz 2008-08-27:
		296	// This hack was necessary, since we are only want to do one transition at a time
		297	// Thus we are only de-exiting the nucleus if the parent is identical as the daughter
		298	// if (daughterExcitation > 0.0)
		299	if (daughterExcitation > 0.0 && parent->GetParticleName() == daughters[0]->GetParticleName())
		300	// <--
296	{	301	{
297	//	302	//
298	// Pop the daughter nucleus off the product vector - we need to retain	303	// Pop the daughter nucleus off the product vector - we need to retain
Lines 314-320 G4DecayProducts *G4NuclearDecayChannel:: Link Here
314	// set the maximum life-time for a level that will be treated. Level with life-time longer than this	319	// set the maximum life-time for a level that will be treated. Level with life-time longer than this
315	// will be outputed as meta-stable isotope	320	// will be outputed as meta-stable isotope
316	//	321	//
317	deexcitation->SetMaxHalfLife(1e-6*second);	322	// --> acz 2008-08-27
		323	deexcitation->SetMaxHalfLife(1e-9*second);
		324	// <--
318	// but in IT mode, we need to force the transition	325	// but in IT mode, we need to force the transition
319	if (decayMode == 0) {	326	if (decayMode == 0) {
320	deexcitation->RDMForced(true);	327	deexcitation->RDMForced(true);
Lines 368-373 G4DecayProducts *G4NuclearDecayChannel:: Link Here
368	// gammas->clearAndDestroy();	375	// gammas->clearAndDestroy();
369	delete gammas;	376	delete gammas;
370	delete deexcitation;	377	delete deexcitation;
		378
		379	// --> acz 2008-09-29:
		380	// Some protection against infinite loops: If the daughter is still the same as the parent, then we are are stuck in an infinite loop
		381	if (products->entries() == 1) {
		382	dynamicDaughter = products->PopProducts();
		383	if (parent->GetParticleName() == dynamicDaughter->GetDefinition()->GetParticleName()) {
		384	delete dynamicDaughter;
		385	return products;
		386	} else {
		387	products->PushProducts (dynamicDaughter);
		388	}
		389	}
		390	// <--
371	}	391	}
372	//	392	//
373	// now we have to take care of the EC product which have to go through the ARM	393	// now we have to take care of the EC product which have to go through the ARM
Lines 401-407 G4DecayProducts *G4NuclearDecayChannel:: Link Here
401	}	421	}
402	}	422	}
403	// now deal with the IT case where ICM may have been applied	423	// now deal with the IT case where ICM may have been applied
404	//
405	if (decayMode == 0) {	424	if (decayMode == 0) {
406	eShell = shellIndex;	425	eShell = shellIndex;
407	}	426	}

Lines 1344-1352 G4VParticleChange* G4RadioactiveDecay::D Link Here

(-)geant4.9.2.b01/source/processes/hadronic/models/radioactive_decay/src/G4RadioactiveDecay.cc (-7 / +75 lines)
1344	return &fParticleChangeForRadDecay;	1344	return &fParticleChangeForRadDecay;
1345	}	1345	}
1346		1346
		1347	// --> acz: 2008-09-25
		1348	// Some problems are related to bad handling of excited states in the current data
		1349	// This determine a flag here, where or not we have an excited nucleus
		1350	bool ExcitedNucleus = false;
		1351	if (theParticleDef->GetParticleType() == "nucleus") {
		1352	if (dynamic_cast<G4Ions*>(theParticleDef)->GetExcitationEnergy() > 0.0) {
		1353	ExcitedNucleus = true;
		1354	}
		1355	}
		1356	// <--
		1357
1347	// now check is the particle is valid for RDM	1358	// now check is the particle is valid for RDM
1348	//	1359	//
1349	if (!(IsApplicable(*theParticleDef)))	1360	if (!(IsApplicable(*theParticleDef)) && ExcitedNucleus == false)
1350	{	1361	{
1351	//	1362	//
1352	// The particle is not a Ion or outside the nucleuslimits for decay	1363	// The particle is not a Ion or outside the nucleuslimits for decay
Lines 1376-1382 G4VParticleChange* G4RadioactiveDecay::D Link Here
1376	}	1387	}
1377	G4DecayTable *theDecayTable = theParticleDef->GetDecayTable();	1388	G4DecayTable *theDecayTable = theParticleDef->GetDecayTable();
1378		1389
1379	if (theDecayTable == 0 \|\| theDecayTable->entries() == 0 )	1390	// --> acz: The decay table might be empty for some excited states
		1391	// if ((theDecayTable == 0 \|\| theDecayTable->entries() == 0))
		1392	if ((theDecayTable == 0 \|\| theDecayTable->entries() == 0) && ExcitedNucleus == false)
		1393	// <--
1380	{	1394	{
1381	//	1395	//
1382	//	1396	//
Lines 1424-1429 G4VParticleChange* G4RadioactiveDecay::D Link Here
1424	#endif	1438	#endif
1425	//	1439	//
1426	G4DecayProducts products = DoDecay(theParticleDef);	1440	G4DecayProducts products = DoDecay(theParticleDef);
		1441
		1442	// --> acz: Prevent infinite loop
		1443	// If there are no secondaries, kill the track:
		1444	if (products->entries() == 0) {
		1445	fParticleChangeForRadDecay.SetNumberOfSecondaries(0);
		1446	fParticleChangeForRadDecay.ProposeTrackStatus( fStopAndKill ) ;
		1447	fParticleChangeForRadDecay.ProposeLocalEnergyDeposit(0.0);
		1448	ClearNumberOfInteractionLengthLeft();
		1449	return &fParticleChangeForRadDecay;
		1450	}
		1451	// <--
		1452
1427	//	1453	//
1428	//	1454	//
1429	// Get parent particle information and boost the decay products to the	1455	// Get parent particle information and boost the decay products to the
Lines 1442-1448 G4VParticleChange* G4RadioactiveDecay::D Link Here
1442	// time lapsed between the particle come to rest and the actual decay. This time	1468	// time lapsed between the particle come to rest and the actual decay. This time
1443	// is simply sampled with the mean-life of the particle.	1469	// is simply sampled with the mean-life of the particle.
1444	//	1470	//
		1471	// --> acz 2008-09-17: Bug: PDGLifeTime can be zero for some excited states!
		1472	if (theParticleDef->GetPDGLifeTime() <= 0.0 && ExcitedNucleus == true) {
		1473	if (dynamic_cast<G4Ions*>(theParticleDef)->GetExcitationEnergy() > 0.0) {
		1474	G4NuclearLevelManager* M =
		1475	G4NuclearLevelStore::GetInstance()->GetManager(dynamic_cast<G4Ions*>(theParticleDef)->GetAtomicNumber(),
		1476	dynamic_cast<G4Ions*>(theParticleDef)->GetAtomicMass());
		1477	if (M->IsValid() == true) {
		1478	const G4NuclearLevel* Level = M->NearestLevel(dynamic_cast<G4Ions*>(theParticleDef)->GetExcitationEnergy());
		1479	if (Level != 0) {
		1480	theParticleDef->SetPDGLifeTime(Level->HalfLife()/log(2));
		1481	}
		1482	}
		1483	}
		1484	}
		1485	if (theParticleDef->GetPDGLifeTime() > 0.0) {
1445	finalGlobalTime += -std::log( G4UniformRand()) * theParticleDef->GetPDGLifeTime() ;	1486	finalGlobalTime += -std::log( G4UniformRand()) * theParticleDef->GetPDGLifeTime() ;
		1487	}
		1488	// <-- acz
1446	energyDeposit += theParticle->GetKineticEnergy();	1489	energyDeposit += theParticle->GetKineticEnergy();
1447	}	1490	}
1448	else	1491	else
Lines 1731-1736 G4DecayProducts* G4RadioactiveDecay::DoD Link Here
1731	G4cout <<"Selecte a channel..."<<G4endl;	1774	G4cout <<"Selecte a channel..."<<G4endl;
1732	}	1775	}
1733	#endif	1776	#endif
		1777
		1778	// --> acz 2008-09-26: Some bug somewhere: No decay table for some excited states, thus we have to fake one here
		1779	if (theDecayTable == 0) { //
		1780	bool ExcitedNucleus = false;
		1781	if (theParticleDef.GetParticleType() == "nucleus") {
		1782	if (dynamic_cast<G4Ions*>(&theParticleDef)->GetExcitationEnergy() > 0.0) {
		1783	ExcitedNucleus = true;
		1784	}
		1785	}
		1786
		1787	if (ExcitedNucleus == true) {
		1788	G4NuclearDecayChannel* Channel = new G4NuclearDecayChannel(IT,
		1789	GetVerboseLevel(),
		1790	&theParticleDef,
		1791	1.0,
		1792	0.0,
		1793	theParticleDef.GetAtomicMass(),
		1794	theParticleDef.GetAtomicNumber(),
		1795	dynamic_cast<G4Ions*>(&theParticleDef)->GetExcitationEnergy());
		1796
		1797	G4double tempmass = theParticleDef.GetPDGMass();
		1798	products = Channel->DecayIt(tempmass);
		1799
		1800	delete Channel;
		1801	}
		1802
		1803	} else {
1734	G4VDecayChannel *theDecayChannel = theDecayTable->SelectADecayChannel();	1804	G4VDecayChannel *theDecayChannel = theDecayTable->SelectADecayChannel();
1735	if (theDecayChannel == 0)	1805	if (theDecayChannel == 0)
1736	{	1806	{
Lines 1738-1744 G4DecayProducts* G4RadioactiveDecay::DoD Link Here
1738	//	1808	//
1739	G4cerr <<"G4RadioactiveDecay::DoIt : can not determine decay channel";	1809	G4cerr <<"G4RadioactiveDecay::DoIt : can not determine decay channel";
1740	G4cerr <<G4endl;	1810	G4cerr <<G4endl;
1741	theDecayTable ->DumpInfo();
1742	}	1811	}
1743	else	1812	else
1744	{	1813	{
Lines 1754-1768 G4DecayProducts* G4RadioactiveDecay::DoD Link Here
1754	#endif	1823	#endif
1755		1824
1756	G4double tempmass = theParticleDef.GetPDGMass();	1825	G4double tempmass = theParticleDef.GetPDGMass();
1757	//
1758
1759	products = theDecayChannel->DecayIt(tempmass);	1826	products = theDecayChannel->DecayIt(tempmass);
1760		1827
1761	}	1828	}
1762	return products;
1763
1764	}	1829	}
		1830	// <--
1765		1831
		1832	return products;
		1833	}
1766		1834
1767		1835
1768		1836

Return to problem 1001