ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4SimpleIntegration.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4SimpleIntegration.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 //
27 //
28 // Implementation file for simple integration methods
29 //
30 
31 #include "globals.hh"
32 #include "G4SimpleIntegration.hh"
33 
34 
36  : fFunction(pFunction),
37  fTolerance(.0001),
38  fMaxDepth(100)
39 {
40 }
41 
43  G4double pTolerance)
44  : fFunction(pFunction),
45  fTolerance(pTolerance),
46  fMaxDepth(100)
47 {
48 }
49 
50 
52 {
53 }
54 
55  // Simple integration methods
56 
59  G4double xFinal,
60  G4int iterationNumber )
61 {
62  G4double Step = (xFinal - xInitial)/iterationNumber ;
63  G4double mean = (fFunction(xInitial) + fFunction(xFinal))*0.5 ;
64  G4double x = xInitial ;
65  for(G4int i=1;i<iterationNumber;i++)
66  {
67  x += Step ;
68  mean += fFunction(x) ;
69  }
70  return mean*Step ;
71 }
72 
73 G4double
75  G4double xFinal,
76  G4int iterationNumber )
77 {
78  G4double Step = (xFinal - xInitial)/iterationNumber ;
79  G4double x = xInitial + 0.5*Step;
80  G4double mean = fFunction(x) ;
81  for(G4int i=1;i<iterationNumber;i++)
82  {
83  x += Step ;
84  mean += fFunction(x) ;
85  }
86  return mean*Step ;
87 }
88 
89 G4double
91  G4double xFinal,
92  G4int iterationNumber )
93 {
94  G4double x=0.;
95  static const G4double root = 1.0/std::sqrt(3.0) ;
96  G4double Step = (xFinal - xInitial)/(2.0*iterationNumber) ;
97  G4double delta = Step*root ;
98  G4double mean = 0.0 ;
99  for(G4int i=0;i<iterationNumber;i++)
100  {
101  x = (2*i + 1)*Step ;
102  mean += (fFunction(x+delta) + fFunction(x-delta)) ;
103  }
104  return mean*Step ;
105 }
106 
107 G4double
109  G4double xFinal,
110  G4int iterationNumber )
111 {
112  G4double Step = (xFinal - xInitial)/iterationNumber ;
113  G4double x = xInitial ;
114  G4double xPlus = xInitial + 0.5*Step ;
115  G4double mean = (fFunction(xInitial) + fFunction(xFinal))*0.5 ;
116  G4double sum = fFunction(xPlus) ;
117  for(G4int i=1;i<iterationNumber;i++)
118  {
119  x += Step ;
120  xPlus += Step ;
121  mean += fFunction(x) ;
122  sum += fFunction(xPlus) ;
123  }
124  mean += 2.0*sum ;
125  return mean*Step/3.0 ;
126 }
127 
128 
129 
130  // Adaptive Gauss integration
131 
132 G4double
134  G4double xFinal )
135 {
136  G4int depth = 0 ;
137  G4double sum = 0.0 ;
138  AdaptGauss(xInitial,xFinal,sum,depth) ;
139  return sum ;
140 }
141 
142 
143 G4double
145  G4double xFinal )
146 {
147  static const G4double root = 1.0/std::sqrt(3.0) ;
148 
149  G4double xMean = (xInitial + xFinal)/2.0 ;
150  G4double Step = (xFinal - xInitial)/2.0 ;
151  G4double delta = Step*root ;
152  G4double sum = (fFunction(xMean + delta) + fFunction(xMean - delta)) ;
153 
154  return sum*Step ;
155 }
156 
157 
158 void
160  G4double xFinal,
161  G4double& sum,
162  G4int& depth )
163 {
164  if(depth >fMaxDepth)
165  {
166  G4Exception("G4SimpleIntegration::AdaptGauss()", "Error",
167  FatalException, "Function varies too rapidly !") ;
168  }
169  G4double xMean = (xInitial + xFinal)/2.0 ;
170  G4double leftHalf = Gauss(xInitial,xMean) ;
171  G4double rightHalf = Gauss(xMean,xFinal) ;
172  G4double full = Gauss(xInitial,xFinal) ;
173  if(std::fabs(leftHalf+rightHalf-full) < fTolerance)
174  {
175  sum += full ;
176  }
177  else
178  {
179  depth++ ;
180  AdaptGauss(xInitial,xMean,sum,depth) ;
181  AdaptGauss(xMean,xFinal,sum,depth) ;
182  }
183 }