BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSIntegratorG4RK4MinStep.cc
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2023.
4
5This file is part of BDSIM.
6
7BDSIM is free software: you can redistribute it and/or modify
8it under the terms of the GNU General Public License as published
9by the Free Software Foundation version 3 of the License.
10
11BDSIM is distributed in the hope that it will be useful, but
12WITHOUT ANY WARRANTY; without even the implied warranty of
13MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14GNU General Public License for more details.
15
16You should have received a copy of the GNU General Public License
17along with BDSIM. If not, see <http://www.gnu.org/licenses/>.
18*/
19#include "BDSIntegratorG4RK4MinStep.hh"
20#include "BDSUtilities.hh"
21
22#include "G4ClassicalRK4.hh"
23#include "G4Mag_EqRhs.hh"
24#include "G4MagIntegratorStepper.hh"
25
27 G4double minimumStepSizeIn):
28 BDSIntegratorMag(eqOfMIn, 6),
29 g4integrator(new G4ClassicalRK4(eqOfMIn, 6)),
30 minimumStepSize(minimumStepSizeIn)
31{;}
32
33BDSIntegratorG4RK4MinStep::~BDSIntegratorG4RK4MinStep()
34{
35 delete g4integrator;
36}
37
38void BDSIntegratorG4RK4MinStep::Stepper(const G4double yIn[],
39 const G4double dydx[],
40 const G4double h,
41 G4double yOut[],
42 G4double yErr[])
43{
44 // in case of zero field or neutral particles do a linear step
45 const G4double fcof = eqOfM->FCof();
46 if (!BDS::IsFiniteStrength(fcof) || h < minimumStepSize)
47 {
48 AdvanceDriftMag(yIn,h,yOut,yErr);
49 SetDistChord(0);
50 return;
51 }
52 else
53 {
54 g4integrator->Stepper(yIn, dydx, h, yOut, yErr);
55 return;
56 }
57}
58
59
void AdvanceDriftMag(const G4double yIn[], const G4double h, G4double yOut[], G4double yErr[]) const
Error array [6] all 0.
virtual void Stepper(const G4double y[], const G4double dydx[], const G4double h, G4double yOut[], G4double yErr[])
BDSIntegratorG4RK4MinStep()
Private default constructor to enforce use of supplied constructor.
Common functionality to BDSIM integrators.
G4Mag_EqRhs * eqOfM
void SetDistChord(G4double distChordIn)
Setter for distChord to private member.
G4bool IsFiniteStrength(G4double variable)