| Geant4 Cross Reference |
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 class implemenation << 26 // $Id: G4Pow.cc,v 1.6 2011-01-08 15:39:03 vnivanch Exp $
>> 27 // GEANT4 tag $Name: geant4-09-04-patch-02 $
27 // 28 //
28 // Author: Vladimir Ivanchenko, 23.05.2009 << 29 // -------------------------------------------------------------------
29 // Revisions: << 30 //
30 // - 23.03.2018, M.Novak: increased accuracy o << 31 // GEANT4 Class file
31 // [1/4,4] interval by << 32 //
32 // ------------------------------------------- << 33 //
>> 34 // File name: G4Pow
>> 35 //
>> 36 // Author: Vladimir Ivanchenko
>> 37 //
>> 38 // Creation date: 23.05.2009
>> 39 //
>> 40 // Modifications:
>> 41 // 08.01.2011 V.Ivanchenko extended maxZ from 256 to 512
>> 42 //
>> 43 // -------------------------------------------------------------------
33 44
34 #include "G4Pow.hh" 45 #include "G4Pow.hh"
35 #include "G4Threading.hh" <<
36 46
37 G4Pow* G4Pow::fpInstance = nullptr; << 47 G4Pow* G4Pow::fInstance = 0;
38 48
39 // ------------------------------------------- 49 // -------------------------------------------------------------------
40 50
41 G4Pow* G4Pow::GetInstance() 51 G4Pow* G4Pow::GetInstance()
42 { 52 {
43 if(fpInstance == nullptr) << 53 if (fInstance == 0)
44 { 54 {
45 static G4Pow geant4pow; << 55 static G4Pow manager;
46 fpInstance = &geant4pow; << 56 fInstance = &manager;
47 } 57 }
48 return fpInstance; << 58 return fInstance;
49 } 59 }
50 60
51 // ------------------------------------------- 61 // -------------------------------------------------------------------
52 62
53 G4Pow::G4Pow() 63 G4Pow::G4Pow()
54 { 64 {
55 #ifdef G4MULTITHREADED << 65 const G4int maxZ = 512;
56 if(G4Threading::IsWorkerThread()) <<
57 { <<
58 G4Exception("G4Pow::G4Pow()", "InvalidSetu <<
59 "Attempt to instantiate G4Pow <<
60 } <<
61 #endif <<
62 const G4int maxZ = 512; <<
63 const G4int maxZfact = 170; <<
64 const G4int numLowA = 17; <<
65 <<
66 maxA = -0.6 + maxZ; <<
67 maxLowA = 4.0; <<
68 maxA2 = 1.25 + max2 * 0.2; <<
69 maxAexp = -0.76 + maxZfact * 0.5; <<
70 <<
71 ener.resize(max2 + 1, 1.0); <<
72 logen.resize(max2 + 1, 0.0); <<
73 lz2.resize(max2 + 1, 0.0); <<
74 pz13.resize(maxZ, 0.0); <<
75 lowa13.resize(numLowA, 0.0); <<
76 lz.resize(maxZ, 0.0); <<
77 fexp.resize(maxZfact, 0.0); <<
78 fact.resize(maxZfact, 0.0); <<
79 logfact.resize(maxZ, 0.0); <<
80 <<
81 G4double f = 1.0; <<
82 G4double logf = 0.0; <<
83 fact[0] = 1.0; <<
84 fexp[0] = 1.0; <<
85 <<
86 for(G4int i = 1; i <= max2; ++i) <<
87 { <<
88 ener[i] = powN(500., i); <<
89 logen[i] = G4Log(ener[i]); <<
90 lz2[i] = G4Log(1.0 + i * 0.2); <<
91 } <<
92 <<
93 for(G4int i = 1; i < maxZ; ++i) <<
94 { <<
95 auto x = G4double(i); <<
96 pz13[i] = std::pow(x, onethird); <<
97 lz[i] = G4Log(x); <<
98 if(i < maxZfact) <<
99 { <<
100 f *= x; <<
101 fact[i] = f; <<
102 fexp[i] = G4Exp(0.5 * i); <<
103 } <<
104 logf += lz[i]; <<
105 logfact[i] = logf; <<
106 } <<
107 <<
108 for(G4int i = 4; i < numLowA; ++i) <<
109 { <<
110 lowa13[i] = std::pow(0.25 * i, onethird); <<
111 } <<
112 } <<
113 <<
114 // ------------------------------------------- <<
115 <<
116 G4double G4Pow::A13(G4double A) const <<
117 { <<
118 G4double res = 0.; <<
119 if(A > 0.) <<
120 { <<
121 const bool invert = (A < 1.); <<
122 const G4double a = invert ? 1. / A : A; <<
123 res = (a < maxLowA) ? A13Low <<
124 } <<
125 return res; <<
126 } <<
127 <<
128 // ------------------------------------------- <<
129 66
130 G4double G4Pow::A13High(const G4double a, cons << 67 pz13.resize(maxZ,0.0);
131 { << 68 lz.resize(maxZ,0.0);
132 G4double res; << 69 fact.resize(maxZ,0.0);
133 if(a < maxA) << 70 logfact.resize(maxZ,0.0);
134 { << 71
135 const auto i = static_cast<G4int>(a + << 72 onethird = 1.0/3.0;
136 const G4double x = (a / i - 1.) * onethird << 73 G4double f = 1.0;
137 res = pz13[i] * (1. + x - x * << 74
>> 75 for(G4int i=1; i<maxZ; ++i)
>> 76 {
>> 77 G4double x = G4double(i);
>> 78 pz13[i] = std::pow(x,onethird);
>> 79 lz[i] = std::log(x);
>> 80 f *= x;
>> 81 fact[i] = f;
>> 82 logfact[i] = std::log(f);
138 } 83 }
139 else << 84 fact[0] = 1.0;
140 { << 85 logfact[0] = 0.0;
141 res = G4Exp(G4Log(a) * onethird); <<
142 } <<
143 res = invert ? 1. / res : res; <<
144 return res; <<
145 } 86 }
146 87
147 // ------------------------------------------- 88 // -------------------------------------------------------------------
148 89
149 G4double G4Pow::A13Low(const G4double a, const << 90 G4Pow::~G4Pow()
150 { << 91 {}
151 G4double res; <<
152 const auto i = static_cast<G4int>(4. * ( <<
153 const G4double y = 0.25 * i; <<
154 const G4double x = (a / y - 1.) * onethird; <<
155 res = lowa13[i] * (1. + x - x * <<
156 res = invert ? 1. / res : res; <<
157 return res; <<
158 } <<
159 92
160 // ------------------------------------------- 93 // -------------------------------------------------------------------
161 94
162 G4double G4Pow::powN(G4double x, G4int n) cons << 95 G4double G4Pow::powN(G4double x, G4int n)
163 { 96 {
164 if(0.0 == x) <<
165 { <<
166 return 0.0; <<
167 } <<
168 if(std::abs(n) > 8) <<
169 { <<
170 return std::pow(x, G4double(n)); <<
171 } <<
172 G4double res = 1.0; 97 G4double res = 1.0;
173 if(n >= 0) << 98 if(n >= 0) { for(G4int i=0; i<n; ++i) { res *= x; } }
174 { << 99 else if((n < 0) && (x != 0.0))
175 for(G4int i = 0; i < n; ++i) <<
176 { <<
177 res *= x; <<
178 } <<
179 } <<
180 else if(n < 0) <<
181 { 100 {
182 G4double y = 1.0 / x; << 101 G4double y = 1.0/x;
183 G4int nn = -n; << 102 G4int nn = -n;
184 for(G4int i = 0; i < nn; ++i) << 103 for(G4int i=0; i<nn; ++i) { res *= y; }
185 { <<
186 res *= y; <<
187 } <<
188 } 104 }
189 return res; 105 return res;
190 } 106 }
191 107