Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4NeutronHPLegendreStore Class Reference

#include <G4NeutronHPLegendreStore.hh>

Public Member Functions

 G4NeutronHPLegendreStore (G4int n)
 
 ~G4NeutronHPLegendreStore ()
 
void Init (G4int i, G4double e, G4int n)
 
void SetNPoints (G4int n)
 
void SetEnergy (G4int i, G4double energy)
 
void SetTemperature (G4int i, G4double temp)
 
void SetCoeff (G4int i, G4int l, G4double coeff)
 
void SetCoeff (G4int i, G4NeutronHPLegendreTable *theTable)
 
G4double GetCoeff (G4int i, G4int l)
 
G4double GetEnergy (G4int i)
 
G4double GetTemperature (G4int i)
 
G4int GetNumberOfPoly (G4int i)
 
G4double SampleDiscreteTwoBody (G4double anEnergy)
 
G4double SampleElastic (G4double anEnergy)
 
G4double Sample (G4double energy)
 
G4double SampleMax (G4double energy)
 
G4double Integrate (G4int k, G4double costh)
 
void InitInterpolation (std::ifstream &aDataFile)
 
void SetManager (G4InterpolationManager &aManager)
 

Detailed Description

Definition at line 37 of file G4NeutronHPLegendreStore.hh.

Constructor & Destructor Documentation

◆ G4NeutronHPLegendreStore()

G4NeutronHPLegendreStore::G4NeutronHPLegendreStore ( G4int  n)
inline

Definition at line 41 of file G4NeutronHPLegendreStore.hh.

◆ ~G4NeutronHPLegendreStore()

G4NeutronHPLegendreStore::~G4NeutronHPLegendreStore ( )
inline

Definition at line 47 of file G4NeutronHPLegendreStore.hh.

48 {
49 delete [] theCoeff;
50 }

Member Function Documentation

◆ GetCoeff()

G4double G4NeutronHPLegendreStore::GetCoeff ( G4int  i,
G4int  l 
)
inline

Definition at line 68 of file G4NeutronHPLegendreStore.hh.

68{return theCoeff[i].GetCoeff(l);}

Referenced by SampleDiscreteTwoBody(), SampleElastic(), and SampleMax().

◆ GetEnergy()

G4double G4NeutronHPLegendreStore::GetEnergy ( G4int  i)
inline

Definition at line 69 of file G4NeutronHPLegendreStore.hh.

69{return theCoeff[i].GetEnergy();}

Referenced by Sample(), SampleDiscreteTwoBody(), SampleElastic(), and SampleMax().

◆ GetNumberOfPoly()

G4int G4NeutronHPLegendreStore::GetNumberOfPoly ( G4int  i)
inline

Definition at line 71 of file G4NeutronHPLegendreStore.hh.

71{return theCoeff[i].GetNumberOfPoly();}

◆ GetTemperature()

G4double G4NeutronHPLegendreStore::GetTemperature ( G4int  i)
inline

Definition at line 70 of file G4NeutronHPLegendreStore.hh.

70{return theCoeff[i].GetTemperature();}

◆ Init()

void G4NeutronHPLegendreStore::Init ( G4int  i,
G4double  e,
G4int  n 
)
inline

Definition at line 52 of file G4NeutronHPLegendreStore.hh.

53 {
54 theCoeff[i].Init(e, n);
55 }
void Init(std::ifstream &aDataFile)

Referenced by G4NeutronHPElasticFS::Init(), G4NeutronHPAngular::Init(), and G4NeutronHPContAngularPar::Sample().

◆ InitInterpolation()

void G4NeutronHPLegendreStore::InitInterpolation ( std::ifstream &  aDataFile)
inline

Definition at line 79 of file G4NeutronHPLegendreStore.hh.

80 {
81 theManager.Init(aDataFile);
82 }
void Init(G4int aScheme, G4int aRange)

Referenced by G4NeutronHPElasticFS::Init(), and G4NeutronHPAngular::Init().

◆ Integrate()

G4double G4NeutronHPLegendreStore::Integrate ( G4int  k,
G4double  costh 
)

Definition at line 301 of file G4NeutronHPLegendreStore.cc.

302{
303 G4double result=0;
305// G4cout <<"the COEFFS "<<k<<" ";
306// G4cout <<theCoeff[k].GetNumberOfPoly()<<" ";
307 for(G4int l=0; l<theCoeff[k].GetNumberOfPoly() ; l++)
308 {
309 result += theCoeff[k].GetCoeff(l)*theLeg.Integrate(l, costh);
310// G4cout << theCoeff[k].GetCoeff(l)<<" ";
311 }
312// G4cout <<G4endl;
313 return result;
314}
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
G4double Integrate(G4int l, G4double costh)

Referenced by Sample().

◆ Sample()

G4double G4NeutronHPLegendreStore::Sample ( G4double  energy)

Definition at line 248 of file G4NeutronHPLegendreStore.cc.

249{
250 G4int i0;
251 G4int low(0), high(0);
252// G4cout << "G4NeutronHPLegendreStore::Sample "<<energy<<" "<<energy<<" "<<nEnergy<<G4endl;
253 for (i0=0; i0<nEnergy; i0++)
254 {
255// G4cout <<"theCoeff["<<i0<<"].GetEnergy() = "<<theCoeff[i0].GetEnergy()<<G4endl;
256 high = i0;
257 if(theCoeff[i0].GetEnergy()>energy) break;
258 }
259 low = std::max(0, high-1);
260// G4cout << "G4NeutronHPLegendreStore::Sample high, low: "<<high<<", "<<low<<G4endl;
261 G4NeutronHPVector theBuffer;
263 G4double x1, x2, y1, y2, y;
264 x1 = theCoeff[low].GetEnergy();
265 x2 = theCoeff[high].GetEnergy();
266// G4cout << "the xes "<<x1<<" "<<x2<<G4endl;
267 G4double costh=0;
268 for(i0=0; i0<601; i0++)
269 {
270 costh = G4double(i0-300)/300.;
271 y1 = Integrate(low, costh);
272 y2 = Integrate(high, costh);
273 y = theInt.Interpolate(theManager.GetScheme(high), energy, x1, x2, y1, y2);
274 theBuffer.SetData(i0, costh, y);
275// G4cout << "Integration "<<low<<" "<<costh<<" "<<y1<<" "<<y2<<" "<<y<<G4endl;
276 }
277 G4double rand = G4UniformRand();
278 G4int it;
279 for (i0=1; i0<601; i0++)
280 {
281 it = i0;
282 if(rand < theBuffer.GetY(i0)/theBuffer.GetY(600)) break;
283// G4cout <<"sampling now "<<i0<<" "
284// << theBuffer.GetY(i0)<<" "
285// << theBuffer.GetY(600)<<" "
286// << rand<<" "
287// << theBuffer.GetY(i0)/theBuffer.GetY(600)<<G4endl;;
288 }
289 if(it==601) it=600;
290// G4cout << "G4NeutronHPLegendreStore::Sample it "<<rand<<" "<<it<<G4endl;
291 G4double norm = theBuffer.GetY(600);
292 if(norm==0) return -DBL_MAX;
293 x1 = theBuffer.GetY(it)/norm;
294 x2 = theBuffer.GetY(it-1)/norm;
295 y1 = theBuffer.GetX(it);
296 y2 = theBuffer.GetX(it-1);
297// G4cout << "G4NeutronHPLegendreStore::Sample x y "<<x1<<" "<<y1<<" "<<x2<<" "<<y2<<G4endl;
298 return theInt.Interpolate(theManager.GetScheme(high), rand, x1, x2, y1, y2);
299}
#define G4UniformRand()
Definition: Randomize.hh:53
G4InterpolationScheme GetScheme(G4int index) const
G4double Interpolate(G4InterpolationScheme aScheme, G4double x, G4double x1, G4double x2, G4double y1, G4double y2) const
G4double Integrate(G4int k, G4double costh)
G4double GetX(G4int i) const
G4double GetY(G4double x)
void SetData(G4int i, G4double x, G4double y)
#define DBL_MAX
Definition: templates.hh:83

◆ SampleDiscreteTwoBody()

G4double G4NeutronHPLegendreStore::SampleDiscreteTwoBody ( G4double  anEnergy)

Definition at line 42 of file G4NeutronHPLegendreStore.cc.

43{
44 G4double result;
45
46 G4int i0;
47 G4int low(0), high(0);
49 for (i0=0; i0<nEnergy; i0++)
50 {
51 high = i0;
52 if(theCoeff[i0].GetEnergy()>anEnergy) break;
53 }
54 low = std::max(0, high-1);
56 G4double x, x1, x2;
57 x = anEnergy;
58 x1 = theCoeff[low].GetEnergy();
59 x2 = theCoeff[high].GetEnergy();
60 G4double theNorm = 0;
61 G4double try01=0, try02=0;
62 G4double max1, max2, costh;
63 max1 = 0; max2 = 0;
64 G4int l,m_tmp;
65 for(i0=0; i0<601; i0++)
66 {
67 costh = G4double(i0-300)/300.;
68 try01 = 0.5;
69 for(m_tmp=0; m_tmp<theCoeff[low].GetNumberOfPoly() ; m_tmp++)
70 {
71 l=m_tmp+1;
72 try01 += (2.*l+1)/2.*theCoeff[low].GetCoeff(m_tmp)*theLeg.Evaluate(l, costh);
73 }
74 if(try01>max1) max1=try01;
75 try02 = 0.5;
76 for(m_tmp=0; m_tmp<theCoeff[high].GetNumberOfPoly() ; m_tmp++)
77 {
78 l=m_tmp+1;
79 try02 += (2.*l+1)/2.*theCoeff[high].GetCoeff(m_tmp)*theLeg.Evaluate(l, costh);
80 }
81 if(try02>max2) max2=try02;
82 }
83 theNorm = theInt.Interpolate(theManager.GetScheme(high), x, x1, x2, max1, max2);
84
85 G4double value, random;
86 G4double v1, v2;
87 do
88 {
89 v1 = 0.5;
90 v2 = 0.5;
91 result = 2.*G4UniformRand()-1.;
92 for(m_tmp=0; m_tmp<theCoeff[low].GetNumberOfPoly() ; m_tmp++)
93 {
94 l=m_tmp+1;
95 G4double legend = theLeg.Evaluate(l, result); // @@@ done to avoid optimization error on SUN
96 v1 += (2.*l+1)/2.*theCoeff[low].GetCoeff(m_tmp)*legend;
97 }
98 for(m_tmp=0; m_tmp<theCoeff[high].GetNumberOfPoly() ; m_tmp++)
99 {
100 l=m_tmp+1;
101 G4double legend = theLeg.Evaluate(l, result); // @@@ done to avoid optimization error on SUN
102 v2 += (2.*l+1)/2.*theCoeff[high].GetCoeff(m_tmp)*legend;
103 }
104 // v1 = std::max(0.,v1); // Workaround in case one of the distributions is fully non-physical.
105 // v2 = std::max(0.,v2);
106 value = theInt.Interpolate(theManager.GetScheme(high), x, x1, x2, v1, v2);
107 random = G4UniformRand();
108 if(0>=theNorm) break; // Workaround for negative cross-section values. @@@@ 31 May 2000
109 }
110 while(random>value/theNorm);
111
112 return result;
113}
G4double Evaluate(G4int l, G4double costh)
G4double GetCoeff(G4int i, G4int l)

Referenced by G4NeutronHPDiscreteTwoBody::Sample().

◆ SampleElastic()

G4double G4NeutronHPLegendreStore::SampleElastic ( G4double  anEnergy)

Definition at line 187 of file G4NeutronHPLegendreStore.cc.

188{
189 G4double result;
190
191 G4int i0;
192 G4int low(0), high(0);
194 for (i0=0; i0<nEnergy; i0++)
195 {
196 high = i0;
197 if(theCoeff[i0].GetEnergy()>anEnergy) break;
198 }
199 low = std::max(0, high-1);
201 G4double x, x1, x2;
202 x = anEnergy;
203 x1 = theCoeff[low].GetEnergy();
204 x2 = theCoeff[high].GetEnergy();
205 G4double theNorm = 0;
206 G4double try01=0, try02=0, try11=0, try12=0;
207 G4double try1, try2;
208 G4int l;
209 for(l=0; l<theCoeff[low].GetNumberOfPoly(); l++)
210 {
211 try01 += (2.*l+1)/2.*theCoeff[low].GetCoeff(l)*theLeg.Evaluate(l, -1.);
212 try11 += (2.*l+1)/2.*theCoeff[low].GetCoeff(l)*theLeg.Evaluate(l, +1.);
213 }
214 for(l=0; l<theCoeff[high].GetNumberOfPoly(); l++)
215 {
216 try02 += (2.*l+1)/2.*theCoeff[high].GetCoeff(l)*theLeg.Evaluate(l, -1.);
217 try12 += (2.*l+1)/2.*theCoeff[high].GetCoeff(l)*theLeg.Evaluate(l, +1.);
218 }
219 try1 = theInt.Interpolate(theManager.GetScheme(high), x, x1, x2, try01, try02);
220 try2 = theInt.Interpolate(theManager.GetScheme(high), x, x1, x2, try11, try12);
221 theNorm = std::max(try1, try2);
222
223 G4double value, random;
224 G4double v1, v2;
225 do
226 {
227 v1 = 0;
228 v2 = 0;
229 result = 2.*G4UniformRand()-1.;
230 for(l=0; l<theCoeff[low].GetNumberOfPoly() ; l++)
231 {
232 G4double legend = theLeg.Evaluate(l, result); // @@@ done to avoid optimization error on SUN
233 v1 += (2.*l+1)/2.*theCoeff[low].GetCoeff(l)*legend;
234 }
235 for(l=0; l<theCoeff[high].GetNumberOfPoly() ; l++)
236 {
237 G4double legend = theLeg.Evaluate(l, result); // @@@ done to avoid optimization error on SUN
238 v2 += (2.*l+1)/2.*theCoeff[high].GetCoeff(l)*legend;
239 }
240 value = theInt.Interpolate(theManager.GetScheme(high), x, x1, x2, v1, v2);
241 random = G4UniformRand();
242 }
243 while(random>value/theNorm);
244
245 return result;
246}

Referenced by G4NeutronHPElasticFS::ApplyYourself().

◆ SampleMax()

G4double G4NeutronHPLegendreStore::SampleMax ( G4double  energy)

Definition at line 117 of file G4NeutronHPLegendreStore.cc.

118{
119 G4double result;
120
121 G4int i0;
122 G4int low(0), high(0);
124 for (i0=0; i0<nEnergy; i0++)
125 {
126 high = i0;
127 if(theCoeff[i0].GetEnergy()>anEnergy) break;
128 }
129 low = std::max(0, high-1);
131 G4double x, x1, x2;
132 x = anEnergy;
133 x1 = theCoeff[low].GetEnergy();
134 x2 = theCoeff[high].GetEnergy();
135 G4double theNorm = 0;
136 G4double try01=0, try02=0;
137 G4double max1, max2, costh;
138 max1 = 0; max2 = 0;
139 G4int l;
140 for(i0=0; i0<601; i0++)
141 {
142 costh = G4double(i0-300)/300.;
143 try01 = 0;
144 for(l=0; l<theCoeff[low].GetNumberOfPoly() ; l++)
145 {
146 try01 += (2.*l+1)/2.*theCoeff[low].GetCoeff(l)*theLeg.Evaluate(l, costh);
147 }
148 if(try01>max1) max1=try01;
149 try02 = 0;
150 for(l=0; l<theCoeff[high].GetNumberOfPoly() ; l++)
151 {
152 try02 += (2.*l+1)/2.*theCoeff[high].GetCoeff(l)*theLeg.Evaluate(l, costh);
153 }
154 if(try02>max2) max2=try02;
155 }
156 theNorm = theInt.Interpolate(theManager.GetScheme(high), x, x1, x2, max1, max2);
157
158 G4double value, random;
159 G4double v1, v2;
160 do
161 {
162 v1 = 0;
163 v2 = 0;
164 result = 2.*G4UniformRand()-1.;
165 for(l=0; l<theCoeff[low].GetNumberOfPoly() ; l++)
166 {
167 G4double legend = theLeg.Evaluate(l, result); // @@@ done to avoid optimization error on SUN
168 v1 += (2.*l+1)/2.*theCoeff[low].GetCoeff(l)*legend;
169 }
170 for(l=0; l<theCoeff[high].GetNumberOfPoly() ; l++)
171 {
172 G4double legend = theLeg.Evaluate(l, result); // @@@ done to avoid optimization error on SUN
173 v2 += (2.*l+1)/2.*theCoeff[high].GetCoeff(l)*legend;
174 }
175 v1 = std::max(0.,v1); // Workaround in case one of the distributions is fully non-physical.
176 v2 = std::max(0.,v2);
177 value = theInt.Interpolate(theManager.GetScheme(high), x, x1, x2, v1, v2);
178 random = G4UniformRand();
179 if(0>=theNorm) break; // Workaround for negative cross-section values. @@@@ 31 May 2000
180 }
181 while(random>value/theNorm);
182
183 return result;
184}

Referenced by G4NeutronHPPhotonDist::GetPhotons(), G4NeutronHPContAngularPar::Sample(), and G4NeutronHPAngular::SampleAndUpdate().

◆ SetCoeff() [1/2]

void G4NeutronHPLegendreStore::SetCoeff ( G4int  i,
G4int  l,
G4double  coeff 
)
inline

◆ SetCoeff() [2/2]

void G4NeutronHPLegendreStore::SetCoeff ( G4int  i,
G4NeutronHPLegendreTable theTable 
)
inline

Definition at line 60 of file G4NeutronHPLegendreStore.hh.

61 {
62 if(i>nEnergy) throw G4HadronicException(__FILE__, __LINE__, "LegendreTableIndex out of range");
63 theCoeff[i] = *theTable;
64// not here -- see G4NeutronHPPhotonDist.cc line 275
65// delete theTable;
66 }

◆ SetEnergy()

void G4NeutronHPLegendreStore::SetEnergy ( G4int  i,
G4double  energy 
)
inline

Definition at line 57 of file G4NeutronHPLegendreStore.hh.

57{ theCoeff[i].SetEnergy(energy); }

◆ SetManager()

void G4NeutronHPLegendreStore::SetManager ( G4InterpolationManager aManager)
inline

Definition at line 84 of file G4NeutronHPLegendreStore.hh.

85 {
86 theManager = aManager;
87 }

Referenced by G4NeutronHPDiscreteTwoBody::Sample(), and G4NeutronHPContAngularPar::Sample().

◆ SetNPoints()

void G4NeutronHPLegendreStore::SetNPoints ( G4int  n)
inline

Definition at line 56 of file G4NeutronHPLegendreStore.hh.

56{ nEnergy = n; }

◆ SetTemperature()

void G4NeutronHPLegendreStore::SetTemperature ( G4int  i,
G4double  temp 
)
inline

Definition at line 58 of file G4NeutronHPLegendreStore.hh.

58{ theCoeff[i].SetTemperature(temp); }

Referenced by G4NeutronHPElasticFS::Init(), and G4NeutronHPAngular::Init().


The documentation for this class was generated from the following files: