BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
|
Exact helix through pure dipole field. More...
#include <BDSIntegratorDipoleRodrigues2.hh>
Public Member Functions | |
BDSIntegratorDipoleRodrigues2 (G4Mag_EqRhs *eqOfMIn, G4double minimumRadiusOfCurvature) | |
virtual void | DumbStepper (const G4double yIn[6], G4ThreeVector field, G4double stepLength, G4double yOut[6]) |
void | SingleStep (const G4double yIn[6], const G4double &h, G4double yOut[6]) |
virtual void | Stepper (const G4double yIn[6], const G4double dydx[], G4double h, G4double yOut[6], G4double yErr[6]) |
Output error array. | |
virtual G4int | IntegratorOrder () const |
void | AdvanceHelixForSpiralling (const G4double yIn[6], const G4ThreeVector &field, const G4double &h, G4double yOut[6], G4double yErr[6]) |
G4double | RadiusOfHelix () const |
Public accessor for protected variable in base class. | |
![]() | |
void | AdvanceDriftMag (const G4double yIn[], const G4double h, G4double yOut[], G4double yErr[]) const |
Error array [6] all 0. | |
void | AdvanceDriftMag (const G4double yIn[], const G4ThreeVector &unitMomentum, const G4double h, G4double yOut[], G4double yErr[]) const |
Protected Member Functions | |
void | FudgeDistChordToZero () |
Protected Attributes | |
G4Mag_EqRhs * | eqOfM |
Cache of equation of motion. This class does not own it. | |
Private Member Functions | |
BDSIntegratorDipoleRodrigues2 ()=delete | |
Private default constructor to force use of provided one. | |
Private Attributes | |
G4double | minimumRadiusOfCurvature |
Exact helix through pure dipole field.
Integrator based on G4MagHelicalStepper that provides required tracking througha pure magnetic field, but is a pure virtual class.
This class implements the required Stepper method that uses the provided AdvanceHelix method from the base class.
Unlike other BDSIM integrators this works in global coordinates and does not require the transforms to curvilinear coordinates.
For low momentum particles (|p| < 40 MeV) we use the G4ClassicalRK4 integrator instead of the one for the pure magnetic field as it is more robust.
Note, some things that could be passed by reference are not purely for correct overload of Geant4 function.
Definition at line 53 of file BDSIntegratorDipoleRodrigues2.hh.
BDSIntegratorDipoleRodrigues2::BDSIntegratorDipoleRodrigues2 | ( | G4Mag_EqRhs * | eqOfMIn, |
G4double | minimumRadiusOfCurvature | ||
) |
Definition at line 32 of file BDSIntegratorDipoleRodrigues2.cc.
|
inlinevirtual |
Definition at line 59 of file BDSIntegratorDipoleRodrigues2.hh.
void BDSIntegratorDipoleRodrigues2::AdvanceHelixForSpiralling | ( | const G4double | yIn[6], |
const G4ThreeVector & | field, | ||
const G4double & | h, | ||
G4double | yOut[6], | ||
G4double | yErr[6] | ||
) |
Variation of AdvanceHelix specifically to deal with particles that are likely to be spiralling in the magnetic field.
Definition at line 156 of file BDSIntegratorDipoleRodrigues2.cc.
References BDSIntegratorDrift::AdvanceDriftMag(), FudgeDistChordToZero(), and BDS::IsFinite().
Referenced by Stepper().
|
virtual |
Required to be provided by base class, but apparently should never be called by the driver. Simply calls AdvanceHelix.
Definition at line 39 of file BDSIntegratorDipoleRodrigues2.cc.
|
protected |
DistChord() is non-virtual function in base class so set Ang and RadHelix appropriately such that DistChord() will return 0.
Definition at line 196 of file BDSIntegratorDipoleRodrigues2.cc.
Referenced by AdvanceHelixForSpiralling(), BDSIntegratorDipoleFringe::BaseStepper(), and Stepper().
|
inlinevirtual |
Definition at line 84 of file BDSIntegratorDipoleRodrigues2.hh.
|
inline |
Public accessor for protected variable in base class.
Definition at line 95 of file BDSIntegratorDipoleRodrigues2.hh.
Referenced by BDSIntegratorDipoleQuadrupole::Stepper().
void BDSIntegratorDipoleRodrigues2::SingleStep | ( | const G4double | yIn[6], |
const G4double & | h, | ||
G4double | yOut[6] | ||
) |
Perform a single step. Useful for other integrators that want to use this one to estimate the radius of curvature in the field, in which case, no need for usual two step calculation, nor error calculation.
Definition at line 47 of file BDSIntegratorDipoleRodrigues2.cc.
References eqOfM.
Referenced by BDSIntegratorDipoleQuadrupole::Stepper().
|
virtual |
Output error array.
Calculate output coordinates. Decide if particle is spiralling or not. Nominally calculate two half steps and compare to one full step for error estimation.
yIn | Input coordinates x,y,z,px,py,pz,t |
dydx | Partial differentials at yInput |
h | Length of trajectory to calculate |
yOut | Output array |
Definition at line 59 of file BDSIntegratorDipoleRodrigues2.cc.
References BDSIntegratorDrift::AdvanceDriftMag(), AdvanceHelixForSpiralling(), eqOfM, FudgeDistChordToZero(), BDS::IsFiniteStrength(), and minimumRadiusOfCurvature.
Referenced by BDSIntegratorDipoleFringe::BaseStepper(), and BDSIntegratorDipoleQuadrupole::Stepper().
|
protected |
Cache of equation of motion. This class does not own it.
Definition at line 103 of file BDSIntegratorDipoleRodrigues2.hh.
Referenced by SingleStep(), BDSIntegratorDipoleFringe::Stepper(), Stepper(), and BDSIntegratorDipoleFringeScaling::Stepper().
|
private |
The minimum tolerable radius of curvature before we decide the particle is spiralling and should be treated differently.
Definition at line 111 of file BDSIntegratorDipoleRodrigues2.hh.
Referenced by Stepper().