Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4Pow.cc
이 파일의 문서화 페이지로 가기
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4Pow.cc 109086 2018-03-26 08:20:25Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4Pow
34 //
35 // Author: Vladimir Ivanchenko
36 //
37 // Creation date: 23.05.2009
38 //
39 // Modifications:
40 // 08.01.2011 V.Ivanchenko extended maxZ from 256 to 512
41 // 02.05.2013 V.Ivanchenko added expA and logX methods,
42 // revised A13, logA, powZ, powA to improved accuracy
43 // 23.03.2018 M.Novak increased accuracy of A13 on the most critical
44 // [1/4,4] interval by introducing a denser(0.25) grid
45 //
46 // -------------------------------------------------------------------
47 
48 #include "G4Pow.hh"
49 #ifdef G4MULTITHREADED
50 #include "G4Threading.hh"
51 #endif
52 
54 
55 // -------------------------------------------------------------------
56 
58 {
59  if (fpInstance == 0)
60  {
61  static G4Pow geant4pow;
62  fpInstance = &geant4pow;
63  }
64  return fpInstance;
65 }
66 
67 // -------------------------------------------------------------------
68 
70  : onethird(1.0/3.0), max2(5)
71 {
72 #ifdef G4MULTITHREADED
74  {
75  G4Exception ("G4Pow::G4Pow()", "InvalidSetup", FatalException,
76  "Attempt to instantiate G4Pow in worker thread!");
77  }
78 #endif
79  const G4int maxZ = 512;
80  const G4int maxZfact = 170;
81  const G4int numLowA = 17;
82 
83  maxA = -0.6 + maxZ;
84  maxLowA = 4.0;
85  maxA2 = 1.25 + max2*0.2;
86  maxAexp = -0.76+ maxZfact*0.5;
87 
88  ener.resize(max2+1,1.0);
89  logen.resize(max2+1,0.0);
90  lz2.resize(max2+1,0.0);
91  pz13.resize(maxZ,0.0);
92  lowa13.resize(numLowA,0.0);
93  lz.resize(maxZ,0.0);
94  fexp.resize(maxZfact,0.0);
95  fact.resize(maxZfact,0.0);
96  logfact.resize(maxZ,0.0);
97 
98  G4double f = 1.0;
99  G4double logf = 0.0;
100  fact[0] = 1.0;
101  fexp[0] = 1.0;
102 
103  for(G4int i=1; i<=max2; ++i)
104  {
105  ener[i] = powN(500.,i);
106  logen[i]= G4Log(ener[i]);
107  lz2[i] = G4Log(1.0 + i*0.2);
108  }
109 
110  for(G4int i=1; i<maxZ; ++i)
111  {
112  G4double x = G4double(i);
113  pz13[i] = std::pow(x,onethird);
114  lz[i] = G4Log(x);
115  if(i < maxZfact)
116  {
117  f *= x;
118  fact[i] = f;
119  fexp[i] = G4Exp(0.5*i);
120  }
121  logf += lz[i];
122  logfact[i] = logf;
123  }
124 
125  for (G4int i=4; i<numLowA; ++i) {
126  lowa13[i] = std::pow(0.25*i,onethird);
127  }
128 
129 }
130 
131 // -------------------------------------------------------------------
132 
134 {}
135 
136 // -------------------------------------------------------------------
137 
139  G4double res = 0.;
140  if (A>0.) {
141  const bool invert = (A<1.);
142  const G4double a = invert ? 1./A : A;
143  res = (a<maxLowA) ? A13Low(a, invert) : A13High(a, invert);
144  }
145  return res;
146 }
147 
148 // -------------------------------------------------------------------
149 
150 G4double G4Pow::A13High(const G4double a, const bool invert) const {
151  G4double res;
152  if (a<maxA) {
153  const G4int i = static_cast<G4int>(a+0.5);
154  const G4double x = (a/i-1.)*onethird;
155  res = pz13[i]*(1.+x-x*x*(1.-1.666667*x));
156  } else {
157  res = G4Exp(G4Log(a)*onethird);
158  }
159  res = invert ? 1./res : res;
160  return res;
161 }
162 
163 // -------------------------------------------------------------------
164 
165 G4double G4Pow::A13Low(const G4double a, const bool invert) const {
166  G4double res;
167  const G4int i = static_cast<G4int>(4.*(a+0.125));
168  const G4double y = 0.25*i;
169  const G4double x = (a/y-1.)*onethird;
170  res = lowa13[i]*(1.+x-x*x*(1.-1.666667*x));
171  res = invert ? 1./res : res;
172  return res;
173 }
174 
175 // -------------------------------------------------------------------
176 
178 {
179  if(0.0 == x) { return 0.0; }
180  if(std::abs(n) > 8) { return std::pow(x, G4double(n)); }
181  G4double res = 1.0;
182  if(n >= 0) { for(G4int i=0; i<n; ++i) { res *= x; } }
183  else if(n < 0)
184  {
185  G4double y = 1.0/x;
186  G4int nn = -n;
187  for(G4int i=0; i<nn; ++i) { res *= y; }
188  }
189  return res;
190 }
G4DataVector fexp
Definition: G4Pow.hh:125
Float_t x
Definition: compare.C:6
G4DataVector logfact
Definition: G4Pow.hh:127
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
G4double A13High(const G4double, const bool) const
Definition: G4Pow.cc:150
G4DataVector lz2
Definition: G4Pow.hh:124
G4DataVector fact
Definition: G4Pow.hh:126
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
G4double A13(G4double A) const
Definition: G4Pow.cc:138
G4DataVector ener
Definition: G4Pow.hh:119
static G4Pow * fpInstance
Definition: G4Pow.hh:109
Float_t y
Definition: compare.C:6
const G4double onethird
Definition: G4Pow.hh:111
~G4Pow()
Definition: G4Pow.cc:133
const G4int max2
Definition: G4Pow.hh:112
G4double maxA
Definition: G4Pow.hh:114
G4double maxLowA
Definition: G4Pow.hh:115
G4DataVector lz
Definition: G4Pow.hh:123
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double powN(G4double x, G4int n) const
Definition: G4Pow.cc:177
static const G4int maxZ
double G4double
Definition: G4Types.hh:76
static G4Pow * GetInstance()
Definition: G4Pow.cc:57
G4bool IsWorkerThread()
Definition: G4Threading.cc:129
G4DataVector pz13
Definition: G4Pow.hh:121
G4DataVector logen
Definition: G4Pow.hh:120
double A(double temperature)
G4Pow()
Definition: G4Pow.cc:69
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
int G4int
Definition: G4Types.hh:78
Definition: G4Pow.hh:56
Char_t n[5]
G4double A13Low(const G4double, const bool) const
Definition: G4Pow.cc:165
G4double maxA2
Definition: G4Pow.hh:116
G4double maxAexp
Definition: G4Pow.hh:117
G4DataVector lowa13
Definition: G4Pow.hh:122