G4ClassicalRK4 Class Reference

#include <G4ClassicalRK4.hh>

Inheritance diagram for G4ClassicalRK4:

G4MagErrorStepper G4MagIntegratorStepper

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

Detailed Description

Definition at line 44 of file G4ClassicalRK4.hh.


Constructor & Destructor Documentation

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

Definition at line 39 of file G4ClassicalRK4.cc.

00040   : G4MagErrorStepper(EqRhs, numberOfVariables)
00041 {
00042    unsigned int noVariables= std::max(numberOfVariables,8); // For Time .. 7+1
00043  
00044    dydxm = new G4double[noVariables];
00045    dydxt = new G4double[noVariables]; 
00046    yt    = new G4double[noVariables]; 
00047 }

G4ClassicalRK4::~G4ClassicalRK4 (  ) 

Definition at line 53 of file G4ClassicalRK4.cc.

00054 {
00055   delete[] dydxm;
00056   delete[] dydxt;
00057   delete[] yt;
00058 }


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.

References G4MagIntegratorStepper::GetNumberOfVariables(), G4MagIntegratorStepper::NormalisePolarizationVector(), and G4MagIntegratorStepper::RightHandSide().

00075 {
00076   const G4int nvar = this->GetNumberOfVariables();   //  fNumberOfVariables(); 
00077   G4int i;
00078   G4double  hh = h*0.5 , h6 = h/6.0  ;
00079 
00080   // Initialise time to t0, needed when it is not updated by the integration.
00081   //        [ Note: Only for time dependent fields (usually electric) 
00082   //                  is it neccessary to integrate the time.] 
00083   yt[7]   = yIn[7]; 
00084   yOut[7] = yIn[7];
00085 
00086   for(i=0;i<nvar;i++)
00087   {
00088     yt[i] = yIn[i] + hh*dydx[i] ;             // 1st Step K1=h*dydx
00089   }
00090   RightHandSide(yt,dydxt) ;                   // 2nd Step K2=h*dydxt
00091 
00092   for(i=0;i<nvar;i++)
00093   { 
00094     yt[i] = yIn[i] + hh*dydxt[i] ;
00095   }
00096   RightHandSide(yt,dydxm) ;                   // 3rd Step K3=h*dydxm
00097 
00098   for(i=0;i<nvar;i++)
00099   {
00100     yt[i]   = yIn[i] + h*dydxm[i] ;
00101     dydxm[i] += dydxt[i] ;                    // now dydxm=(K2+K3)/h
00102   }
00103   RightHandSide(yt,dydxt) ;                   // 4th Step K4=h*dydxt
00104  
00105   for(i=0;i<nvar;i++)    // Final RK4 output
00106   {
00107     yOut[i] = yIn[i]+h6*(dydx[i]+dydxt[i]+2.0*dydxm[i]); //+K1/6+K4/6+(K2+K3)/3
00108   }
00109   if ( nvar == 12 )  { NormalisePolarizationVector ( yOut ); }
00110   
00111 }  // end of DumbStepper ....................................................

G4int G4ClassicalRK4::IntegratorOrder (  )  const [inline, virtual]

Implements G4MagIntegratorStepper.

Definition at line 69 of file G4ClassicalRK4.hh.

00069 { return 4; }


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:51:38 2013 for Geant4 by  doxygen 1.4.7