G4MagHelicalStepper Class Reference

#include <G4MagHelicalStepper.hh>

Inheritance diagram for G4MagHelicalStepper:

Collaboration diagram for G4MagHelicalStepper:

Public Member Functions
	G4MagHelicalStepper (G4Mag_EqRhs *EqRhs)
virtual	~G4MagHelicalStepper ()
virtual void	Stepper (const G4double y[], const G4double dydx[], G4double h, G4double yout[], G4double yerr[])
virtual void	DumbStepper (const G4double y[], G4ThreeVector Bfld, G4double h, G4double yout[])=0
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
virtual G4int	IntegratorOrder () const =0
G4EquationOfMotion *	GetEquationOfMotion ()
void	SetEquationOfMotion (G4EquationOfMotion *newEquation)

Protected Member Functions
void	LinearStep (const G4double yIn[], G4double h, G4double yHelix[]) const
void	AdvanceHelix (const G4double yIn[], G4ThreeVector Bfld, G4double h, G4double yHelix[], G4double yHelix2[]=0)
void	MagFieldEvaluate (const G4double y[], G4ThreeVector &Bfield)
G4double	GetInverseCurve (const G4double Momentum, const G4double Bmag)
void	SetAngCurve (const G4double Ang)
G4double	GetAngCurve () const
void	SetCurve (const G4double Curve)
G4double	GetCurve () const
void	SetRadHelix (const G4double Rad)
G4double	GetRadHelix () const

Private Member Functions
	G4MagHelicalStepper (const G4MagHelicalStepper &)
G4MagHelicalStepper &	operator= (const G4MagHelicalStepper &)

Private Attributes
G4Mag_EqRhs *	fPtrMagEqOfMot
G4double	fAngCurve
G4double	frCurve
G4double	frHelix
G4ThreeVector	yInitial
G4ThreeVector	yMidPoint
G4ThreeVector	yFinal

Static Private Attributes
static const G4double	fUnitConstant = 0.299792458*(GeV/(tesla*m))

Detailed Description

Definition at line 58 of file G4MagHelicalStepper.hh.

Constructor & Destructor Documentation

G4MagHelicalStepper() [1/2]

G4MagHelicalStepper::G4MagHelicalStepper

(

G4Mag_EqRhs *

EqRhs

)

Definition at line 47 of file G4MagHelicalStepper.cc.

   : G4MagIntegratorStepper(EqRhs, 6), // integrate over 6 variables only !!
                                       // position & velocity
     fPtrMagEqOfMot(EqRhs), fAngCurve(0.), frCurve(0.), frHelix(0.)
{
}

~G4MagHelicalStepper()

G4MagHelicalStepper::~G4MagHelicalStepper ( )

virtual

Definition at line 54 of file G4MagHelicalStepper.cc.

{
}

G4MagHelicalStepper() [2/2]

G4MagHelicalStepper::G4MagHelicalStepper ( const G4MagHelicalStepper &  )

private

Member Function Documentation

AdvanceHelix()

void G4MagHelicalStepper::AdvanceHelix ( const G4double  yIn[], G4ThreeVector  Bfld, G4double  h, G4double  yHelix[], G4double  yHelix2[] = 0  )

protected

Definition at line 59 of file G4MagHelicalStepper.cc.

{
  // const G4int    nvar = 6;
 
  // OLD  const G4double approc_limit = 0.05;
  // OLD  approc_limit = 0.05 gives max.error=x^5/5!=(0.05)^5/5!=2.6*e-9
  // NEW  approc_limit = 0.005 gives max.error=x^5/5!=2.6*e-14

  const G4double approc_limit = 0.005;
  G4ThreeVector  Bnorm, B_x_P, vperp, vpar;

  G4double B_d_P;
  G4double B_v_P;
  G4double Theta;
  G4double R_1;
  G4double R_Helix;
  G4double CosT2, SinT2, CosT, SinT;
  G4ThreeVector positionMove, endTangent;

  G4double Bmag = Bfld.mag();
  const G4double *pIn = yIn+3;
  G4ThreeVector initVelocity= G4ThreeVector( pIn[0], pIn[1], pIn[2]);
  G4double      velocityVal = initVelocity.mag();
  G4ThreeVector initTangent = (1.0/velocityVal) * initVelocity;
  
  R_1=GetInverseCurve(velocityVal,Bmag);

  // for too small magnetic fields there is no curvature
  // (include momentum here) FIXME

  if( (std::fabs(R_1) < 1e-10)||(Bmag<1e-12) )
  {
    LinearStep( yIn, h, yHelix );

    // Store and/or calculate parameters for chord distance

    SetAngCurve(1.);     
    SetCurve(h);
    SetRadHelix(0.);
  }
  else
  {
    Bnorm = (1.0/Bmag)*Bfld;

    // calculate the direction of the force
    
    B_x_P = Bnorm.cross(initTangent);

    // parallel and perp vectors

    B_d_P = Bnorm.dot(initTangent); // this is the fraction of P parallel to B

    vpar = B_d_P * Bnorm;       // the component parallel      to B
    vperp= initTangent - vpar;  // the component perpendicular to B
    
    B_v_P  = std::sqrt( 1 - B_d_P * B_d_P); // Fraction of P perp to B

    // calculate  the stepping angle
  
    Theta   = R_1 * h; // * B_v_P;

    // Trigonometrix
      
    if( std::fabs(Theta) > approc_limit )
    {
       SinT     = std::sin(Theta);
       CosT     = std::cos(Theta);
    }
    else
    {
      G4double Theta2 = Theta*Theta;
      G4double Theta3 = Theta2 * Theta;
      G4double Theta4 = Theta2 * Theta2;
      SinT     = Theta - 1.0/6.0 * Theta3;
      CosT     = 1 - 0.5 * Theta2 + 1.0/24.0 * Theta4;
    }

    // the actual "rotation"

    G4double R = 1.0 / R_1;

    positionMove  = R * ( SinT * vperp + (1-CosT) * B_x_P) + h * vpar;
    endTangent    = CosT * vperp + SinT * B_x_P + vpar;

    // Store the resulting position and tangent

    yHelix[0]   = yIn[0] + positionMove.x(); 
    yHelix[1]   = yIn[1] + positionMove.y(); 
    yHelix[2]   = yIn[2] + positionMove.z();
    yHelix[3] = velocityVal * endTangent.x();
    yHelix[4] = velocityVal * endTangent.y();
    yHelix[5] = velocityVal * endTangent.z();

    // Store 2*h step Helix if exist

    if(yHelix2)
    {
      SinT2     = 2.0 * SinT * CosT;
      CosT2     = 1.0 - 2.0 * SinT * SinT;
      endTangent    = (CosT2 * vperp + SinT2 * B_x_P + vpar);
      positionMove  = R * ( SinT2 * vperp + (1-CosT2) * B_x_P) + h*2 * vpar;
      
      yHelix2[0]   = yIn[0] + positionMove.x(); 
      yHelix2[1]   = yIn[1] + positionMove.y(); 
      yHelix2[2]   = yIn[2] + positionMove.z(); 
      yHelix2[3] = velocityVal * endTangent.x();
      yHelix2[4] = velocityVal * endTangent.y();
      yHelix2[5] = velocityVal * endTangent.z();
    }

    // Store and/or calculate parameters for chord distance

    G4double ptan=velocityVal*B_v_P;

    G4double particleCharge = fPtrMagEqOfMot->FCof() / (eplus*c_light); 
    R_Helix =std::abs( ptan/(fUnitConstant  * particleCharge*Bmag));
       
    SetAngCurve(std::abs(Theta));
    SetCurve(std::abs(R));
    SetRadHelix(R_Helix);
  }
}

Here is the call graph for this function:

Here is the caller graph for this function:

DistChord()

G4double G4MagHelicalStepper::DistChord ( ) const

virtual

Implements G4MagIntegratorStepper.

Definition at line 238 of file G4MagHelicalStepper.cc.

{
  // Check whether h/R >  pi  !!
  // Method DistLine is good only for <  pi

  G4double Ang=GetAngCurve();
  if(Ang<=pi)
  {
    return GetRadHelix()*(1-std::cos(0.5*Ang));
  }
  else
  {
    if(Ang<twopi)
    {
      return GetRadHelix()*(1+std::cos(0.5*(twopi-Ang)));
    }
    else  // return Diameter of projected circle
    {
      return 2*GetRadHelix();
    }
  }
}

Here is the call graph for this function:

DumbStepper()

virtual void G4MagHelicalStepper::DumbStepper ( const G4double  y[], G4ThreeVector  Bfld, G4double  h, G4double  yout[]  )

pure virtual

Implemented in G4HelixMixedStepper, G4ExactHelixStepper, G4HelixExplicitEuler, G4HelixImplicitEuler, G4HelixHeum, and G4HelixSimpleRunge.

Here is the caller graph for this function:

GetAngCurve()

G4double G4MagHelicalStepper::GetAngCurve ( ) const

inlineprotected

Here is the caller graph for this function:

GetCurve()

G4double G4MagHelicalStepper::GetCurve ( ) const

inlineprotected

GetInverseCurve()

G4double G4MagHelicalStepper::GetInverseCurve ( const G4double  Momentum, const G4double  Bmag  )

inlineprotected

Here is the caller graph for this function:

GetRadHelix()

G4double G4MagHelicalStepper::GetRadHelix ( ) const

inlineprotected

Here is the caller graph for this function:

LinearStep()

void G4MagHelicalStepper::LinearStep ( const G4double  yIn[], G4double  h, G4double  yHelix[]  ) const

inlineprotected

Here is the caller graph for this function:

MagFieldEvaluate()

void G4MagHelicalStepper::MagFieldEvaluate ( const G4double  y[], G4ThreeVector &  Bfield  )

inlineprotected

Here is the caller graph for this function:

operator=()

G4MagHelicalStepper& G4MagHelicalStepper::operator= ( const G4MagHelicalStepper &  )

private

SetAngCurve()

void G4MagHelicalStepper::SetAngCurve ( const G4double  Ang )

inlineprotected

Here is the caller graph for this function:

SetCurve()

void G4MagHelicalStepper::SetCurve ( const G4double  Curve )

inlineprotected

Here is the caller graph for this function:

SetRadHelix()

void G4MagHelicalStepper::SetRadHelix ( const G4double  Rad )

inlineprotected

Here is the caller graph for this function:

Stepper()

void G4MagHelicalStepper::Stepper ( const G4double  y[], const G4double  dydx[], G4double  h, G4double  yout[], G4double  yerr[]  )

virtual

Implements G4MagIntegratorStepper.

Reimplemented in G4HelixMixedStepper, and G4ExactHelixStepper.

Definition at line 192 of file G4MagHelicalStepper.cc.

{  
   const G4int nvar = 6;

   G4int i;

   // correction for Richardson Extrapolation.
   // G4double  correction = 1. / ( (1 << IntegratorOrder()) -1 );
   
   G4double      yTemp[7], yIn[7] ;
   G4ThreeVector Bfld_initial, Bfld_midpoint;
   
   //  Saving yInput because yInput and yOut can be aliases for same array

   for(i=0;i<nvar;i++) { yIn[i]=yInput[i]; }

   G4double h = hstep * 0.5; 

   MagFieldEvaluate(yIn, Bfld_initial) ;      

   // Do two half steps

   DumbStepper(yIn,   Bfld_initial,  h, yTemp);
   MagFieldEvaluate(yTemp, Bfld_midpoint) ;     
   DumbStepper(yTemp, Bfld_midpoint, h, yOut); 

   // Do a full Step

   h = hstep ;
   DumbStepper(yIn, Bfld_initial, h, yTemp);

   // Error estimation

   for(i=0;i<nvar;i++)
   {
     yErr[i] = yOut[i] - yTemp[i] ;
   }
   
   return;
}

Here is the call graph for this function:

Member Data Documentation

fAngCurve

G4double G4MagHelicalStepper::fAngCurve

private

Definition at line 133 of file G4MagHelicalStepper.hh.

fPtrMagEqOfMot

G4Mag_EqRhs* G4MagHelicalStepper::fPtrMagEqOfMot

private

Definition at line 130 of file G4MagHelicalStepper.hh.

frCurve

G4double G4MagHelicalStepper::frCurve

private

Definition at line 134 of file G4MagHelicalStepper.hh.

frHelix

G4double G4MagHelicalStepper::frHelix

private

Definition at line 135 of file G4MagHelicalStepper.hh.

fUnitConstant

const G4double G4MagHelicalStepper::fUnitConstant = 0.299792458*(GeV/(tesla*m))

staticprivate

Definition at line 127 of file G4MagHelicalStepper.hh.

yFinal

G4ThreeVector G4MagHelicalStepper::yFinal

private

Definition at line 137 of file G4MagHelicalStepper.hh.

yInitial

G4ThreeVector G4MagHelicalStepper::yInitial

private

Definition at line 137 of file G4MagHelicalStepper.hh.

yMidPoint

G4ThreeVector G4MagHelicalStepper::yMidPoint

private

Definition at line 137 of file G4MagHelicalStepper.hh.

---

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

• G4MagHelicalStepper.hh
• G4MagHelicalStepper.cc