Geant4 11.4.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4RK547FEq2 Class Reference

G4RK547FEq2 is an implementation of the 7 stage embedded Runge-Kutta 4,5 pair. More...

#include <G4RK547FEq2.hh>

Inheritance diagram for G4RK547FEq2:

Public Member Functions

 G4RK547FEq2 (G4EquationOfMotion *EqRhs, G4int integrationVariables=6)
 ~G4RK547FEq2 () override=default
 G4RK547FEq2 (const G4RK547FEq2 &)=delete
G4RK547FEq2operator= (const G4RK547FEq2 &)=delete
void Stepper (const G4double yInput[], const G4double dydx[], G4double hstep, G4double yOutput[], G4double yError[]) override
void Stepper (const G4double yInput[], const G4double dydx[], G4double hstep, G4double yOutput[], G4double yError[], G4double dydxOutput[])
G4double DistChord () const override
G4int IntegratorOrder () const override
G4StepperType StepperType () const override
Public Member Functions inherited from G4MagIntegratorStepper
 G4MagIntegratorStepper (G4EquationOfMotion *Equation, G4int numIntegrationVariables, G4int numStateVariables=12, G4bool isFSAL=false)
virtual ~G4MagIntegratorStepper ()=default
 G4MagIntegratorStepper (const G4MagIntegratorStepper &)=delete
G4MagIntegratorStepperoperator= (const G4MagIntegratorStepper &)=delete
void NormaliseTangentVector (G4double vec[6])
void NormalisePolarizationVector (G4double vec[12])
void RightHandSide (const G4double y[], G4double dydx[]) const
void RightHandSide (const G4double y[], G4double dydx[], G4double field[]) const
G4int GetNumberOfVariables () const
G4int GetNumberOfStateVariables () const
G4int IntegrationOrder ()
G4EquationOfMotionGetEquationOfMotion ()
const G4EquationOfMotionGetEquationOfMotion () const
void SetEquationOfMotion (G4EquationOfMotion *newEquation)
unsigned long GetfNoRHSCalls ()
void ResetfNORHSCalls ()
G4bool IsFSAL () const
G4bool isQSS () const
void SetIsQSS (G4bool val)

Additional Inherited Members

Protected Member Functions inherited from G4MagIntegratorStepper
void SetIntegrationOrder (G4int order)
void SetFSAL (G4bool flag=true)

Detailed Description

G4RK547FEq2 is an implementation of the 7 stage embedded Runge-Kutta 4,5 pair.

Definition at line 49 of file G4RK547FEq2.hh.

Constructor & Destructor Documentation

◆ G4RK547FEq2() [1/2]

G4RK547FEq2::G4RK547FEq2 ( G4EquationOfMotion * EqRhs,
G4int integrationVariables = 6 )

Constructor for G4RK547FEq2.

Parameters
[in]EqRhsPointer to the provided equation of motion.
[in]integrationVariablesThe number of integration variables.

Definition at line 51 of file G4RK547FEq2.cc.

52 : G4MagIntegratorStepper(EqRhs, integrationVariables)
53{
54}
G4MagIntegratorStepper(G4EquationOfMotion *Equation, G4int numIntegrationVariables, G4int numStateVariables=12, G4bool isFSAL=false)

Referenced by G4RK547FEq2(), and operator=().

◆ ~G4RK547FEq2()

G4RK547FEq2::~G4RK547FEq2 ( )
overridedefault

Default Destructor.

◆ G4RK547FEq2() [2/2]

G4RK547FEq2::G4RK547FEq2 ( const G4RK547FEq2 & )
delete

Copy constructor and assignment operator not allowed.

Member Function Documentation

◆ DistChord()

G4double G4RK547FEq2::DistChord ( ) const
overridevirtual

Returns the distance from chord line.

Implements G4MagIntegratorStepper.

Definition at line 180 of file G4RK547FEq2.cc.

181{
183 makeStep(fyIn, fdydx, fhstep / 2., yMid);
184
185 const G4ThreeVector begin = makeVector(fyIn, Value3D::Position);
186 const G4ThreeVector mid = makeVector(yMid, Value3D::Position);
187 const G4ThreeVector end = makeVector(fyOut, Value3D::Position);
188
189 return G4LineSection::Distline(mid, begin, end);
190}
CLHEP::Hep3Vector G4ThreeVector
double G4double
Definition G4Types.hh:83
static G4double Distline(const G4ThreeVector &OtherPnt, const G4ThreeVector &LinePntA, const G4ThreeVector &LinePntB)
G4ThreeVector makeVector(const ArrayType &array, Value3D value)

◆ IntegratorOrder()

G4int G4RK547FEq2::IntegratorOrder ( ) const
inlineoverridevirtual

Returns the order, 4, of integration.

Implements G4MagIntegratorStepper.

Definition at line 113 of file G4RK547FEq2.hh.

113{ return 4; }

◆ operator=()

G4RK547FEq2 & G4RK547FEq2::operator= ( const G4RK547FEq2 & )
delete

◆ Stepper() [1/2]

void G4RK547FEq2::Stepper ( const G4double yInput[],
const G4double dydx[],
G4double hstep,
G4double yOutput[],
G4double yError[] )
overridevirtual

The stepper for the Runge Kutta integration. The stepsize is fixed, with the step size given by 'hstep'. Integrates ODE starting values yInput[0 to 6]. Outputs yOutput[] and its estimated error yError[].

Parameters
[in]yInputStarting values array of integration variables.
[in]dydxDerivatives array.
[in]hstepThe given step size.
[out]yOutputIntegration output.
[out]yErrorThe estimated error.

Implements G4MagIntegratorStepper.

Definition at line 148 of file G4RK547FEq2.cc.

153{
154 copy(fyIn, yInput);
155 copy(fdydx, dydx);
156 fhstep = hstep;
157
158 makeStep(fyIn, fdydx, fhstep, fyOut, fdydxOut, yError);
159
160 copy(yOutput, fyOut);
161}
void copy(G4double dst[], const G4double src[], std::size_t size=G4FieldTrack::ncompSVEC)

◆ Stepper() [2/2]

void G4RK547FEq2::Stepper ( const G4double yInput[],
const G4double dydx[],
G4double hstep,
G4double yOutput[],
G4double yError[],
G4double dydxOutput[] )

Same as the Stepper() function above, with dydx also in ouput.

Parameters
[in]yInputStarting values array of integration variables.
[in]dydxDerivatives array.
[in]hstepThe given step size.
[out]yOutputIntegration output.
[out]yErrorThe estimated error.
[out]dydxOutputdydx in output.

Definition at line 163 of file G4RK547FEq2.cc.

169{
170 copy(fyIn, yInput);
171 copy(fdydx, dydx);
172 fhstep = hstep;
173
174 makeStep(fyIn, fdydx, fhstep, fyOut, fdydxOut, yError);
175
176 copy(yOutput, fyOut);
177 copy(dydxOutput, fdydxOut);
178}

◆ StepperType()

G4StepperType G4RK547FEq2::StepperType ( ) const
inlineoverridevirtual

Returns the stepper type-ID, "kRK547FEq2".

Reimplemented from G4MagIntegratorStepper.

Definition at line 118 of file G4RK547FEq2.hh.

118{ return kRK547FEq2; }
@ kRK547FEq2
G4RK547FEq2.

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