Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/de_excitation/multifragmentation/src/G4StatMFMacroCanonical.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/G4StatMFMacroCanonical.cc (Version 11.3.0) and /processes/hadronic/models/de_excitation/multifragmentation/src/G4StatMFMacroCanonical.cc (Version 6.2)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                    <<   3 // * DISCLAIMER                                                       *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th <<   5 // * The following disclaimer summarizes all the specific disclaimers *
  6 // * the Geant4 Collaboration.  It is provided <<   6 // * of contributors to this software. The specific disclaimers,which *
  7 // * conditions of the Geant4 Software License <<   7 // * govern, are listed with their locations in:                      *
  8 // * LICENSE and available at  http://cern.ch/ <<   8 // *   http://cern.ch/geant4/license                                  *
  9 // * include a list of copyright holders.      << 
 10 // *                                                9 // *                                                                  *
 11 // * Neither the authors of this software syst     10 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     11 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     12 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     13 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file  <<  14 // * use.                                                             *
 16 // * for the full disclaimer and the limitatio << 
 17 // *                                               15 // *                                                                  *
 18 // * This  code  implementation is the result  <<  16 // * This  code  implementation is the  intellectual property  of the *
 19 // * technical work of the GEANT4 collaboratio <<  17 // * GEANT4 collaboration.                                            *
 20 // * By using,  copying,  modifying or  distri <<  18 // * By copying,  distributing  or modifying the Program (or any work *
 21 // * any work based  on the software)  you  ag <<  19 // * based  on  the Program)  you indicate  your  acceptance of  this *
 22 // * use  in  resulting  scientific  publicati <<  20 // * statement, and all its terms.                                    *
 23 // * acceptance of all terms of the Geant4 Sof << 
 24 // *******************************************     21 // ********************************************************************
 25 //                                                 22 //
 26 //                                                 23 //
                                                   >>  24 // $Id: G4StatMFMacroCanonical.cc,v 1.3 2003/11/06 18:11:44 lara Exp $
                                                   >>  25 // GEANT4 tag $Name: geant4-06-00-patch-01 $
 27 //                                                 26 //
 28 // by V. Lara                                      27 // by V. Lara
 29 // -------------------------------------------     28 // --------------------------------------------------------------------
 30 //                                             << 
 31 // Modified:                                   << 
 32 // 25.07.08 I.Pshenichnov (in collaboration wi << 
 33 //          Mishustin (FIAS, Frankfurt, INR, M << 
 34 //          Moscow, pshenich@fias.uni-frankfur << 
 35 //          a fagment with Z=A; fixed memory l << 
 36                                                    29 
 37 #include "G4StatMFMacroCanonical.hh"               30 #include "G4StatMFMacroCanonical.hh"
 38 #include "G4PhysicalConstants.hh"              <<  31 
 39 #include "G4SystemOfUnits.hh"                  << 
 40 #include "G4Pow.hh"                            << 
 41                                                    32 
 42 // constructor                                     33 // constructor
 43 G4StatMFMacroCanonical::G4StatMFMacroCanonical     34 G4StatMFMacroCanonical::G4StatMFMacroCanonical(const G4Fragment & theFragment) 
 44 {                                                  35 {
 45                                                    36 
 46   // Get memory for clusters                       37   // Get memory for clusters
 47   _theClusters.push_back(new G4StatMFMacroNucl     38   _theClusters.push_back(new G4StatMFMacroNucleon);              // Size 1
 48   _theClusters.push_back(new G4StatMFMacroBiNu     39   _theClusters.push_back(new G4StatMFMacroBiNucleon);            // Size 2
 49   _theClusters.push_back(new G4StatMFMacroTriN     40   _theClusters.push_back(new G4StatMFMacroTriNucleon);           // Size 3
 50   _theClusters.push_back(new G4StatMFMacroTetr     41   _theClusters.push_back(new G4StatMFMacroTetraNucleon);         // Size 4
 51   for (G4int i = 4; i < theFragment.GetA_asInt <<  42   for (G4int i = 4; i < theFragment.GetA(); i++)   
 52     _theClusters.push_back(new G4StatMFMacroMu     43     _theClusters.push_back(new G4StatMFMacroMultiNucleon(i+1)); // Size 5 ... A
 53                                                    44   
 54   // Perform class initialization                  45   // Perform class initialization
 55   Initialize(theFragment);                         46   Initialize(theFragment);
 56                                                    47     
 57 }                                                  48 }
 58                                                    49 
                                                   >>  50 
 59 // destructor                                      51 // destructor
 60 G4StatMFMacroCanonical::~G4StatMFMacroCanonica     52 G4StatMFMacroCanonical::~G4StatMFMacroCanonical() 
 61 {                                                  53 {
 62   // garbage collection                            54   // garbage collection
 63   if (!_theClusters.empty())                       55   if (!_theClusters.empty()) 
 64     {                                              56     {
 65       std::for_each(_theClusters.begin(),_theC     57       std::for_each(_theClusters.begin(),_theClusters.end(),DeleteFragment());
 66     }                                              58     }
 67 }                                                  59 }
 68                                                    60 
                                                   >>  61 // operators definitions
                                                   >>  62 G4StatMFMacroCanonical & 
                                                   >>  63 G4StatMFMacroCanonical::operator=(const G4StatMFMacroCanonical & ) 
                                                   >>  64 {
                                                   >>  65   throw G4HadronicException(__FILE__, __LINE__, "G4StatMFMacroCanonical::operator= meant to not be accessable");
                                                   >>  66   return *this;
                                                   >>  67 }
                                                   >>  68 
                                                   >>  69 G4bool G4StatMFMacroCanonical::operator==(const G4StatMFMacroCanonical & ) const 
                                                   >>  70 {
                                                   >>  71   return false;
                                                   >>  72 }
                                                   >>  73 
                                                   >>  74 
                                                   >>  75 G4bool G4StatMFMacroCanonical::operator!=(const G4StatMFMacroCanonical & ) const 
                                                   >>  76 {
                                                   >>  77   return true;
                                                   >>  78 }
                                                   >>  79 
                                                   >>  80 
 69 // Initialization method                           81 // Initialization method
                                                   >>  82 
                                                   >>  83 
 70 void G4StatMFMacroCanonical::Initialize(const      84 void G4StatMFMacroCanonical::Initialize(const G4Fragment & theFragment) 
 71 {                                                  85 {
 72                                                    86   
 73   G4int A = theFragment.GetA_asInt();          <<  87   G4double A = theFragment.GetA();
 74   G4int Z = theFragment.GetZ_asInt();          <<  88   G4double Z = theFragment.GetZ();
 75   G4double x = 1.0 - 2.0*Z/G4double(A);        << 
 76   G4Pow* g4calc = G4Pow::GetInstance();        << 
 77                                                    89   
 78   // Free Internal energy at T = 0                 90   // Free Internal energy at T = 0
 79   __FreeInternalE0 = A*( -G4StatMFParameters:: <<  91   __FreeInternalE0 = A*( -G4StatMFParameters::GetE0() +          // Volume term (for T = 0)
 80        G4StatMFParameters::GetGamma0()*x*x) // <<  92        G4StatMFParameters::GetGamma0()*        // Symmetry term
 81     + G4StatMFParameters::GetBeta0()*g4calc->Z <<  93        (1.0-2.0*Z/A)*(1.0-2.0*Z/A) ) +
 82     0.6*elm_coupling*Z*Z/(G4StatMFParameters:: <<  94     G4StatMFParameters::GetBeta0()*pow(A,2.0/3.0) +              // Surface term (for T = 0)
 83         g4calc->Z13(A));                       <<  95     (3.0/5.0)*elm_coupling*Z*Z/(G4StatMFParameters::Getr0()*     // Coulomb term 
                                                   >>  96         pow(A,1.0/3.0));
                                                   >>  97   
                                                   >>  98   
 84                                                    99   
 85   CalculateTemperature(theFragment);              100   CalculateTemperature(theFragment);
                                                   >> 101   
 86   return;                                         102   return;
 87 }                                                 103 }
 88                                                   104 
                                                   >> 105 
                                                   >> 106 
                                                   >> 107 
                                                   >> 108 
 89 void G4StatMFMacroCanonical::CalculateTemperat    109 void G4StatMFMacroCanonical::CalculateTemperature(const G4Fragment & theFragment)
 90 {                                                 110 {
 91   // Excitation Energy                            111   // Excitation Energy
 92   G4double U = theFragment.GetExcitationEnergy    112   G4double U = theFragment.GetExcitationEnergy();
 93                                                   113   
 94   G4int A = theFragment.GetA_asInt();          << 114   G4double A = theFragment.GetA();
 95   G4int Z = theFragment.GetZ_asInt();          << 115   G4double Z = theFragment.GetZ();
 96                                                   116   
 97   // Fragment Multiplicity                        117   // Fragment Multiplicity
 98   G4double FragMult = std::max((1.0+(2.31/MeV)    118   G4double FragMult = std::max((1.0+(2.31/MeV)*(U/A - 3.5*MeV))*A/100.0, 2.0);
 99                                                   119 
                                                   >> 120 
100   // Parameter Kappa                              121   // Parameter Kappa
101   G4Pow* g4calc = G4Pow::GetInstance();        << 122   _Kappa = (1.0+elm_coupling*(pow(FragMult,1./3.)-1)/
102   _Kappa = (1.0+elm_coupling*(g4calc->A13(Frag << 123       (G4StatMFParameters::Getr0()*pow(A,1./3.)));
103       (G4StatMFParameters::Getr0()*g4calc->Z13 << 
104   _Kappa = _Kappa*_Kappa*_Kappa - 1.0;            124   _Kappa = _Kappa*_Kappa*_Kappa - 1.0;
                                                   >> 125 
105                                                   126   
106   G4StatMFMacroTemperature * theTemp = new        127   G4StatMFMacroTemperature * theTemp = new  
107     G4StatMFMacroTemperature(A,Z,U,__FreeInter    128     G4StatMFMacroTemperature(A,Z,U,__FreeInternalE0,_Kappa,&_theClusters);
108                                                   129   
109   __MeanTemperature = theTemp->CalcTemperature    130   __MeanTemperature = theTemp->CalcTemperature();
110   _ChemPotentialNu = theTemp->GetChemicalPoten    131   _ChemPotentialNu = theTemp->GetChemicalPotentialNu();
111   _ChemPotentialMu = theTemp->GetChemicalPoten    132   _ChemPotentialMu = theTemp->GetChemicalPotentialMu();
112   __MeanMultiplicity = theTemp->GetMeanMultipl    133   __MeanMultiplicity = theTemp->GetMeanMultiplicity();
113   __MeanEntropy = theTemp->GetEntropy();          134   __MeanEntropy = theTemp->GetEntropy();
114                                                   135   
115   delete theTemp;                                 136   delete theTemp;     
116                                                   137   
117   return;                                         138   return;
118 }                                                 139 }
119                                                   140 
                                                   >> 141 
120 // -------------------------------------------    142 // --------------------------------------------------------------------------
121                                                   143 
122 G4StatMFChannel * G4StatMFMacroCanonical::Choo    144 G4StatMFChannel * G4StatMFMacroCanonical::ChooseAandZ(const G4Fragment &theFragment)
123 // Calculate total fragments multiplicity, fra << 145     // Calculate total fragments multiplicity, fragment atomic numbers and charges
124 {                                                 146 {
125   G4int A = theFragment.GetA_asInt();          << 147   G4double A = theFragment.GetA();
126   G4int Z = theFragment.GetZ_asInt();          << 148   G4double Z = theFragment.GetZ();
127                                                   149   
128   std::vector<G4int> ANumbers(A);              << 150   std::vector<G4double> ANumbers(static_cast<G4int>(A));
129                                                   151   
130   G4double Multiplicity = ChooseA(A,ANumbers);    152   G4double Multiplicity = ChooseA(A,ANumbers);
131                                                   153   
132   std::vector<G4int> FragmentsA;               << 
133                                                   154   
134   G4int i = 0;                                 << 155   std::vector<G4double> FragmentsA;
135   for (i = 0; i < A; i++)                      << 
136     {                                          << 
137       for (G4int j = 0; j < ANumbers[i]; j++)  << 
138     }                                          << 
139                                                   156   
140   // Sort fragments in decreasing order        << 157   G4int i = 0;  
141   G4int im = 0;                                << 158     for (i = 0; i < A; i++) 
142   for (G4int j = 0; j < Multiplicity; j++)     << 159       {
143     {                                          << 160   for (G4int j = 0; j < ANumbers[i]; j++) FragmentsA.push_back(i+1);
144       G4int FragmentsAMax = 0;                 << 161       }
145       im = j;                                  << 162 
146       for (i = j; i < Multiplicity; i++)       << 163     
147   {                                            << 164     // Sort fragments in decreasing order
148     if (FragmentsA[i] <= FragmentsAMax) { cont << 165     G4int im = 0;
149     else                                       << 166     for (G4int j = 0; j < Multiplicity; j++) 
150       {                                        << 167       {
151         im = i;                                << 168   G4double FragmentsAMax = 0.0;
152         FragmentsAMax = FragmentsA[im];        << 169   im = j;
153       }                                        << 170   for (i = j; i < Multiplicity; i++) 
154   }                                            << 171     {
155       if (im != j)                             << 172       if (FragmentsA[i] <= FragmentsAMax) continue;
156   {                                            << 173       else 
157     FragmentsA[im] = FragmentsA[j];            << 174         {
158     FragmentsA[j]  = FragmentsAMax;            << 175     im = i;
159   }                                            << 176     FragmentsAMax = FragmentsA[im];
160     }                                          << 177         }
161   return ChooseZ(Z,FragmentsA);                << 178     }
                                                   >> 179   
                                                   >> 180   if (im != j) 
                                                   >> 181     {
                                                   >> 182       FragmentsA[im] = FragmentsA[j];
                                                   >> 183       FragmentsA[j] = FragmentsAMax;
                                                   >> 184     }
                                                   >> 185       }
                                                   >> 186     
                                                   >> 187     return ChooseZ(static_cast<G4int>(Z),FragmentsA);
162 }                                                 188 }
163                                                   189 
164 G4double G4StatMFMacroCanonical::ChooseA(G4int << 190 
                                                   >> 191 
                                                   >> 192 G4double G4StatMFMacroCanonical::ChooseA(const G4double A, std::vector<G4double> & ANumbers)
165   // Determines fragments multiplicities and c    193   // Determines fragments multiplicities and compute total fragment multiplicity
166 {                                                 194 {
167   G4double multiplicity = 0.0;                    195   G4double multiplicity = 0.0;
168   G4int i;                                        196   G4int i;
169                                                << 197   
                                                   >> 198   
170   std::vector<G4double> AcumMultiplicity;         199   std::vector<G4double> AcumMultiplicity;
171   AcumMultiplicity.reserve(A);                 << 200   AcumMultiplicity.reserve(static_cast<G4int>(A));
172                                                   201 
173   AcumMultiplicity.push_back((*(_theClusters.b    202   AcumMultiplicity.push_back((*(_theClusters.begin()))->GetMeanMultiplicity());
174   for (std::vector<G4VStatMFMacroCluster*>::it    203   for (std::vector<G4VStatMFMacroCluster*>::iterator it = _theClusters.begin()+1;
175        it != _theClusters.end(); ++it)            204        it != _theClusters.end(); ++it)
176     {                                             205     {
177       AcumMultiplicity.push_back((*it)->GetMea    206       AcumMultiplicity.push_back((*it)->GetMeanMultiplicity()+AcumMultiplicity.back());
178     }                                             207     }
179                                                   208   
180   G4int CheckA;                                   209   G4int CheckA;
181   do {                                            210   do {
182     CheckA = -1;                                  211     CheckA = -1;
183     G4int SumA = 0;                               212     G4int SumA = 0;
184     G4int ThisOne = 0;                            213     G4int ThisOne = 0;
185     multiplicity = 0.0;                           214     multiplicity = 0.0;
186     for (i = 0; i < A; i++) ANumbers[i] = 0;   << 215     for (i = 0; i < A; i++) ANumbers[i] = 0.0;
187     do {                                          216     do {
188       G4double RandNumber = G4UniformRand()*__    217       G4double RandNumber = G4UniformRand()*__MeanMultiplicity;
189       for (i = 0; i < A; i++) {                   218       for (i = 0; i < A; i++) {
190   if (RandNumber < AcumMultiplicity[i]) {         219   if (RandNumber < AcumMultiplicity[i]) {
191     ThisOne = i;                                  220     ThisOne = i;
192     break;                                        221     break;
193   }                                               222   }
194       }                                           223       }
195       multiplicity++;                             224       multiplicity++;
196       ANumbers[ThisOne] = ANumbers[ThisOne]+1;    225       ANumbers[ThisOne] = ANumbers[ThisOne]+1;
197       SumA += ThisOne+1;                          226       SumA += ThisOne+1;
198       CheckA = A - SumA;                       << 227       CheckA = static_cast<G4int>(A) - SumA;
199                                                   228       
200       // Loop checking, 05-Aug-2015, Vladimir  << 
201     } while (CheckA > 0);                         229     } while (CheckA > 0);
202                                                   230     
203     // Loop checking, 05-Aug-2015, Vladimir Iv << 231   } while (CheckA < 0 || abs(__MeanMultiplicity - multiplicity) > sqrt(__MeanMultiplicity) + 1./2.);
204   } while (CheckA < 0 || std::abs(__MeanMultip << 
205                                                   232   
206   return multiplicity;                            233   return multiplicity;
207 }                                                 234 }
208                                                   235 
209 G4StatMFChannel * G4StatMFMacroCanonical::Choo << 236 
210               std::vector<G4int> & FragmentsA) << 237 G4StatMFChannel * G4StatMFMacroCanonical::ChooseZ(const G4int & Z, 
                                                   >> 238               std::vector<G4double> & FragmentsA)
211     //                                            239     // 
212 {                                                 240 {
213   G4Pow* g4calc = G4Pow::GetInstance();        << 241   std::vector<G4double> FragmentsZ;
214   std::vector<G4int> FragmentsZ;               << 
215                                                   242   
216   G4int DeltaZ = 0;                            << 243   G4double DeltaZ = 0.0;
217   G4double CP =  G4StatMFParameters::GetCoulom << 244   G4double CP = (3./5.)*(elm_coupling/G4StatMFParameters::Getr0())*
218   G4int multiplicity = (G4int)FragmentsA.size( << 245     (1.0 - 1.0/pow(1.0+G4StatMFParameters::GetKappaCoulomb(),1./3.));
219                                                   246   
220   do {                                         << 247   G4int multiplicity = FragmentsA.size();
221     FragmentsZ.clear();                        << 248   
222     G4int SumZ = 0;                            << 249   do 
223     for (G4int i = 0; i < multiplicity; ++i)   << 250     {
224       {                                        << 251       FragmentsZ.clear();
225   G4int A = FragmentsA[i];                     << 252       G4int SumZ = 0;
226   if (A <= 1)                                  << 253       for (G4int i = 0; i < multiplicity; i++) 
227     {                                          << 254   {
228       G4double RandNumber = G4UniformRand();   << 255     G4double A = FragmentsA[i];
229       if (RandNumber < (*_theClusters.begin()) << 256     if (A <= 1.0) 
230         {                                      << 257       {
231     FragmentsZ.push_back(1);                   << 258         G4double RandNumber = G4UniformRand();
232     SumZ += FragmentsZ[i];                     << 259         if (RandNumber < (*_theClusters.begin())->GetZARatio()) 
233         }                                      << 260     {
234       else FragmentsZ.push_back(0);            << 261       FragmentsZ.push_back(1.0);
235     }                                          << 262       SumZ += static_cast<G4int>(FragmentsZ[i]);
236   else                                         << 263     } 
237     {                                          << 264         else FragmentsZ.push_back(0.0);
238       G4double RandZ;                          << 265       } 
239       G4double CC = 8.0*G4StatMFParameters::Ge << 266     else 
240         + 2*CP*g4calc->Z23(FragmentsA[i]);     << 267       {
241       G4double ZMean;                          << 268         G4double RandZ;
242       if (FragmentsA[i] > 1 && FragmentsA[i] < << 269         G4double CC = 8.0*G4StatMFParameters::GetGamma0()+2.0*CP*pow(FragmentsA[i],2./3.);
243       else {                                   << 270         G4double ZMean;
244         ZMean = FragmentsA[i]*(4.0*G4StatMFPar << 271         if (FragmentsA[i] > 1.5 && FragmentsA[i] < 4.5) ZMean = 0.5*FragmentsA[i];
245              + _ChemPotentialNu)/CC;           << 272         else ZMean = FragmentsA[i]*(4.0*G4StatMFParameters::GetGamma0()+_ChemPotentialNu)/CC;
                                                   >> 273         G4double ZDispersion = sqrt(FragmentsA[i]*__MeanTemperature/CC);
                                                   >> 274         G4int z;
                                                   >> 275         do 
                                                   >> 276     {
                                                   >> 277       RandZ = G4RandGauss::shoot(ZMean,ZDispersion);
                                                   >> 278       z = static_cast<G4int>(RandZ+0.5);
                                                   >> 279     } while (z < 0 || z > A);
                                                   >> 280         FragmentsZ.push_back(z);
                                                   >> 281         SumZ += z;
246       }                                           282       }
247       G4double ZDispersion = std::sqrt(Fragmen << 283   }
248       G4int z;                                 << 284       DeltaZ = Z - SumZ;
249       do                                       << 285     }
250         {                                      << 286   while (abs(DeltaZ) > 1.1);
251     RandZ = G4RandGauss::shoot(ZMean,ZDispersi << 
252     z = G4lrint(RandZ+0.5);                    << 
253     // Loop checking, 05-Aug-2015, Vladimir Iv << 
254         } while (z < 0 || z > A);              << 
255       FragmentsZ.push_back(z);                 << 
256       SumZ += z;                               << 
257     }                                          << 
258       }                                        << 
259     DeltaZ = Z - SumZ;                         << 
260   // Loop checking, 05-Aug-2015, Vladimir Ivan << 
261   } while (std::abs(DeltaZ) > 1);              << 
262                                                   287     
263   // DeltaZ can be 0, 1 or -1                     288   // DeltaZ can be 0, 1 or -1
264   G4int idx = 0;                                  289   G4int idx = 0;
265   if (DeltaZ < 0.0)                               290   if (DeltaZ < 0.0)
266     {                                             291     {
267       while (FragmentsZ[idx] < 1) { ++idx; }   << 292       while (FragmentsZ[idx] < 0.5) ++idx;
268     }                                             293     }
269   FragmentsZ[idx] += DeltaZ;                      294   FragmentsZ[idx] += DeltaZ;
270                                                   295   
271   G4StatMFChannel * theChannel = new G4StatMFC    296   G4StatMFChannel * theChannel = new G4StatMFChannel;
272   for (G4int i = multiplicity-1; i >= 0; --i)  << 297   for (G4int i = multiplicity-1; i >= 0; i--) 
273     {                                             298     {
274       theChannel->CreateFragment(FragmentsA[i]    299       theChannel->CreateFragment(FragmentsA[i],FragmentsZ[i]);
275     }                                             300     }
276                                                << 301   
277   return theChannel;                           << 302     
                                                   >> 303     return theChannel;
278 }                                                 304 }
279                                                   305