Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/global/management/include/G4Pow.hh

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 /global/management/include/G4Pow.hh (Version 11.3.0) and /global/management/include/G4Pow.hh (Version 10.6)


  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 // G4Pow                                       <<  26 //
                                                   >>  27 //
                                                   >>  28 // -------------------------------------------------------------------
                                                   >>  29 //
                                                   >>  30 // Class G4Pow
 27 //                                                 31 //
 28 // Class description:                              32 // Class description:
 29 //                                                 33 //
 30 // Utility singleton class for the fast comput     34 // Utility singleton class for the fast computation of log and pow
 31 // functions. Integer argument should be in th <<  35 // functions. Integer argument should in the interval 0-512, no
 32 // check is performed inside these methods for     36 // check is performed inside these methods for performance reasons.
 33 // For factorial integer argument should be in     37 // For factorial integer argument should be in the interval 0-170
 34 // Computations with double arguments are fast     38 // Computations with double arguments are fast for the interval
 35 // 0.002-511.5 for all functions except expone <<  39 // 0.002-511.5 for all functions except exponent, which is computed 
 36 // for the interval 0-84.4, standard library i     40 // for the interval 0-84.4, standard library is used in the opposite case
 37                                                    41 
 38 // Author: Vladimir Ivanchenko, 23.05.2009     <<  42 // Author: Vladimir Ivanchenko
 39 // ------------------------------------------- <<  43 //
 40 #ifndef G4Pow_hh                               <<  44 // Creation date: 23.05.2009
 41 #define G4Pow_hh 1                             <<  45 // -------------------------------------------------------------------
                                                   >>  46 
                                                   >>  47 #ifndef G4Pow_h
                                                   >>  48 #define G4Pow_h 1
 42                                                    49 
 43 #include "G4DataVector.hh"                     << 
 44 #include "G4Exp.hh"                            << 
 45 #include "G4Log.hh"                            << 
 46 #include "globals.hh"                              50 #include "globals.hh"
                                                   >>  51 #include "G4Log.hh"
                                                   >>  52 #include "G4Exp.hh"
                                                   >>  53 #include "G4DataVector.hh"
 47                                                    54 
 48 class G4Pow                                        55 class G4Pow
 49 {                                                  56 {
 50  public:                                       <<  57 
 51   static G4Pow* GetInstance();                 <<  58   public:
 52   ~G4Pow() = default;                          <<  59 
 53                                                <<  60     static G4Pow* GetInstance();
 54   // Fast computation of Z^1/3                 <<  61    ~G4Pow();
 55   //                                           <<  62 
 56   inline G4double Z13(G4int Z) const;          <<  63     // Fast computation of Z^1/3
 57   G4double A13(G4double A) const;              <<  64     //
 58                                                <<  65     inline G4double Z13(G4int Z) const;
 59   // Fast computation of Z^2/3                 <<  66            G4double A13(G4double A) const;
 60   //                                           <<  67 
 61   inline G4double Z23(G4int Z) const;          <<  68     // Fast computation of Z^2/3
 62   inline G4double A23(G4double A) const;       <<  69     //
 63                                                <<  70     inline G4double Z23(G4int Z) const;
 64   // Fast computation of log(Z)                <<  71     inline G4double A23(G4double A) const;
 65   //                                           <<  72 
 66   inline G4double logZ(G4int Z) const;         <<  73     // Fast computation of log(Z)
 67   inline G4double logA(G4double A) const;      <<  74     //
 68   inline G4double logX(G4double x) const;      <<  75     inline G4double logZ(G4int Z) const;
 69                                                <<  76     inline G4double logA(G4double A) const;
 70   // Fast computation of log10(Z)              <<  77     inline G4double logX(G4double x) const;
 71   //                                           <<  78 
 72   inline G4double log10Z(G4int Z) const;       <<  79     // Fast computation of log10(Z)
 73   inline G4double log10A(G4double A) const;    <<  80     //
 74                                                <<  81     inline G4double log10Z(G4int Z) const;
 75   // Fast computation of exp(X)                <<  82     inline G4double log10A(G4double A) const;
 76   //                                           <<  83 
 77   inline G4double expA(G4double A) const;      <<  84     // Fast computation of exp(X)
 78                                                <<  85     //
 79   // Fast computation of pow(Z,X)              <<  86     inline G4double expA(G4double A) const;
 80   //                                           <<  87 
 81   inline G4double powZ(G4int Z, G4double y) co <<  88     // Fast computation of pow(Z,X)
 82   inline G4double powA(G4double A, G4double y) <<  89     //
 83   G4double powN(G4double x, G4int n) const;    <<  90     inline G4double powZ(G4int Z, G4double y) const;
 84                                                <<  91     inline G4double powA(G4double A, G4double y) const;
 85   // Fast factorial                            <<  92            G4double powN(G4double x, G4int n) const;
 86   //                                           <<  93 
 87   inline G4double factorial(G4int Z) const;    <<  94     // Fast factorial
 88   inline G4double logfactorial(G4int Z) const; <<  95     //
 89                                                <<  96     inline G4double factorial(G4int Z) const;
 90  private:                                      <<  97     inline G4double logfactorial(G4int Z) const;
 91   G4Pow();                                     <<  98 
 92                                                <<  99   private:
 93   G4double A13Low(const G4double, const G4bool << 100 
 94   G4double A13High(const G4double, const G4boo << 101     G4Pow();
 95                                                << 102 
 96   inline G4double logBase(G4double x) const;   << 103     G4double A13Low(const G4double, const bool) const;
 97                                                << 104     G4double A13High(const G4double, const bool) const;
 98   static G4Pow* fpInstance;                    << 105 
 99                                                << 106     inline G4double logBase(G4double x) const;
100   const G4double onethird = 1.0 / 3.0;         << 107 
101   const G4int max2        = 5;                 << 108     static G4Pow* fpInstance;
102                                                << 109 
103   G4double maxA;                               << 110     const G4double onethird;
104   G4double maxLowA;                            << 111     const G4int    max2;
105   G4double maxA2;                              << 112 
106   G4double maxAexp;                            << 113     G4double maxA;
107                                                << 114     G4double maxLowA;
108   G4DataVector ener;                           << 115     G4double maxA2;
109   G4DataVector logen;                          << 116     G4double maxAexp;
110   G4DataVector pz13;                           << 117 
111   G4DataVector lowa13;                         << 118     G4DataVector ener;
112   G4DataVector lz;                             << 119     G4DataVector logen;
113   G4DataVector lz2;                            << 120     G4DataVector pz13;
114   G4DataVector fexp;                           << 121     G4DataVector lowa13;
115   G4DataVector fact;                           << 122     G4DataVector lz;
116   G4DataVector logfact;                        << 123     G4DataVector lz2;
                                                   >> 124     G4DataVector fexp;
                                                   >> 125     G4DataVector fact;
                                                   >> 126     G4DataVector logfact;
117 };                                                127 };
118                                                   128 
119 // -----------------------------               << 129 // -------------------------------------------------------------------
120 // Inline methods implementation               << 
121 // -----------------------------               << 
122                                                   130 
123 inline G4double G4Pow::Z13(G4int Z) const { re << 131 inline G4double G4Pow::Z13(G4int Z) const
                                                   >> 132 {
                                                   >> 133   return pz13[Z];
                                                   >> 134 }
124                                                   135 
125 inline G4double G4Pow::Z23(G4int Z) const         136 inline G4double G4Pow::Z23(G4int Z) const
126 {                                                 137 {
127   G4double x = Z13(Z);                            138   G4double x = Z13(Z);
128   return x * x;                                << 139   return x*x;
129 }                                                 140 }
130                                                   141 
131 inline G4double G4Pow::A23(G4double A) const      142 inline G4double G4Pow::A23(G4double A) const
132 {                                                 143 {
133   G4double x = A13(A);                            144   G4double x = A13(A);
134   return x * x;                                << 145   return x*x;
135 }                                                 146 }
136                                                   147 
137 inline G4double G4Pow::logZ(G4int Z) const { r << 148 inline G4double G4Pow::logZ(G4int Z) const
                                                   >> 149 {
                                                   >> 150   return lz[Z];
                                                   >> 151 }
138                                                   152 
139 inline G4double G4Pow::logBase(G4double a) con    153 inline G4double G4Pow::logBase(G4double a) const
140 {                                                 154 {
141   G4double res;                                   155   G4double res;
142   if(a <= maxA2)                               << 156   if(a <= maxA2) 
143   {                                               157   {
144     G4int i = G4int(max2 * (a - 1) + 0.5);     << 158     G4int i = G4int(max2*(a - 1) + 0.5);
145     if(i > max2)                               << 159     if(i > max2) { i = max2; }
146     {                                          << 160     G4double x = a/(G4double(i)/max2 + 1) - 1;
147       i = max2;                                << 161     res = lz2[i] + x*(1.0 - (0.5 - onethird*x)*x);
148     }                                          << 
149     G4double x = a / (G4double(i) / max2 + 1)  << 
150     res        = lz2[i] + x * (1.0 - (0.5 - on << 
151   }                                               162   }
152   else if(a <= maxA)                              163   else if(a <= maxA)
153   {                                               164   {
154     G4int i    = G4int(a + 0.5);               << 165     G4int i = G4int(a + 0.5);
155     G4double x = a / G4double(i) - 1;          << 166     G4double x = a/G4double(i) - 1;
156     res        = lz[i] + x * (1.0 - (0.5 - one << 167     res = lz[i] + x*(1.0 - (0.5 - onethird*x)*x);
157   }                                               168   }
158   else                                            169   else
159   {                                               170   {
160     res = G4Log(a);                               171     res = G4Log(a);
161   }                                               172   }
162   return res;                                     173   return res;
163 }                                                 174 }
164                                                   175 
165 inline G4double G4Pow::logA(G4double A) const     176 inline G4double G4Pow::logA(G4double A) const
166 {                                                 177 {
167   return (1.0 <= A ? logBase(A) : -logBase(1.  << 178   return (1.0 <= A ? logBase(A) : -logBase(1./A));
168 }                                                 179 }
169                                                   180 
170 inline G4double G4Pow::logX(G4double x) const     181 inline G4double G4Pow::logX(G4double x) const
171 {                                                 182 {
172   G4double res = 0.0;                             183   G4double res = 0.0;
173   G4double a   = (1.0 <= x) ? x : 1.0 / x;     << 184   G4double a = (1.0 <= x) ? x : 1.0/x;
174                                                   185 
175   if(a <= maxA)                                << 186   if(a <= maxA) 
176   {                                               187   {
177     res = logBase(a);                             188     res = logBase(a);
178   }                                               189   }
179   else if(a <= ener[2])                           190   else if(a <= ener[2])
180   {                                               191   {
181     res = logen[1] + logBase(a / ener[1]);     << 192     res = logen[1] + logBase(a/ener[1]);
182   }                                               193   }
183   else if(a <= ener[3])                           194   else if(a <= ener[3])
184   {                                               195   {
185     res = logen[2] + logBase(a / ener[2]);     << 196     res = logen[2] + logBase(a/ener[2]);
186   }                                               197   }
187   else                                            198   else
188   {                                               199   {
189     res = G4Log(a);                               200     res = G4Log(a);
190   }                                               201   }
191                                                   202 
192   if(1.0 > x)                                  << 203   if(1.0 > x) { res = -res; }
193   {                                            << 
194     res = -res;                                << 
195   }                                            << 
196   return res;                                     204   return res;
197 }                                                 205 }
198                                                   206 
199 inline G4double G4Pow::log10Z(G4int Z) const { << 207 inline G4double G4Pow::log10Z(G4int Z) const
                                                   >> 208 {
                                                   >> 209   return lz[Z]/lz[10];
                                                   >> 210 }
200                                                   211 
201 inline G4double G4Pow::log10A(G4double A) cons << 212 inline G4double G4Pow::log10A(G4double A) const
                                                   >> 213 {
                                                   >> 214   return logX(A)/lz[10];
                                                   >> 215 }
202                                                   216 
203 inline G4double G4Pow::expA(G4double A) const     217 inline G4double G4Pow::expA(G4double A) const
204 {                                                 218 {
205   G4double res;                                   219   G4double res;
206   G4double a = (0.0 <= A) ? A : -A;               220   G4double a = (0.0 <= A) ? A : -A;
207                                                   221 
208   if(a <= maxAexp)                                222   if(a <= maxAexp)
209   {                                               223   {
210     G4int i    = G4int(2 * a + 0.5);           << 224     G4int i = G4int(2*a + 0.5);
211     G4double x = a - i * 0.5;                  << 225     G4double x = a - i*0.5;
212     res        = fexp[i] * (1.0 + x * (1.0 + 0 << 226     res = fexp[i]*(1.0 + x*(1.0 + 0.5*(1.0 + onethird*x)*x));
213   }                                               227   }
214   else                                            228   else
215   {                                               229   {
216     res = G4Exp(a);                               230     res = G4Exp(a);
217   }                                               231   }
218   if(0.0 > A)                                  << 232   if(0.0 > A) { res = 1.0/res; }
219   {                                            << 
220     res = 1.0 / res;                           << 
221   }                                            << 
222   return res;                                     233   return res;
223 }                                                 234 }
224                                                   235 
225 inline G4double G4Pow::powZ(G4int Z, G4double     236 inline G4double G4Pow::powZ(G4int Z, G4double y) const
226 {                                                 237 {
227   return expA(y * lz[Z]);                      << 238   return expA(y*lz[Z]);
228 }                                                 239 }
229                                                   240 
230 inline G4double G4Pow::powA(G4double A, G4doub    241 inline G4double G4Pow::powA(G4double A, G4double y) const
231 {                                                 242 {
232   return (0.0 == A ? 0.0 : expA(y * logX(A))); << 243   return (0.0 == A ? 0.0 : expA(y*logX(A)));
233 }                                                 244 }
234                                                   245 
235 inline G4double G4Pow::factorial(G4int Z) cons << 246 inline G4double G4Pow::factorial(G4int Z) const
                                                   >> 247 {
                                                   >> 248   return fact[Z];
                                                   >> 249 }
                                                   >> 250 
                                                   >> 251 inline G4double G4Pow::logfactorial(G4int Z) const
                                                   >> 252 {
                                                   >> 253   return logfact[Z];
                                                   >> 254 }
236                                                   255 
237 inline G4double G4Pow::logfactorial(G4int Z) c << 256 // -------------------------------------------------------------------
238                                                   257 
239 #endif                                            258 #endif
240                                                   259