Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4DoLoMcPriRK34 Class Reference

#include <G4DoLoMcPriRK34.hh>

Inheritance diagram for G4DoLoMcPriRK34:
Collaboration diagram for G4DoLoMcPriRK34:

Public Member Functions

 G4DoLoMcPriRK34 (G4EquationOfMotion *EqRhs, G4int numberOfVariables=6, G4bool primary=true)
 
 ~G4DoLoMcPriRK34 ()
 
void Stepper (const G4double y[], const G4double dydx[], G4double h, G4double yout[], G4double yerr[])
 
void SetupInterpolation ()
 
void SetupInterpolate (const G4double yInput[], const G4double dydx[], const G4double Step)
 
void Interpolate (const G4double yInput[], const G4double dydx[], const G4double Step, G4double yOut[], G4double tau)
 
void Interpolate (G4double tau, G4double yOut[])
 
void interpolate (const G4double yInput[], const G4double dydx[], G4double yOut[], G4double Step, G4double tau)
 
G4double DistChord () const
 
G4int IntegratorOrder () const
 
- Public Member Functions inherited from G4MagIntegratorStepper
 G4MagIntegratorStepper (G4EquationOfMotion *Equation, G4int numIntegrationVariables, G4int numStateVariables=12, bool isFSAL=false)
 
virtual ~G4MagIntegratorStepper ()
 
virtual void ComputeRightHandSide (const G4double y[], G4double dydx[])
 
void NormaliseTangentVector (G4double vec[6])
 
void NormalisePolarizationVector (G4double vec[12])
 
void RightHandSide (const double y[], double dydx[])
 
G4int GetNumberOfVariables () const
 
G4int GetNumberOfStateVariables () const
 
G4int IntegrationOrder ()
 
G4EquationOfMotionGetEquationOfMotion ()
 
void SetEquationOfMotion (G4EquationOfMotion *newEquation)
 
unsigned long GetfNoRHSCalls ()
 
void ResetfNORHSCalls ()
 
bool IsFSAL ()
 

Additional Inherited Members

- Protected Member Functions inherited from G4MagIntegratorStepper
void SetIntegrationOrder (int order)
 
void SetFSAL (bool flag=true)
 

Detailed Description

Definition at line 55 of file G4DoLoMcPriRK34.hh.

Constructor & Destructor Documentation

G4DoLoMcPriRK34::G4DoLoMcPriRK34 ( G4EquationOfMotion EqRhs,
G4int  numberOfVariables = 6,
G4bool  primary = true 
)

Definition at line 53 of file G4DoLoMcPriRK34.cc.

56  : G4MagIntegratorStepper(EqRhs, noIntegrationVariables),
57  fAuxStepper(0)
58 {
59  const G4int numberOfVariables = noIntegrationVariables;
60 
61  //New Chunk of memory being created for use by the Stepper
62 
63  //aki - for storing intermediate RHS
64  ak2 = new G4double[numberOfVariables];
65  ak3 = new G4double[numberOfVariables];
66  ak4 = new G4double[numberOfVariables];
67  ak5 = new G4double[numberOfVariables];
68  ak6 = new G4double[numberOfVariables];
69 
70  yTemp = new G4double[numberOfVariables] ;
71  yIn = new G4double[numberOfVariables] ;
72 
73  fLastInitialVector = new G4double[numberOfVariables] ;
74  fLastFinalVector = new G4double[numberOfVariables] ;
75  fLastDyDx = new G4double[numberOfVariables];
76 
77  fMidVector = new G4double[numberOfVariables];
78  fMidError = new G4double[numberOfVariables];
79  if( primary )
80  {
81  fAuxStepper = new G4DoLoMcPriRK34(EqRhs, numberOfVariables,
82  !primary);
83  }
84 }
int G4int
Definition: G4Types.hh:78
G4MagIntegratorStepper(G4EquationOfMotion *Equation, G4int numIntegrationVariables, G4int numStateVariables=12, bool isFSAL=false)
double G4double
Definition: G4Types.hh:76
G4DoLoMcPriRK34(G4EquationOfMotion *EqRhs, G4int numberOfVariables=6, G4bool primary=true)
G4DoLoMcPriRK34::~G4DoLoMcPriRK34 ( )

Definition at line 88 of file G4DoLoMcPriRK34.cc.

89 {
90  //clear all previously allocated memory for Stepper and DistChord
91  delete[] ak2;
92  delete[] ak3;
93  delete[] ak4;
94  delete[] ak5;
95  delete[] ak6;
96 
97  delete[] yTemp;
98  delete[] yIn;
99 
100  delete[] fLastInitialVector;
101  delete[] fLastFinalVector;
102  delete[] fLastDyDx;
103  delete[] fMidVector;
104  delete[] fMidError;
105 
106  delete fAuxStepper;
107 
108 
109 }

Member Function Documentation

G4double G4DoLoMcPriRK34::DistChord ( ) const
virtual

Implements G4MagIntegratorStepper.

Definition at line 233 of file G4DoLoMcPriRK34.cc.

234 {
235  G4double distLine, distChord;
236  G4ThreeVector initialPoint, finalPoint, midPoint;
237 
238  // Store last initial and final points (they will be overwritten in self-Stepper call!)
239  initialPoint = G4ThreeVector( fLastInitialVector[0],
240  fLastInitialVector[1], fLastInitialVector[2]);
241  finalPoint = G4ThreeVector( fLastFinalVector[0],
242  fLastFinalVector[1], fLastFinalVector[2]);
243 
244  // Do half a Step using StepNoErr
245 
246  fAuxStepper->Stepper( fLastInitialVector, fLastDyDx, 0.5 * fLastStepLength,
247  fMidVector, fMidError );
248 
249  midPoint = G4ThreeVector( fMidVector[0], fMidVector[1], fMidVector[2]);
250 
251  // Use stored values of Initial and Endpoint + new Midpoint to evaluate
252  // distance of Chord
253 
254 
255  if (initialPoint != finalPoint)
256  {
257  distLine = G4LineSection::Distline( midPoint, initialPoint, finalPoint );
258  distChord = distLine;
259  }
260  else
261  {
262  distChord = (midPoint-initialPoint).mag();
263  }
264  return distChord;
265 }
CLHEP::Hep3Vector G4ThreeVector
static G4double Distline(const G4ThreeVector &OtherPnt, const G4ThreeVector &LinePntA, const G4ThreeVector &LinePntB)
void Stepper(const G4double y[], const G4double dydx[], G4double h, G4double yout[], G4double yerr[])
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4int G4DoLoMcPriRK34::IntegratorOrder ( ) const
inlinevirtual

Implements G4MagIntegratorStepper.

Definition at line 93 of file G4DoLoMcPriRK34.hh.

93 {return 3; }
void G4DoLoMcPriRK34::Interpolate ( const G4double  yInput[],
const G4double  dydx[],
const G4double  Step,
G4double  yOut[],
G4double  tau 
)

Definition at line 285 of file G4DoLoMcPriRK34.cc.

289  {
290  G4double
291  bf1, bf2, bf3, bf4, bf5, bf6;
292 
293 
294  const G4int numberOfVariables= this->GetNumberOfVariables();
295 
296  for(int i=0;i<numberOfVariables;i++)
297  {
298  yIn[i]=yInput[i];
299  }
300 
301  G4double
302  tau_2 = tau*tau ,
303  tau_3 = tau*tau_2;
304 
305  //Calculating the polynomials (coefficients for the respective stages)
306  bf1 = -(162.0*tau_3 - 504.0*tau_2 + 551.0*tau - 238.0)/238.0 ,
307  bf2 = 0.0 ,
308  bf3 = 27.0*tau*(27.0*tau_2 - 70.0*tau + 51.0 )/385.0 ,
309  bf4 = -27*tau*(27.0*tau_2 - 50.0*tau + 21.0)/85.0 ,
310  bf5 = 7.0*tau*(2232.0*tau_2 - 4166.0*tau + 1785.0 )/3278.0 ,
311  bf6 = tau*(tau - 1.0)*(387.0*tau - 238.0)/149.0 ;
312 
313  for( int i=0; i<numberOfVariables; i++){
314  yOut[i] = yIn[i] + Step*tau*(bf1*dydx[i] + bf2*ak2[i] + bf3*ak3[i] +
315  bf4*ak4[i] + bf5*ak5[i] + bf6*ak6[i] ) ;
316  }
317 
318 
319 
320 }
int G4int
Definition: G4Types.hh:78
G4int GetNumberOfVariables() const
Definition: Step.hh:41
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4DoLoMcPriRK34::Interpolate ( G4double  tau,
G4double  yOut[] 
)

Definition at line 278 of file G4DoLoMcPriRK34.cc.

280 {
281  Interpolate( fLastInitialVector, fLastDyDx, fLastStepLength, yOut, tau );
282 }
void Interpolate(const G4double yInput[], const G4double dydx[], const G4double Step, G4double yOut[], G4double tau)

Here is the call graph for this function:

void G4DoLoMcPriRK34::interpolate ( const G4double  yInput[],
const G4double  dydx[],
G4double  yOut[],
G4double  Step,
G4double  tau 
)
void G4DoLoMcPriRK34::SetupInterpolate ( const G4double  yInput[],
const G4double  dydx[],
const G4double  Step 
)

Definition at line 270 of file G4DoLoMcPriRK34.cc.

273 {
274  //Do Nothing
275 }
void G4DoLoMcPriRK34::SetupInterpolation ( )

Definition at line 267 of file G4DoLoMcPriRK34.cc.

268 {}
void G4DoLoMcPriRK34::Stepper ( const G4double  y[],
const G4double  dydx[],
G4double  h,
G4double  yout[],
G4double  yerr[] 
)
virtual

Implements G4MagIntegratorStepper.

Definition at line 117 of file G4DoLoMcPriRK34.cc.

122 {
123  G4int i;
124 
125  //The various constants defined on the basis of butcher tableu
126  const G4double //G4double - only once
127 
128  b21 = 7.0/27.0 ,
129 
130 
131  b31 = 7.0/72.0 ,
132  b32 = 7.0/24.0 ,
133 
134  b41 = 3043.0/3528.0 ,
135  b42 = -3757.0/1176.0 ,
136  b43 = 1445.0/441.0,
137 
138  b51 = 17617.0/11662.0 ,
139  b52 = -4023.0/686.0 ,
140  b53 = 9372.0/1715.0 ,
141  b54 = -66.0/595.0 ,
142 
143  b61 = 29.0/238.0 ,
144  b62 = 0.0 ,
145  b63 = 216.0/385.0 ,
146  b64 = 54.0/85.0 ,
147  b65 = -7.0/22.0 ,
148 
149 
150 
151  dc1 = 363.0/2975.0 - b61 ,
152  dc2 = 0.0 - b62 ,
153  dc3 = 981.0/1750.0 - b63,
154  dc4 = 2709.0/4250.0 - b64 ,
155  dc5 = -3.0/10.0 - b65 ,
156  dc6 = -1.0/50.0 ; //end of declaration
157 
158 
159  const G4int numberOfVariables= this->GetNumberOfVariables();
160 
161  // The number of variables to be integrated over
162  yOut[7] = yTemp[7] = yIn[7];
163  // Saving yInput because yInput and yOut can be aliases for same array
164 
165  for(i=0;i<numberOfVariables;i++)
166  {
167  yIn[i]=yInput[i];
168  }
169 
170 
171 
172  // RightHandSide(yIn, DyDx) ;
173  // 1st stage - Not doing, getting passed
174 
175  for(i=0;i<numberOfVariables;i++)
176  {
177  yTemp[i] = yIn[i] + b21*Step*DyDx[i] ;
178  }
179  RightHandSide(yTemp, ak2) ; // 2nd stage
180 
181  for(i=0;i<numberOfVariables;i++)
182  {
183  yTemp[i] = yIn[i] + Step*(b31*DyDx[i] + b32*ak2[i]) ;
184  }
185  RightHandSide(yTemp, ak3) ; // 3rd stage
186 
187  for(i=0;i<numberOfVariables;i++)
188  {
189  yTemp[i] = yIn[i] + Step*(b41*DyDx[i] + b42*ak2[i] + b43*ak3[i]) ;
190  }
191  RightHandSide(yTemp, ak4) ; // 4th stage
192 
193  for(i=0;i<numberOfVariables;i++)
194  {
195  yTemp[i] = yIn[i] + Step*(b51*DyDx[i] + b52*ak2[i] + b53*ak3[i] +
196  b54*ak4[i]) ;
197  }
198  RightHandSide(yTemp, ak5) ; // 5th stage
199 
200  for(i=0;i<numberOfVariables;i++)
201  {
202  yOut[i] = yIn[i] + Step*(b61*DyDx[i] + b62*ak2[i] + b63*ak3[i] +
203  b64*ak4[i] + b65*ak5[i]) ;
204  }
205  RightHandSide(yOut, ak6) ; // 6th and Final stage
206 
207 
208 
209  for(i=0;i<numberOfVariables;i++)
210  {
211 
212  yErr[i] = Step*(dc1*DyDx[i] + dc2*ak2[i] + dc3*ak3[i] + dc4*ak4[i] +
213  dc5*ak5[i] + dc6*ak6[i] ) ;
214 
215 
216  // Store Input and Final values, for possible use in calculating chord
217  fLastInitialVector[i] = yIn[i] ;
218  fLastFinalVector[i] = yOut[i];
219  fLastDyDx[i] = DyDx[i];
220 
221 
222  }
223 
224  fLastStepLength = Step;
225 
226  return ;
227 }
int G4int
Definition: G4Types.hh:78
G4int GetNumberOfVariables() const
Definition: Step.hh:41
void RightHandSide(const double y[], double dydx[])
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:


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