Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/de_excitation/multifragmentation/src/G4StatMFMicroCanonical.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /processes/hadronic/models/de_excitation/multifragmentation/src/G4StatMFMicroCanonical.cc (Version 11.3.0) and /processes/hadronic/models/de_excitation/multifragmentation/src/G4StatMFMicroCanonical.cc (Version 9.6.p3)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                         3 // * License and Disclaimer                                           *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th      5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided      6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License      7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/      8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.           9 // * include a list of copyright holders.                             *
 10 // *                                               10 // *                                                                  *
 11 // * Neither the authors of this software syst     11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file      15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitatio     16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                               17 // *                                                                  *
 18 // * This  code  implementation is the result      18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboratio     19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distri     20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  ag     21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publicati     22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Sof     23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // *******************************************     24 // ********************************************************************
 25 //                                                 25 //
 26 //                                                 26 //
                                                   >>  27 // $Id$
 27 //                                                 28 //
 28 // Hadronic Process: Nuclear De-excitations        29 // Hadronic Process: Nuclear De-excitations
 29 // by V. Lara                                      30 // by V. Lara
 30                                                    31 
 31 #include <numeric>                                 32 #include <numeric>
 32                                                    33 
 33 #include "G4StatMFMicroCanonical.hh"               34 #include "G4StatMFMicroCanonical.hh"
 34 #include "G4PhysicalConstants.hh"                  35 #include "G4PhysicalConstants.hh"
 35 #include "G4SystemOfUnits.hh"                      36 #include "G4SystemOfUnits.hh"
 36 #include "G4HadronicException.hh"                  37 #include "G4HadronicException.hh"
 37 #include "G4Pow.hh"                                38 #include "G4Pow.hh"
 38                                                    39 
 39 // constructor                                     40 // constructor
 40 G4StatMFMicroCanonical::G4StatMFMicroCanonical     41 G4StatMFMicroCanonical::G4StatMFMicroCanonical(G4Fragment const & theFragment) 
 41 {                                                  42 {
 42   // Perform class initialization              <<  43     // Perform class initialization
 43   Initialize(theFragment);                     <<  44     Initialize(theFragment);
                                                   >>  45 
 44 }                                                  46 }
 45                                                    47 
                                                   >>  48 
 46 // destructor                                      49 // destructor
 47 G4StatMFMicroCanonical::~G4StatMFMicroCanonica     50 G4StatMFMicroCanonical::~G4StatMFMicroCanonical() 
 48 {                                                  51 {
 49   // garbage collection                            52   // garbage collection
 50   if (!_ThePartitionManagerVector.empty()) {       53   if (!_ThePartitionManagerVector.empty()) {
 51     std::for_each(_ThePartitionManagerVector.b     54     std::for_each(_ThePartitionManagerVector.begin(),
 52         _ThePartitionManagerVector.end(),          55         _ThePartitionManagerVector.end(),
 53         DeleteFragment());                         56         DeleteFragment());
 54   }                                                57   }
 55 }                                                  58 }
 56                                                    59 
                                                   >>  60 
                                                   >>  61 
                                                   >>  62 // Initialization method
                                                   >>  63 
 57 void G4StatMFMicroCanonical::Initialize(const      64 void G4StatMFMicroCanonical::Initialize(const G4Fragment & theFragment) 
 58 {                                                  65 {
 59                                                    66   
 60   std::vector<G4StatMFMicroManager*>::iterator     67   std::vector<G4StatMFMicroManager*>::iterator it;
 61                                                    68 
 62   // Excitation Energy                             69   // Excitation Energy 
 63   G4double U = theFragment.GetExcitationEnergy     70   G4double U = theFragment.GetExcitationEnergy();
 64                                                    71 
 65   G4int A = theFragment.GetA_asInt();              72   G4int A = theFragment.GetA_asInt();
 66   G4int Z = theFragment.GetZ_asInt();              73   G4int Z = theFragment.GetZ_asInt();
 67   G4double x = 1.0 - 2.0*Z/G4double(A);            74   G4double x = 1.0 - 2.0*Z/G4double(A);
 68   G4Pow* g4calc = G4Pow::GetInstance();        <<  75   G4Pow* g4pow = G4Pow::GetInstance();
 69                                                    76     
 70   // Configuration temperature                     77   // Configuration temperature
 71   G4double TConfiguration = std::sqrt(8.0*U/G4     78   G4double TConfiguration = std::sqrt(8.0*U/G4double(A));
 72                                                    79   
 73   // Free internal energy at Temperature T = 0     80   // Free internal energy at Temperature T = 0
 74   __FreeInternalE0 = A*(                           81   __FreeInternalE0 = A*( 
 75       // Volume term (for T = 0)                   82       // Volume term (for T = 0)
 76       -G4StatMFParameters::GetE0() +               83       -G4StatMFParameters::GetE0() +  
 77       // Symmetry term                             84       // Symmetry term
 78       G4StatMFParameters::GetGamma0()*x*x          85       G4StatMFParameters::GetGamma0()*x*x 
 79       ) +                                          86       ) + 
 80     // Surface term (for T = 0)                    87     // Surface term (for T = 0)
 81     G4StatMFParameters::GetBeta0()*g4calc->Z23 <<  88     G4StatMFParameters::GetBeta0()*g4pow->Z23(A) + 
 82     // Coulomb term                                89     // Coulomb term 
 83     elm_coupling*0.6*Z*Z/(G4StatMFParameters:: <<  90     elm_coupling*(3.0/5.0)*Z*Z/(G4StatMFParameters::Getr0()*g4pow->Z13(A));
 84                                                    91   
 85   // Statistical weight                            92   // Statistical weight
 86   G4double W = 0.0;                                93   G4double W = 0.0;
 87                                                    94   
                                                   >>  95   
 88   // Mean breakup multiplicity                     96   // Mean breakup multiplicity
 89   __MeanMultiplicity = 0.0;                        97   __MeanMultiplicity = 0.0;
 90                                                    98   
 91   // Mean channel temperature                      99   // Mean channel temperature
 92   __MeanTemperature = 0.0;                        100   __MeanTemperature = 0.0;
 93                                                   101   
 94   // Mean channel entropy                         102   // Mean channel entropy
 95   __MeanEntropy = 0.0;                            103   __MeanEntropy = 0.0;
 96                                                   104   
 97   // Calculate entropy of compound nucleus        105   // Calculate entropy of compound nucleus
 98   G4double SCompoundNucleus = CalcEntropyOfCom    106   G4double SCompoundNucleus = CalcEntropyOfCompoundNucleus(theFragment,TConfiguration);
 99                                                   107   
100   // Statistical weight of compound nucleus       108   // Statistical weight of compound nucleus
101   _WCompoundNucleus = 1.0;                     << 109   _WCompoundNucleus = 1.0; // std::exp(SCompoundNucleus - SCompoundNucleus);
102                                                   110   
103   W += _WCompoundNucleus;                         111   W += _WCompoundNucleus;
104                                                << 112   
                                                   >> 113   
                                                   >> 114   
105   // Maximal fragment multiplicity allowed in     115   // Maximal fragment multiplicity allowed in direct simulation
106   G4int MaxMult = G4StatMFMicroCanonical::MaxA    116   G4int MaxMult = G4StatMFMicroCanonical::MaxAllowedMultiplicity;
107   if (A > 110) MaxMult -= 1;                      117   if (A > 110) MaxMult -= 1;
108                                                   118   
                                                   >> 119   
                                                   >> 120   
109   for (G4int im = 2; im <= MaxMult; im++) {       121   for (G4int im = 2; im <= MaxMult; im++) {
110     G4StatMFMicroManager * aMicroManager =        122     G4StatMFMicroManager * aMicroManager = 
111       new G4StatMFMicroManager(theFragment,im,    123       new G4StatMFMicroManager(theFragment,im,__FreeInternalE0,SCompoundNucleus);
112     _ThePartitionManagerVector.push_back(aMicr    124     _ThePartitionManagerVector.push_back(aMicroManager);
113   }                                               125   }
114                                                   126   
115   // W is the total probability                   127   // W is the total probability
116   W = std::accumulate(_ThePartitionManagerVect    128   W = std::accumulate(_ThePartitionManagerVector.begin(),
117           _ThePartitionManagerVector.end(),    << 129       _ThePartitionManagerVector.end(),
118           W, [](const G4double& running_total, << 130       W,SumProbabilities());
119                             G4StatMFMicroManag << 
120                          {                     << 
121                return running_total + manager- << 
122              } );                              << 
123                                                   131   
124   // Normalization of statistical weights         132   // Normalization of statistical weights
125   for (it =  _ThePartitionManagerVector.begin(    133   for (it =  _ThePartitionManagerVector.begin(); it !=  _ThePartitionManagerVector.end(); ++it) 
126     {                                             134     {
127       (*it)->Normalize(W);                        135       (*it)->Normalize(W);
128     }                                             136     }
129                                                   137 
130   _WCompoundNucleus /= W;                         138   _WCompoundNucleus /= W;
131                                                   139   
132   __MeanMultiplicity += 1.0 * _WCompoundNucleu    140   __MeanMultiplicity += 1.0 * _WCompoundNucleus;
133   __MeanTemperature += TConfiguration * _WComp    141   __MeanTemperature += TConfiguration * _WCompoundNucleus;
134   __MeanEntropy += SCompoundNucleus * _WCompou    142   __MeanEntropy += SCompoundNucleus * _WCompoundNucleus;
135                                                   143   
                                                   >> 144 
136   for (it =  _ThePartitionManagerVector.begin(    145   for (it =  _ThePartitionManagerVector.begin(); it !=  _ThePartitionManagerVector.end(); ++it) 
137     {                                             146     {
138       __MeanMultiplicity += (*it)->GetMeanMult    147       __MeanMultiplicity += (*it)->GetMeanMultiplicity();
139       __MeanTemperature += (*it)->GetMeanTempe    148       __MeanTemperature += (*it)->GetMeanTemperature();
140       __MeanEntropy += (*it)->GetMeanEntropy()    149       __MeanEntropy += (*it)->GetMeanEntropy();
141     }                                             150     }
142                                                   151   
143   return;                                         152   return;
144 }                                                 153 }
145                                                   154 
                                                   >> 155 
146 G4double G4StatMFMicroCanonical::CalcFreeInter    156 G4double G4StatMFMicroCanonical::CalcFreeInternalEnergy(const G4Fragment & theFragment, 
147               G4double T)                         157               G4double T)
148 {                                                 158 {
149   G4int A = theFragment.GetA_asInt();             159   G4int A = theFragment.GetA_asInt();
150   G4int Z = theFragment.GetZ_asInt();             160   G4int Z = theFragment.GetZ_asInt();
151   G4double A13 = G4Pow::GetInstance()->Z13(A);    161   G4double A13 = G4Pow::GetInstance()->Z13(A);
152                                                   162   
153   G4double InvLevelDensityPar = G4StatMFParame << 163   G4double InvLevelDensityPar = G4StatMFParameters::GetEpsilon0()*(1.0 + 3.0/G4double(A-1));
154     *(1.0 + 3.0/G4double(A-1));                << 
155                                                   164   
156   G4double VolumeTerm = (-G4StatMFParameters::    165   G4double VolumeTerm = (-G4StatMFParameters::GetE0()+T*T/InvLevelDensityPar)*A;
157                                                   166   
158   G4double SymmetryTerm = G4StatMFParameters:: << 167   G4double SymmetryTerm = G4StatMFParameters::GetGamma0()*(A - 2*Z)*(A - 2*Z)/G4double(A);
159     *(A - 2*Z)*(A - 2*Z)/G4double(A);          << 
160                                                   168   
161   G4double SurfaceTerm = (G4StatMFParameters:: << 169   G4double SurfaceTerm = (G4StatMFParameters::Beta(T)-T*G4StatMFParameters::DBetaDT(T))*A13*A13;
162         - T*G4StatMFParameters::DBetaDT(T))*A1 << 
163                                                   170   
164   G4double CoulombTerm = elm_coupling*0.6*Z*Z/ << 171   G4double CoulombTerm = elm_coupling*(3.0/5.0)*Z*Z/(G4StatMFParameters::Getr0()*A13);
165                                                   172   
166   return VolumeTerm + SymmetryTerm + SurfaceTe    173   return VolumeTerm + SymmetryTerm + SurfaceTerm + CoulombTerm;
167 }                                                 174 }
168                                                   175 
169 G4double                                          176 G4double 
170 G4StatMFMicroCanonical::CalcEntropyOfCompoundN    177 G4StatMFMicroCanonical::CalcEntropyOfCompoundNucleus(const G4Fragment & theFragment,
171                  G4double & TConf)                178                  G4double & TConf)
172   // Calculates Temperature and Entropy of com    179   // Calculates Temperature and Entropy of compound nucleus
173 {                                                 180 {
174   G4int A = theFragment.GetA_asInt();             181   G4int A = theFragment.GetA_asInt();
                                                   >> 182   //    const G4double Z = theFragment.GetZ();
175   G4double U = theFragment.GetExcitationEnergy    183   G4double U = theFragment.GetExcitationEnergy();
176   G4double A13 = G4Pow::GetInstance()->Z13(A);    184   G4double A13 = G4Pow::GetInstance()->Z13(A);
177                                                   185   
178   G4double Ta = std::max(std::sqrt(U/(0.125*A)    186   G4double Ta = std::max(std::sqrt(U/(0.125*A)),0.0012*MeV); 
179   G4double Tb = Ta;                               187   G4double Tb = Ta;
180                                                   188   
181   G4double ECompoundNucleus = CalcFreeInternal    189   G4double ECompoundNucleus = CalcFreeInternalEnergy(theFragment,Ta);
182   G4double Da = (U+__FreeInternalE0-ECompoundN    190   G4double Da = (U+__FreeInternalE0-ECompoundNucleus)/U;
183   G4double Db = 0.0;                              191   G4double Db = 0.0;
184                                                   192     
185   G4double InvLevelDensity = CalcInvLevelDensi << 193   G4double InvLevelDensity = CalcInvLevelDensity(static_cast<G4int>(A));
186                                                   194   
187   // bracketing the solution                      195   // bracketing the solution
188   if (Da == 0.0) {                                196   if (Da == 0.0) {
189     TConf = Ta;                                   197     TConf = Ta;
190     return 2*Ta*A/InvLevelDensity - G4StatMFPa    198     return 2*Ta*A/InvLevelDensity - G4StatMFParameters::DBetaDT(Ta)*A13*A13;
191   } else if (Da < 0.0) {                          199   } else if (Da < 0.0) {
192     do {                                          200     do {
193       Tb -= 0.5*Tb;                               201       Tb -= 0.5*Tb;
194       ECompoundNucleus = CalcFreeInternalEnerg    202       ECompoundNucleus = CalcFreeInternalEnergy(theFragment,Tb);
195       Db = (U+__FreeInternalE0-ECompoundNucleu    203       Db = (U+__FreeInternalE0-ECompoundNucleus)/U;
196     } while (Db < 0.0);                           204     } while (Db < 0.0);
197   } else {                                        205   } else {
198     do {                                          206     do {
199       Tb += 0.5*Tb;                               207       Tb += 0.5*Tb;
200       ECompoundNucleus = CalcFreeInternalEnerg    208       ECompoundNucleus = CalcFreeInternalEnergy(theFragment,Tb);
201       Db = (U+__FreeInternalE0-ECompoundNucleu    209       Db = (U+__FreeInternalE0-ECompoundNucleus)/U;
202     } while (Db > 0.0);                           210     } while (Db > 0.0);
203   }                                               211   }
204                                                   212   
205   G4double eps = 1.0e-14 * std::abs(Tb-Ta);    << 213   G4double eps = 1.0e-14 * std::fabs(Tb-Ta);
206                                                   214   
207   for (G4int i = 0; i < 1000; i++) {              215   for (G4int i = 0; i < 1000; i++) {
208     G4double Tc = (Ta+Tb)*0.5;                 << 216     G4double Tc = (Ta+Tb)/2.0;
209     if (std::abs(Ta-Tb) <= eps) {                 217     if (std::abs(Ta-Tb) <= eps) {
210       TConf = Tc;                                 218       TConf = Tc;
211       return 2*Tc*A/InvLevelDensity - G4StatMF    219       return 2*Tc*A/InvLevelDensity - G4StatMFParameters::DBetaDT(Tc)*A13*A13;
212     }                                             220     }
213     ECompoundNucleus = CalcFreeInternalEnergy(    221     ECompoundNucleus = CalcFreeInternalEnergy(theFragment,Tc);
214     G4double Dc = (U+__FreeInternalE0-ECompoun    222     G4double Dc = (U+__FreeInternalE0-ECompoundNucleus)/U;
215                                                   223     
216     if (Dc == 0.0) {                              224     if (Dc == 0.0) {
217       TConf = Tc;                                 225       TConf = Tc;
218       return 2*Tc*A/InvLevelDensity - G4StatMF    226       return 2*Tc*A/InvLevelDensity - G4StatMFParameters::DBetaDT(Tc)*A13*A13;
219     }                                             227     }
220                                                   228     
221     if (Da*Dc < 0.0) {                            229     if (Da*Dc < 0.0) {
222       Tb = Tc;                                    230       Tb = Tc;
223       Db = Dc;                                    231       Db = Dc;
224     } else {                                      232     } else {
225       Ta = Tc;                                    233       Ta = Tc;
226       Da = Dc;                                    234       Da = Dc;
227     }                                             235     } 
228   }                                               236   }
229                                                   237   
230   G4cout <<                                    << 238   G4cerr << 
231     "G4StatMFMicrocanoncal::CalcEntropyOfCompo    239     "G4StatMFMicrocanoncal::CalcEntropyOfCompoundNucleus: I can't calculate the temperature" 
232    << G4endl;                                     240    << G4endl;
233                                                   241   
234   return 0.0;                                     242   return 0.0;
235 }                                                 243 }
236                                                   244 
237 G4StatMFChannel *  G4StatMFMicroCanonical::Cho    245 G4StatMFChannel *  G4StatMFMicroCanonical::ChooseAandZ(const G4Fragment & theFragment)
238   // Choice of fragment atomic numbers and cha    246   // Choice of fragment atomic numbers and charges 
239 {                                                 247 {
240   // We choose a multiplicity (1,2,3,...) and  << 248     // We choose a multiplicity (1,2,3,...) and then a channel
241   G4double RandNumber = G4UniformRand();       << 249     G4double RandNumber = G4UniformRand();
242                                                   250 
243   if (RandNumber < _WCompoundNucleus) {        << 251     if (RandNumber < _WCompoundNucleus) { 
244                                                   252   
245     G4StatMFChannel * aChannel = new G4StatMFC << 253   G4StatMFChannel * aChannel = new G4StatMFChannel;
246     aChannel->CreateFragment(theFragment.GetA_ << 254   aChannel->CreateFragment(theFragment.GetA_asInt(),theFragment.GetZ_asInt());
247     return aChannel;                           << 255   return aChannel;
248                                                   256   
249   } else {                                     << 257     } else {
250                                                   258   
251     G4double AccumWeight = _WCompoundNucleus;  << 259   G4double AccumWeight = _WCompoundNucleus;
252     std::vector<G4StatMFMicroManager*>::iterat << 260   std::vector<G4StatMFMicroManager*>::iterator it;
253     for (it = _ThePartitionManagerVector.begin << 261   for (it = _ThePartitionManagerVector.begin(); it != _ThePartitionManagerVector.end(); ++it) {
254       AccumWeight += (*it)->GetProbability();  << 262       AccumWeight += (*it)->GetProbability();
255       if (RandNumber < AccumWeight) {          << 263       if (RandNumber < AccumWeight) {
256   return (*it)->ChooseChannel(theFragment.GetA << 264     return (*it)->ChooseChannel(theFragment.GetA(),theFragment.GetZ(),__MeanTemperature);
257       }                                        << 265       }
                                                   >> 266   }
                                                   >> 267   throw G4HadronicException(__FILE__, __LINE__, "G4StatMFMicroCanonical::ChooseAandZ: wrong normalization!");
258     }                                             268     }
259     throw G4HadronicException(__FILE__, __LINE << 
260   }                                            << 
261                                                   269 
262   return 0;                                    << 270     return 0; 
263 }                                                 271 }
264                                                   272 
265 G4double G4StatMFMicroCanonical::CalcInvLevelD    273 G4double G4StatMFMicroCanonical::CalcInvLevelDensity(G4int anA)
266 {                                                 274 {
267   G4double res = 0.0;                          << 275     // Calculate Inverse Density Level
268   if (anA > 1) {                               << 276     // Epsilon0*(1 + 3 /(Af - 1))
269     res = G4StatMFParameters::GetEpsilon0()*(1 << 277     if (anA == 1) return 0.0;
270   }                                            << 278     else return
271   return res;                                  << 279        G4StatMFParameters::GetEpsilon0()*(1.0+3.0/(anA - 1.0));
272 }                                                 280 }
273                                                   281