Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4ClassicalRK4 Class Reference

#include <G4ClassicalRK4.hh>

Inheritance diagram for G4ClassicalRK4:
Collaboration diagram for G4ClassicalRK4:

Public Member Functions

 G4ClassicalRK4 (G4EquationOfMotion *EquationMotion, G4int numberOfVariables=6)
 
 ~G4ClassicalRK4 ()
 
void DumbStepper (const G4double yIn[], const G4double dydx[], G4double h, G4double yOut[])
 
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, 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 48 of file G4ClassicalRK4.hh.

Constructor & Destructor Documentation

G4ClassicalRK4::G4ClassicalRK4 ( G4EquationOfMotion EquationMotion,
G4int  numberOfVariables = 6 
)

Definition at line 39 of file G4ClassicalRK4.cc.

40  : G4MagErrorStepper(EqRhs, numberOfVariables)
41 {
42  unsigned int noVariables= std::max(numberOfVariables,8); // For Time .. 7+1
43 
44  dydxm = new G4double[noVariables];
45  dydxt = new G4double[noVariables];
46  yt = new G4double[noVariables];
47 }
G4MagErrorStepper(G4EquationOfMotion *EqRhs, G4int numberOfVariables, G4int numStateVariables=12)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4ClassicalRK4::~G4ClassicalRK4 ( )

Definition at line 53 of file G4ClassicalRK4.cc.

54 {
55  delete[] dydxm;
56  delete[] dydxt;
57  delete[] yt;
58 }

Member Function Documentation

void G4ClassicalRK4::DumbStepper ( const G4double  yIn[],
const G4double  dydx[],
G4double  h,
G4double  yOut[] 
)
virtual

Implements G4MagErrorStepper.

Definition at line 71 of file G4ClassicalRK4.cc.

75 {
76  const G4int nvar = this->GetNumberOfVariables(); // fNumberOfVariables();
77  G4int i;
78  G4double hh = h*0.5 , h6 = h/6.0 ;
79 
80  // Initialise time to t0, needed when it is not updated by the integration.
81  // [ Note: Only for time dependent fields (usually electric)
82  // is it neccessary to integrate the time.]
83  yt[7] = yIn[7];
84  yOut[7] = yIn[7];
85 
86  for(i=0;i<nvar;i++)
87  {
88  yt[i] = yIn[i] + hh*dydx[i] ; // 1st Step K1=h*dydx
89  }
90  RightHandSide(yt,dydxt) ; // 2nd Step K2=h*dydxt
91 
92  for(i=0;i<nvar;i++)
93  {
94  yt[i] = yIn[i] + hh*dydxt[i] ;
95  }
96  RightHandSide(yt,dydxm) ; // 3rd Step K3=h*dydxm
97 
98  for(i=0;i<nvar;i++)
99  {
100  yt[i] = yIn[i] + h*dydxm[i] ;
101  dydxm[i] += dydxt[i] ; // now dydxm=(K2+K3)/h
102  }
103  RightHandSide(yt,dydxt) ; // 4th Step K4=h*dydxt
104 
105  for(i=0;i<nvar;i++) // Final RK4 output
106  {
107  yOut[i] = yIn[i]+h6*(dydx[i]+dydxt[i]+2.0*dydxm[i]); //+K1/6+K4/6+(K2+K3)/3
108  }
109  if ( nvar == 12 ) { NormalisePolarizationVector ( yOut ); }
110 
111 } // end of DumbStepper ....................................................
int G4int
Definition: G4Types.hh:78
G4int GetNumberOfVariables() const
void RightHandSide(const double y[], double dydx[])
void NormalisePolarizationVector(G4double vec[12])
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4int G4ClassicalRK4::IntegratorOrder ( ) const
inlinevirtual

Implements G4MagIntegratorStepper.

Definition at line 73 of file G4ClassicalRK4.hh.

73 { return 4; }

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