Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions
G4ConstRK4 Class Reference

#include <G4ConstRK4.hh>

Inheritance diagram for G4ConstRK4:
G4MagErrorStepper G4MagIntegratorStepper

Public Member Functions

 G4ConstRK4 (G4Mag_EqRhs *EquationMotion, G4int numberOfStateVariables=8)
 
 ~G4ConstRK4 ()
 
void Stepper (const G4double y[], const G4double dydx[], G4double h, G4double yout[], G4double yerr[])
 
void DumbStepper (const G4double yIn[], const G4double dydx[], G4double h, G4double yOut[])
 
G4double DistChord () const
 
void RightHandSideConst (const G4double y[], G4double dydx[]) const
 
void GetConstField (const G4double y[], G4double Field[])
 
G4int IntegratorOrder () const
 
- Public Member Functions inherited from G4MagErrorStepper
 G4MagErrorStepper (G4EquationOfMotion *EqRhs, G4int numberOfVariables, G4int numStateVariables=12)
 
virtual ~G4MagErrorStepper ()
 
void Stepper (const G4double y[], const G4double dydx[], G4double h, G4double yout[], G4double yerr[])
 
G4double DistChord () const
 
- Public Member Functions inherited from G4MagIntegratorStepper
 G4MagIntegratorStepper (G4EquationOfMotion *Equation, G4int numIntegrationVariables, G4int numStateVariables=12)
 
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
 
G4EquationOfMotion * GetEquationOfMotion ()
 
void SetEquationOfMotion (G4EquationOfMotion *newEquation)
 

Detailed Description

Definition at line 51 of file G4ConstRK4.hh.

Constructor & Destructor Documentation

G4ConstRK4::G4ConstRK4 ( G4Mag_EqRhs *  EquationMotion,
G4int  numberOfStateVariables = 8 
)

Definition at line 42 of file G4ConstRK4.cc.

References FatalException, G4endl, and G4Exception().

43  : G4MagErrorStepper(EqRhs, 6, numStateVariables)
44 {
45  // const G4int numberOfVariables= 6;
46  if( numStateVariables < 8 )
47  {
48  std::ostringstream message;
49  message << "The number of State variables at least 8 " << G4endl
50  << "Instead it is - numStateVariables= " << numStateVariables;
51  G4Exception("G4ConstRK4::G4ConstRK4()", "GeomField0002",
52  FatalException, message, "Use another Stepper!");
53  }
54 
55  fEq = EqRhs;
56  yMiddle = new G4double[8];
57  dydxMid = new G4double[8];
58  yInitial = new G4double[8];
59  yOneStep = new G4double[8];
60 
61  dydxm = new G4double[8];
62  dydxt = new G4double[8];
63  yt = new G4double[8];
64  Field[0]=0.; Field[1]=0.; Field[2]=0.;
65 }
G4MagErrorStepper::G4MagErrorStepper
G4MagErrorStepper(G4EquationOfMotion *EqRhs, G4int numberOfVariables, G4int numStateVariables=12)
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
FatalException
Definition: G4ExceptionSeverity.hh:60
G4endl
#define G4endl
Definition: G4ios.hh:61
G4double
double G4double
Definition: G4Types.hh:76
G4ConstRK4::~G4ConstRK4 ( )

Definition at line 71 of file G4ConstRK4.cc.

72 {
73  delete [] yMiddle;
74  delete [] dydxMid;
75  delete [] yInitial;
76  delete [] yOneStep;
77  delete [] dydxm;
78  delete [] dydxt;
79  delete [] yt;
80 }

Member Function Documentation

G4double G4ConstRK4::DistChord ( ) const
virtual

Implements G4MagIntegratorStepper.

Definition at line 213 of file G4ConstRK4.cc.

References G4LineSection::Distline().

214 {
215  G4double distLine, distChord;
216 
217  if (fInitialPoint != fFinalPoint)
218  {
219  distLine= G4LineSection::Distline( fMidPoint, fInitialPoint, fFinalPoint );
220  // This is a class method that gives distance of Mid
221  // from the Chord between the Initial and Final points
222  distChord = distLine;
223  }
224  else
225  {
226  distChord = (fMidPoint-fInitialPoint).mag();
227  }
228  return distChord;
229 }
G4LineSection::Distline
static G4double Distline(const G4ThreeVector &OtherPnt, const G4ThreeVector &LinePntA, const G4ThreeVector &LinePntB)
Definition: G4LineSection.cc:101
G4double
double G4double
Definition: G4Types.hh:76
void G4ConstRK4::DumbStepper ( const G4double  yIn[],
const G4double  dydx[],
G4double  h,
G4double  yOut[] 
)
virtual

Implements G4MagErrorStepper.

Definition at line 92 of file G4ConstRK4.cc.

References RightHandSideConst().

Referenced by Stepper().

96 {
97  G4double hh = h*0.5 , h6 = h/6.0 ;
98 
99  // 1st Step K1=h*dydx
100  yt[5] = yIn[5] + hh*dydx[5] ;
101  yt[4] = yIn[4] + hh*dydx[4] ;
102  yt[3] = yIn[3] + hh*dydx[3] ;
103  yt[2] = yIn[2] + hh*dydx[2] ;
104  yt[1] = yIn[1] + hh*dydx[1] ;
105  yt[0] = yIn[0] + hh*dydx[0] ;
106  RightHandSideConst(yt,dydxt) ;
107 
108  // 2nd Step K2=h*dydxt
109  yt[5] = yIn[5] + hh*dydxt[5] ;
110  yt[4] = yIn[4] + hh*dydxt[4] ;
111  yt[3] = yIn[3] + hh*dydxt[3] ;
112  yt[2] = yIn[2] + hh*dydxt[2] ;
113  yt[1] = yIn[1] + hh*dydxt[1] ;
114  yt[0] = yIn[0] + hh*dydxt[0] ;
115  RightHandSideConst(yt,dydxm) ;
116 
117  // 3rd Step K3=h*dydxm
118  // now dydxm=(K2+K3)/h
119  yt[5] = yIn[5] + h*dydxm[5] ;
120  dydxm[5] += dydxt[5] ;
121  yt[4] = yIn[4] + h*dydxm[4] ;
122  dydxm[4] += dydxt[4] ;
123  yt[3] = yIn[3] + h*dydxm[3] ;
124  dydxm[3] += dydxt[3] ;
125  yt[2] = yIn[2] + h*dydxm[2] ;
126  dydxm[2] += dydxt[2] ;
127  yt[1] = yIn[1] + h*dydxm[1] ;
128  dydxm[1] += dydxt[1] ;
129  yt[0] = yIn[0] + h*dydxm[0] ;
130  dydxm[0] += dydxt[0] ;
131  RightHandSideConst(yt,dydxt) ;
132 
133  // 4th Step K4=h*dydxt
134  yOut[5] = yIn[5]+h6*(dydx[5]+dydxt[5]+2.0*dydxm[5]);
135  yOut[4] = yIn[4]+h6*(dydx[4]+dydxt[4]+2.0*dydxm[4]);
136  yOut[3] = yIn[3]+h6*(dydx[3]+dydxt[3]+2.0*dydxm[3]);
137  yOut[2] = yIn[2]+h6*(dydx[2]+dydxt[2]+2.0*dydxm[2]);
138  yOut[1] = yIn[1]+h6*(dydx[1]+dydxt[1]+2.0*dydxm[1]);
139  yOut[0] = yIn[0]+h6*(dydx[0]+dydxt[0]+2.0*dydxm[0]);
140 
141 } // end of DumbStepper ....................................................
G4ConstRK4::RightHandSideConst
void RightHandSideConst(const G4double y[], G4double dydx[]) const
Definition: G4ConstRK4.hh:96
G4double
double G4double
Definition: G4Types.hh:76
void G4ConstRK4::GetConstField ( const G4double  y[],
G4double  Field[] 
)
inline

Definition at line 114 of file G4ConstRK4.hh.

Referenced by Stepper().

115 {
116  G4double PositionAndTime[4];
117 
118  PositionAndTime[0] = y[0];
119  PositionAndTime[1] = y[1];
120  PositionAndTime[2] = y[2];
121  // Global Time
122  PositionAndTime[3] = y[7];
123  fEq -> GetFieldValue(PositionAndTime, B) ;
124 }
G4double
double G4double
Definition: G4Types.hh:76
G4int G4ConstRK4::IntegratorOrder ( ) const
inlinevirtual

Implements G4MagIntegratorStepper.

Definition at line 76 of file G4ConstRK4.hh.

Referenced by Stepper().

76 { return 4; }
void G4ConstRK4::RightHandSideConst ( const G4double  y[],
G4double  dydx[] 
) const
inline

Definition at line 96 of file G4ConstRK4.hh.

References G4Mag_EqRhs::FCof().

Referenced by DumbStepper(), and Stepper().

98 {
99 
100  G4double momentum_mag_square = y[3]*y[3] + y[4]*y[4] + y[5]*y[5];
101  G4double inv_momentum_magnitude = 1.0 / std::sqrt( momentum_mag_square );
102 
103  G4double cof =fEq->FCof()*inv_momentum_magnitude;
104 
105  dydx[0] = y[3]*inv_momentum_magnitude; // (d/ds)x = Vx/V
106  dydx[1] = y[4]*inv_momentum_magnitude; // (d/ds)y = Vy/V
107  dydx[2] = y[5]*inv_momentum_magnitude; // (d/ds)z = Vz/V
108 
109  dydx[3] = cof*(y[4]*Field[2] - y[5]*Field[1]) ; // Ax = a*(Vy*Bz - Vz*By)
110  dydx[4] = cof*(y[5]*Field[0] - y[3]*Field[2]) ; // Ay = a*(Vz*Bx - Vx*Bz)
111  dydx[5] = cof*(y[3]*Field[1] - y[4]*Field[0]) ; // Az = a*(Vx*By - Vy*Bx)
112 }
G4Mag_EqRhs::FCof
G4double FCof() const
Definition: G4Mag_EqRhs.hh:84
G4double
double G4double
Definition: G4Types.hh:76
void G4ConstRK4::Stepper ( const G4double  y[],
const G4double  dydx[],
G4double  h,
G4double  yout[],
G4double  yerr[] 
)
virtual

Implements G4MagIntegratorStepper.

Definition at line 148 of file G4ConstRK4.cc.

References DumbStepper(), GetConstField(), G4MagIntegratorStepper::GetNumberOfStateVariables(), IntegratorOrder(), and RightHandSideConst().

153 {
154  const G4int nvar = 6; // number of variables integrated
155  const G4int maxvar= GetNumberOfStateVariables();
156 
157  // Correction for Richardson extrapolation
158  G4double correction = 1. / ( (1 << IntegratorOrder()) -1 );
159 
160  G4int i;
161 
162  // Saving yInput because yInput and yOutput can be aliases for same array
163  for (i=0; i<maxvar; i++) { yInitial[i]= yInput[i]; }
164 
165  // Must copy the part of the state *not* integrated to the output
166  for (i=nvar; i<maxvar; i++) { yOutput[i]= yInput[i]; }
167 
168  // yInitial[7]= yInput[7]; // The time is typically needed
169  yMiddle[7] = yInput[7]; // Copy the time from initial value
170  yOneStep[7] = yInput[7]; // As it contributes to final value of yOutput ?
171  // yOutput[7] = yInput[7]; // -> dumb stepper does it too for RK4
172  yError[7] = 0.0;
173 
174  G4double halfStep = hstep * 0.5;
175 
176  // Do two half steps
177  //
178  GetConstField(yInitial,Field);
179  DumbStepper (yInitial, dydx, halfStep, yMiddle);
180  RightHandSideConst(yMiddle, dydxMid);
181  DumbStepper (yMiddle, dydxMid, halfStep, yOutput);
182 
183  // Store midpoint, chord calculation
184  //
185  fMidPoint = G4ThreeVector( yMiddle[0], yMiddle[1], yMiddle[2]);
186 
187  // Do a full Step
188  //
189  DumbStepper(yInitial, dydx, hstep, yOneStep);
190  for(i=0;i<nvar;i++)
191  {
192  yError [i] = yOutput[i] - yOneStep[i] ;
193  yOutput[i] += yError[i]*correction ;
194  // Provides accuracy increased by 1 order via the
195  // Richardson extrapolation
196  }
197 
198  fInitialPoint = G4ThreeVector( yInitial[0], yInitial[1], yInitial[2]);
199  fFinalPoint = G4ThreeVector( yOutput[0], yOutput[1], yOutput[2]);
200 
201  return;
202 }
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:42
G4int
int G4int
Definition: G4Types.hh:78
G4ConstRK4::GetConstField
void GetConstField(const G4double y[], G4double Field[])
Definition: G4ConstRK4.hh:114
G4ConstRK4::RightHandSideConst
void RightHandSideConst(const G4double y[], G4double dydx[]) const
Definition: G4ConstRK4.hh:96
G4MagIntegratorStepper::GetNumberOfStateVariables
G4int GetNumberOfStateVariables() const
G4ConstRK4::IntegratorOrder
G4int IntegratorOrder() const
Definition: G4ConstRK4.hh:76
G4double
double G4double
Definition: G4Types.hh:76
G4ConstRK4::DumbStepper
void DumbStepper(const G4double yIn[], const G4double dydx[], G4double h, G4double yOut[])
Definition: G4ConstRK4.cc:92
The documentation for this class was generated from the following files:
Generated on Wed Apr 30 2014 15:56:17 for Geant4.10 by   doxygen 1.8.7