G4Solver< Function > Class Template Reference

#include <G4Solver.hh>


Public Types

 DefaultMaxIter = 100
enum  { DefaultMaxIter = 100 }

Public Member Functions

 G4Solver ()
 G4Solver (const G4int iterations, const G4double tol)
 G4Solver (const G4Solver &right)
 ~G4Solver ()
G4Solveroperator= (const G4Solver &right)
G4bool operator== (const G4Solver &right) const
G4bool operator!= (const G4Solver &right) const
G4int GetMaxIterations (void) const
void SetMaxIterations (const G4int iterations)
G4double GetTolerance (void) const
void SetTolerance (const G4double epsilon)
G4double GetIntervalLowerLimit (void) const
G4double GetIntervalUpperLimit (void) const
void SetIntervalLimits (const G4double Limit1, const G4double Limit2)
G4double GetRoot (void) const
G4bool Bisection (Function &theFunction)
G4bool RegulaFalsi (Function &theFunction)
G4bool Brent (Function &theFunction)
G4bool Crenshaw (Function &theFunction)


Detailed Description

template<class Function>
class G4Solver< Function >

Definition at line 41 of file G4Solver.hh.


Member Enumeration Documentation

template<class Function>
anonymous enum

Enumerator:
DefaultMaxIter 

Definition at line 44 of file G4Solver.hh.

00044 {DefaultMaxIter = 100};


Constructor & Destructor Documentation

template<class Function>
G4Solver< Function >::G4Solver (  )  [inline]

Definition at line 47 of file G4Solver.hh.

00047                : MaxIter(DefaultMaxIter), tolerance(DefaultTolerance),
00048                  a(0.0), b(0.0), root(0.0) {};

template<class Function>
G4Solver< Function >::G4Solver ( const G4int  iterations,
const G4double  tol 
) [inline]

Definition at line 50 of file G4Solver.hh.

00050                                                          :
00051         MaxIter(iterations), tolerance(tol),
00052         a(0.0), b(0.0), root(0.0) {};

template<class Function>
G4Solver< Function >::G4Solver ( const G4Solver< Function > &  right  ) 

Definition at line 36 of file G4Solver.cc.

References G4Solver< Function >::a, G4Solver< Function >::b, G4Solver< Function >::MaxIter, G4Solver< Function >::root, and G4Solver< Function >::tolerance.

00037 {
00038     MaxIter = right.MaxIter;
00039     tolerance = right.tolerance;
00040     a = right.a;
00041     b = right.b;
00042     root = right.root;
00043 }

template<class Function>
G4Solver< Function >::~G4Solver (  )  [inline]

Definition at line 58 of file G4Solver.hh.

00058 {};


Member Function Documentation

template<class Function>
G4bool G4Solver< Function >::Bisection ( Function &  theFunction  ) 

Definition at line 27 of file G4Solver.icc.

References G4cerr, and G4endl.

00028 {
00029     // Check the interval before start
00030     if (a > b || std::abs(a-b) <= tolerance) 
00031     {
00032         G4cerr << "G4Solver::Bisection: The interval must be properly set." << G4endl;
00033         return false;
00034     }
00035     G4double fa = theFunction(a);
00036     G4double fb = theFunction(b);
00037     if (fa*fb > 0.0) 
00038     {
00039         G4cerr << "G4Solver::Bisection: The interval must include a root." << G4endl;
00040         return false;
00041     }
00042     
00043     G4double eps=tolerance*(b-a);
00044     
00045     
00046     // Finding the root
00047     for (G4int i = 0; i < MaxIter; i++) 
00048     {
00049         G4double c = (a+b)/2.0;
00050         if ((b-a) < eps) 
00051         {
00052             root = c;
00053             return true;
00054         }
00055         G4double fc = theFunction(c);
00056         if (fc == 0.0) 
00057         {
00058             root = c;
00059             return true;
00060         }
00061         if (fa*fc < 0.0) 
00062         {
00063             a=c;
00064             fa=fc;
00065         } 
00066         else 
00067         {
00068             b=c;
00069             fb=fc;
00070         }
00071     }
00072     G4cerr << "G4Solver::Bisection: Excedded maximum number of iterations whithout convegence." << G4endl;
00073     return false;
00074 }

template<class Function>
G4bool G4Solver< Function >::Brent ( Function &  theFunction  ) 

Definition at line 130 of file G4Solver.icc.

References G4cerr, and G4endl.

Referenced by G4StatMFMacroMultiplicity::CalcChemicalPotentialMu(), G4StatMFMacroChemicalPotential::CalcChemicalPotentialNu(), and G4StatMFMacroTemperature::CalcTemperature().

00131 {
00132     
00133     const G4double precision = 3.0e-8;
00134     
00135     // Check the interval before start
00136     if (a > b || std::abs(a-b) <= tolerance) 
00137     {
00138         G4cerr << "G4Solver::Brent: The interval must be properly set." << G4endl;
00139         return false;
00140     }
00141     G4double fa = theFunction(a);
00142     G4double fb = theFunction(b);
00143     if (fa*fb > 0.0) 
00144     {
00145         G4cerr << "G4Solver::Brent: The interval must include a root." << G4endl;
00146         return false;
00147     }
00148     
00149     G4double c = b;
00150     G4double fc = fb;
00151     G4double d = 0.0;
00152     G4double e = 0.0;
00153     
00154     for (G4int i=0; i < MaxIter; i++) 
00155     {
00156         // Rename a,b,c and adjust bounding interval d
00157         if (fb*fc > 0.0) 
00158         {
00159             c = a;
00160             fc = fa;
00161             d = b - a;
00162             e = d;
00163         }
00164         if (std::abs(fc) < std::abs(fb)) 
00165         {
00166             a = b;
00167             b = c;
00168             c = a;
00169             fa = fb;
00170             fb = fc;
00171             fc = fa;
00172         }
00173         G4double Tol1 = 2.0*precision*std::abs(b) + 0.5*tolerance;
00174         G4double xm = 0.5*(c-b);
00175         if (std::abs(xm) <= Tol1 || fb == 0.0) 
00176         {
00177             root = b;
00178             return true;
00179         }
00180         // Inverse quadratic interpolation
00181         if (std::abs(e) >= Tol1 && std::abs(fa) > std::abs(fb)) 
00182         {
00183             G4double ss = fb/fa;
00184             G4double p = 0.0;
00185             G4double q = 0.0;
00186             if (a == c) 
00187             {
00188                 p = 2.0*xm*ss;
00189                 q = 1.0 - ss;
00190             } 
00191             else 
00192             {
00193                 q = fa/fc;
00194                 G4double r = fb/fc;
00195                 p = ss*(2.0*xm*q*(q-r)-(b-a)*(r-1.0));
00196                 q = (q-1.0)*(r-1.0)*(ss-1.0);
00197             }
00198             // Check bounds
00199             if (p > 0.0) q = -q;
00200             p = std::abs(p);
00201             G4double min1 = 3.0*xm*q-std::abs(Tol1*q);
00202             G4double min2 = std::abs(e*q);
00203             if (2.0*p < std::min(min1,min2)) 
00204             {
00205                 // Interpolation
00206                 e = d;
00207                 d = p/q;
00208             } 
00209             else 
00210             {
00211                 // Bisection
00212                 d = xm;
00213                 e = d;
00214             }
00215         } 
00216         else 
00217         {
00218             // Bounds decreasing too slowly, use bisection
00219             d = xm;
00220             e = d;
00221         }
00222         // Move last guess to a 
00223         a = b;
00224         fa = fb;
00225         if (std::abs(d) > Tol1) b += d;
00226         else 
00227         {
00228             if (xm >= 0.0) b += std::abs(Tol1);
00229             else b -= std::abs(Tol1);
00230         }
00231         fb = theFunction(b);
00232     }
00233     G4cerr << "G4Solver::Brent: Number of iterations exceeded." << G4endl;
00234     return false;
00235 }

template<class Function>
G4bool G4Solver< Function >::Crenshaw ( Function &  theFunction  ) 

Definition at line 240 of file G4Solver.icc.

References G4cerr, and G4endl.

Referenced by G4StatMFMacroTemperature::CalcTemperature().

00241 {
00242     // Check the interval before start
00243     if (a > b || std::abs(a-b) <= tolerance) 
00244     {
00245         G4cerr << "G4Solver::Crenshaw: The interval must be properly set." << G4endl;
00246         return false;
00247     }
00248 
00249     G4double fa = theFunction(a);
00250     if (fa == 0.0) 
00251     {
00252         root = a;
00253         return true;
00254     }
00255 
00256     G4double Mlast = a;
00257 
00258     G4double fb = theFunction(b);
00259     if (fb == 0.0)
00260     {
00261         root = b;
00262         return true;
00263     }
00264 
00265     if (fa*fb > 0.0) 
00266     {
00267         G4cerr << "G4Solver::Crenshaw: The interval must include a root." << G4endl;
00268         return false;
00269     }
00270 
00271     
00272     for (G4int i=0; i < MaxIter; i++) 
00273       {
00274         G4double c = 0.5 * (b + a);
00275         G4double fc = theFunction(c);
00276         if (fc == 0.0 || std::abs(c - a) < tolerance) 
00277           {
00278             root = c;
00279             return true;
00280           }
00281 
00282         if (fc * fa  > 0.0)
00283           {
00284             G4double tmp = a;
00285             a = b;
00286             b = tmp;
00287             tmp = fa;
00288             fa = fb;
00289             fb = tmp;
00290           }
00291         
00292         G4double fc0 = fc - fa;
00293         G4double fb1 = fb - fc;
00294         G4double fb0 = fb - fa;
00295         if (fb * fb0 < 2.0 * fc * fc0)
00296           {
00297             b = c;
00298             fb = fc;
00299           }
00300         else 
00301           {
00302             G4double B = (c - a) / fc0;
00303             G4double C = (fc0 - fb1) / (fb1 * fb0);
00304             G4double M = a - B * fa * (1.0 - C * fc);
00305             G4double fM = theFunction(M);
00306             if (fM == 0.0 || std::abs(M - Mlast) < tolerance)
00307               {
00308                 root = M;
00309                 return true;
00310               }
00311             Mlast = M;
00312             if (fM * fa < 0.0)
00313               {
00314                 b = M;
00315                 fb = fM;
00316               }
00317             else 
00318             {
00319                 a = M;
00320                 fa = fM;
00321                 b = c;
00322                 fb = fc;
00323             }
00324         }
00325     }
00326     return false;
00327 }

template<class Function>
G4double G4Solver< Function >::GetIntervalLowerLimit ( void   )  const [inline]

Definition at line 72 of file G4Solver.hh.

00072 {return a;}

template<class Function>
G4double G4Solver< Function >::GetIntervalUpperLimit ( void   )  const [inline]

Definition at line 73 of file G4Solver.hh.

00073 {return b;}

template<class Function>
G4int G4Solver< Function >::GetMaxIterations ( void   )  const [inline]

Definition at line 65 of file G4Solver.hh.

00065 {return MaxIter;}

template<class Function>
G4double G4Solver< Function >::GetRoot ( void   )  const [inline]

Definition at line 77 of file G4Solver.hh.

Referenced by G4StatMFMacroMultiplicity::CalcChemicalPotentialMu(), G4StatMFMacroChemicalPotential::CalcChemicalPotentialNu(), and G4StatMFMacroTemperature::CalcTemperature().

00077 {return root;}

template<class Function>
G4double G4Solver< Function >::GetTolerance ( void   )  const [inline]

Definition at line 68 of file G4Solver.hh.

00068 {return tolerance;}

template<class Function>
G4bool G4Solver< Function >::operator!= ( const G4Solver< Function > &  right  )  const

Definition at line 65 of file G4Solver.cc.

References G4Solver< Function >::operator==().

00066 {
00067     return !operator==(right);
00068 }

template<class Function>
G4Solver< Function > & G4Solver< Function >::operator= ( const G4Solver< Function > &  right  ) 

Definition at line 47 of file G4Solver.cc.

References G4Solver< Function >::a, G4Solver< Function >::b, G4Solver< Function >::MaxIter, G4Solver< Function >::root, and G4Solver< Function >::tolerance.

00048 {
00049     MaxIter = right.MaxIter;
00050     tolerance = right.tolerance;
00051     a = right.a;
00052     b = right.b;
00053     root = right.root;  
00054     return *this;
00055 }

template<class Function>
G4bool G4Solver< Function >::operator== ( const G4Solver< Function > &  right  )  const

Definition at line 58 of file G4Solver.cc.

Referenced by G4Solver< Function >::operator!=().

00059 {
00060     if (this == &right) return true;
00061     else return false;
00062 }

template<class Function>
G4bool G4Solver< Function >::RegulaFalsi ( Function &  theFunction  ) 

Definition at line 78 of file G4Solver.icc.

References G4cerr, and G4endl.

00079 {
00080     // Check the interval before start
00081     if (a > b || std::abs(a-b) <= tolerance) 
00082     {
00083         G4cerr << "G4Solver::RegulaFalsi: The interval must be properly set." << G4endl;
00084         return false;
00085     }
00086     G4double fa = theFunction(a);
00087     G4double fb = theFunction(b);
00088     if (fa*fb > 0.0) 
00089     {
00090         G4cerr << "G4Solver::RegulaFalsi: The interval must include a root." << G4endl;
00091         return false;
00092     }
00093     
00094     G4double eps=tolerance*(b-a);
00095         
00096         
00097     // Finding the root
00098     for (G4int i = 0; i < MaxIter; i++) 
00099     {
00100         G4double c = (a*fb-b*fa)/(fb-fa);
00101         G4double delta = std::min(std::abs(c-a),std::abs(b-c));
00102         if (delta < eps) 
00103         {
00104             root = c;
00105             return true;
00106         }
00107         G4double fc = theFunction(c);
00108         if (fc == 0.0) 
00109         {
00110             root = c;
00111             return true;
00112         }
00113         if (fa*fc < 0.0) 
00114         {
00115             b=c;
00116             fb=fc;
00117         } 
00118         else 
00119         {
00120             a=c;
00121             fa=fc;
00122         }
00123     }
00124     G4cerr << "G4Solver::Bisection: Excedded maximum number of iterations whithout convegence." << G4endl;
00125     return false;
00126 
00127 }       

template<class Function>
void G4Solver< Function >::SetIntervalLimits ( const G4double  Limit1,
const G4double  Limit2 
)

Definition at line 333 of file G4Solver.icc.

References G4cerr, and G4endl.

Referenced by G4StatMFMacroMultiplicity::CalcChemicalPotentialMu(), G4StatMFMacroChemicalPotential::CalcChemicalPotentialNu(), and G4StatMFMacroTemperature::CalcTemperature().

00334 {
00335     if (std::abs(Limit1-Limit2) <= tolerance) 
00336     {
00337         G4cerr << "G4Solver::SetIntervalLimits: Interval must be wider than tolerance." << G4endl;
00338         return;
00339     }
00340     if (Limit1 < Limit2) 
00341     {
00342         a = Limit1;
00343         b = Limit2;
00344     } 
00345     else 
00346     {
00347         a = Limit2;
00348         b = Limit1;
00349     }
00350     return;
00351 }

template<class Function>
void G4Solver< Function >::SetMaxIterations ( const G4int  iterations  )  [inline]

Definition at line 66 of file G4Solver.hh.

00066 {MaxIter=iterations;}

template<class Function>
void G4Solver< Function >::SetTolerance ( const G4double  epsilon  )  [inline]

Definition at line 69 of file G4Solver.hh.

00069 {tolerance = epsilon;}


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:53:24 2013 for Geant4 by  doxygen 1.4.7