G4RKG3_Stepper Class Reference

#include <G4RKG3_Stepper.hh>

Inheritance diagram for G4RKG3_Stepper:

G4MagIntegratorStepper

Public Member Functions

 G4RKG3_Stepper (G4Mag_EqRhs *EqRhs)
 ~G4RKG3_Stepper ()
void Stepper (const G4double yIn[], const G4double dydx[], G4double h, G4double yOut[], G4double yErr[])
G4double DistChord () const
void StepNoErr (const G4double tIn[8], const G4double dydx[6], G4double Step, G4double tOut[8], G4double B[3])
void StepWithEst (const G4double tIn[8], const G4double dydx[6], G4double Step, G4double tOut[8], G4double &alpha2, G4double &beta2, const G4double B1[3], G4double B2[3])
G4int IntegratorOrder () const

Detailed Description

Definition at line 51 of file G4RKG3_Stepper.hh.


Constructor & Destructor Documentation

G4RKG3_Stepper::G4RKG3_Stepper ( G4Mag_EqRhs EqRhs  ) 

Definition at line 35 of file G4RKG3_Stepper.cc.

00036   : G4MagIntegratorStepper(EqRhs,6), hStep(0.)
00037 {
00038 }

G4RKG3_Stepper::~G4RKG3_Stepper (  ) 

Definition at line 40 of file G4RKG3_Stepper.cc.

00041 {
00042 }


Member Function Documentation

G4double G4RKG3_Stepper::DistChord (  )  const [virtual]

Implements G4MagIntegratorStepper.

Definition at line 197 of file G4RKG3_Stepper.cc.

References G4LineSection::Distline().

00198  {
00199    // Soon: must check whether h/R > 2 pi  !!
00200    //  Method below is good only for < 2 pi
00201    G4double distChord,distLine;
00202    
00203    if (fyInitial != fyFinal) {
00204       distLine= G4LineSection::Distline(fyMidPoint,fyInitial,fyFinal );
00205   
00206         distChord = distLine;
00207    }else{
00208       distChord = (fyMidPoint-fyInitial).mag();
00209    }
00210  
00211   
00212    return distChord;
00213    
00214  }

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

Implements G4MagIntegratorStepper.

Definition at line 95 of file G4RKG3_Stepper.hh.

00095 { return 4; }

void G4RKG3_Stepper::StepNoErr ( const G4double  tIn[8],
const G4double  dydx[6],
G4double  Step,
G4double  tOut[8],
G4double  B[3] 
)

Definition at line 126 of file G4RKG3_Stepper.cc.

References G4EquationOfMotion::EvaluateRhsGivenB(), G4EquationOfMotion::EvaluateRhsReturnB(), and G4MagIntegratorStepper::GetEquationOfMotion().

00132 { 
00133   
00134   //  Copy and edit the routine above, to delete alpha2, beta2, ...
00135    G4double K1[7],K2[7],K3[7],K4[7] ;
00136    G4double tTemp[8], yderiv[6] ;
00137 
00138   // Need Momentum value to give correct values to the coefficients in equation
00139   // Integration on unit velocity, but  tIn[3,4,5] is momentum 
00140    G4double mom,inverse_mom;
00141    G4int i ;
00142    const G4double c1=0.5,c2=0.125,c3=1./6.;
00143   
00144    // GetEquationOfMotion()->EvaluateRhsReturnB(tIn,dydx,B1) ;
00145    // Correction for momentum not a velocity
00146    // Need the protection !!! must be not zero 
00147      mom=std::sqrt(tIn[3]*tIn[3]+tIn[4]*tIn[4]+tIn[5]*tIn[5]); 
00148      inverse_mom=1./mom;    
00149    for(i=0;i<3;i++)
00150    {
00151       K1[i] = Step * dydx[i+3]*inverse_mom;
00152       tTemp[i] = tIn[i] + Step*(c1*tIn[i+3]*inverse_mom + c2*K1[i]) ;
00153       tTemp[i+3] = tIn[i+3] + c1*K1[i]*mom ;
00154      
00155    }
00156     
00157    GetEquationOfMotion()->EvaluateRhsReturnB(tTemp,yderiv,B) ;
00158     
00159       
00160    for(i=0;i<3;i++)
00161    {
00162       K2[i] = Step * yderiv[i+3]*inverse_mom;
00163       tTemp[i+3] = tIn[i+3] + c1*K2[i]*mom ;
00164    }
00165    
00166    //  Given B, calculate yderiv !
00167    GetEquationOfMotion()->EvaluateRhsGivenB(tTemp,B,yderiv) ;  
00168  
00169    for(i=0;i<3;i++)
00170    {
00171       K3[i] = Step * yderiv[i+3]*inverse_mom;
00172       tTemp[i] = tIn[i] + Step*(tIn[i+3]*inverse_mom + c1*K3[i]) ;
00173       tTemp[i+3] = tIn[i+3] + K3[i]*mom ;
00174    }
00175    
00176 
00177    //  Calculates y-deriv(atives) & returns B too!
00178    GetEquationOfMotion()->EvaluateRhsReturnB(tTemp,yderiv,B) ;  
00179     
00180 
00181    for(i=0;i<3;i++)        // Output trajectory vector
00182    {
00183       K4[i] = Step * yderiv[i+3]*inverse_mom;
00184       tOut[i] = tIn[i] + Step*(tIn[i+3]*inverse_mom+ (K1[i] + K2[i] + K3[i])*c3) ;
00185       tOut[i+3] = tIn[i+3] + mom*(K1[i] + 2*K2[i] + 2*K3[i] +K4[i])*c3 ;
00186    }
00187    tOut[6] = tIn[6];
00188    tOut[7] = tIn[7];
00189    // NormaliseTangentVector( tOut );
00190   
00191 
00192 }

void G4RKG3_Stepper::Stepper ( const G4double  yIn[],
const G4double  dydx[],
G4double  h,
G4double  yOut[],
G4double  yErr[] 
) [virtual]

Implements G4MagIntegratorStepper.

void G4RKG3_Stepper::StepWithEst ( const G4double  tIn[8],
const G4double  dydx[6],
G4double  Step,
G4double  tOut[8],
G4double alpha2,
G4double beta2,
const G4double  B1[3],
G4double  B2[3] 
)


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