G4StatMFMacroMultiplicity.cc

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration.                      *
00020 // * By using,  copying,  modifying or  distributing the software (or *
00021 // * any work based  on the software)  you  agree  to acknowledge its *
00022 // * use  in  resulting  scientific  publications,  and indicate your *
00023 // * acceptance of all terms of the Geant4 Software license.          *
00024 // ********************************************************************
00025 //
00026 //
00027 // $Id$
00028 //
00029 // Hadronic Process: Nuclear De-excitations
00030 // by V. Lara
00031 //
00032 // Modified:
00033 // 25.07.08 I.Pshenichnov (in collaboration with Alexander Botvina and Igor 
00034 //          Mishustin (FIAS, Frankfurt, INR, Moscow and Kurchatov Institute, 
00035 //          Moscow, pshenich@fias.uni-frankfurt.de) additional checks in
00036 //          solver of equation for the chemical potential
00037 
00038 #include "G4StatMFMacroMultiplicity.hh"
00039 #include "G4PhysicalConstants.hh"
00040 
00041 // operators definitions
00042 G4StatMFMacroMultiplicity & 
00043 G4StatMFMacroMultiplicity::operator=(const G4StatMFMacroMultiplicity & ) 
00044 {
00045     throw G4HadronicException(__FILE__, __LINE__, "G4StatMFMacroMultiplicity::operator= meant to not be accessable");
00046     return *this;
00047 }
00048 
00049 G4bool G4StatMFMacroMultiplicity::operator==(const G4StatMFMacroMultiplicity & ) const 
00050 {
00051     throw G4HadronicException(__FILE__, __LINE__, "G4StatMFMacroMultiplicity::operator== meant to not be accessable");
00052     return false;
00053 }
00054 
00055 
00056 G4bool G4StatMFMacroMultiplicity::operator!=(const G4StatMFMacroMultiplicity & ) const 
00057 {
00058     throw G4HadronicException(__FILE__, __LINE__, "G4StatMFMacroMultiplicity::operator!= meant to not be accessable");
00059     return true;
00060 }
00061 
00062 
00063 
00064 
00065 G4double G4StatMFMacroMultiplicity::CalcChemicalPotentialMu(void) 
00066     //  Calculate Chemical potential \mu
00067     // For that is necesary to calculate mean multiplicities
00068 {
00069     G4double CP = ((3./5.)*elm_coupling/G4StatMFParameters::Getr0())*
00070         (1.0-1.0/std::pow(1.0+G4StatMFParameters::GetKappaCoulomb(),1.0/3.0));
00071 
00072     // starting value for chemical potential \mu
00073     // it is the derivative of F(T,V)-\nu*Z w.r.t. Af in Af=5
00074     G4double ZA5 = _theClusters->operator[](4)->GetZARatio();
00075     G4double ILD5 = _theClusters->operator[](4)->GetInvLevelDensity();
00076     _ChemPotentialMu = -G4StatMFParameters::GetE0()-
00077         _MeanTemperature*_MeanTemperature/ILD5 -
00078         _ChemPotentialNu*ZA5 + 
00079         G4StatMFParameters::GetGamma0()*(1.0-2.0*ZA5)*(1.0-2.0*ZA5) +
00080         (2.0/3.0)*G4StatMFParameters::Beta(_MeanTemperature)/std::pow(5.,1./3.) +
00081         (5.0/3.0)*CP*ZA5*ZA5*std::pow(5.,2./3.) -
00082         1.5*_MeanTemperature/5.0;
00083                 
00084 
00085 
00086     G4double ChemPa = _ChemPotentialMu;
00087     if (ChemPa/_MeanTemperature > 10.0) ChemPa = 10.0*_MeanTemperature;
00088     G4double ChemPb = ChemPa - 0.5*std::abs(ChemPa);
00089     
00090     
00091     G4double fChemPa = this->operator()(ChemPa); 
00092     G4double fChemPb = this->operator()(ChemPb); 
00093     
00094 
00095     // Set the precision level for locating the root. 
00096     // If the root is inside this interval, then it's done! 
00097     G4double intervalWidth = 1.e-4;
00098 
00099     // bracketing the solution
00100     G4int iterations = 0;
00101     while (fChemPa*fChemPb > 0.0 && iterations < 100) 
00102     {
00103         if (std::abs(fChemPa) <= std::abs(fChemPb)) 
00104         {
00105             ChemPa += 0.6*(ChemPa-ChemPb);
00106             fChemPa = this->operator()(ChemPa);
00107             iterations++;
00108         } 
00109         else 
00110         {
00111             ChemPb += 0.6*(ChemPb-ChemPa);
00112             fChemPb = this->operator()(ChemPb);
00113             iterations++;
00114         }
00115     }
00116 
00117     if (fChemPa*fChemPb > 0.0) // the bracketing failed, complain 
00118     {
00119       G4cerr <<"G4StatMFMacroMultiplicity:"<<" ChemPa="<<ChemPa<<" ChemPb="<<ChemPb<< G4endl;
00120       G4cerr <<"G4StatMFMacroMultiplicity:"<<" fChemPa="<<fChemPa<<" fChemPb="<<fChemPb<< G4endl;
00121         throw G4HadronicException(__FILE__, __LINE__, "G4StatMFMacroMultiplicity::CalcChemicalPotentialMu: I couldn't bracket the root.");
00122     }
00123     else if (fChemPa*fChemPb < 0.0 && std::abs(ChemPa-ChemPb) > intervalWidth) // the bracketing was OK, try to locate the root
00124     {   
00125     G4Solver<G4StatMFMacroMultiplicity> * theSolver = new G4Solver<G4StatMFMacroMultiplicity>(100,intervalWidth);
00126     theSolver->SetIntervalLimits(ChemPa,ChemPb);
00127     //    if (!theSolver->Crenshaw(*this)) 
00128     if (!theSolver->Brent(*this)) 
00129     {
00130       G4cerr <<"G4StatMFMacroMultiplicity:"<<" ChemPa="<<ChemPa<<" ChemPb="<<ChemPb<< G4endl;
00131       G4cerr <<"G4StatMFMacroMultiplicity:"<<" fChemPa="<<fChemPa<<" fChemPb="<<fChemPb<< G4endl;
00132         throw G4HadronicException(__FILE__, __LINE__, "G4StatMFMacroMultiplicity::CalcChemicalPotentialMu: I couldn't find the root.");
00133     }
00134     _ChemPotentialMu = theSolver->GetRoot();
00135     delete theSolver;
00136     }
00137     else // the root is within the interval, which is shorter then the precision level - all done 
00138     {
00139      _ChemPotentialMu = ChemPa;
00140     }
00141 
00142     return _ChemPotentialMu;
00143 }
00144 
00145 
00146 
00147 G4double G4StatMFMacroMultiplicity::CalcMeanA(const G4double mu)
00148 {
00149   G4double r03 = G4StatMFParameters::Getr0(); r03 *= r03*r03;
00150   G4double V0 = (4.0/3.0)*pi*theA*r03;
00151 
00152   G4double MeanA = 0.0;
00153         
00154   _MeanMultiplicity = 0.0;
00155         
00156  
00157   G4int n = 1;
00158  for (std::vector<G4VStatMFMacroCluster*>::iterator i = _theClusters->begin(); 
00159       i != _theClusters->end(); ++i) 
00160    {
00161      G4double multip = (*i)->CalcMeanMultiplicity(V0*_Kappa,mu,_ChemPotentialNu,_MeanTemperature);
00162      MeanA += multip*static_cast<G4double>(n++);
00163      _MeanMultiplicity += multip;
00164    }
00165 
00166   return MeanA;
00167 }

Generated on Mon May 27 17:49:54 2013 for Geant4 by  doxygen 1.4.7