| Geant4 Cross Reference |
1 // ********************************************************************
2 // * License and Disclaimer *
3 // * *
4 // * The Geant4 software is copyright of the Copyright Holders of *
5 // * the Geant4 Collaboration. It is provided under the terms and *
6 // * conditions of the Geant4 Software License, included in the file *
7 // * LICENSE and available at http://cern.ch/geant4/license . These *
8 // * include a list of copyright holders. *
9 // * *
10 // * Neither the authors of this software system, nor their employing *
11 // * institutes,nor the agencies providing financial support for this *
12 // * work make any representation or warranty, express or implied, *
13 // * regarding this software system or assume any liability for its *
14 // * use. Please see the license in the file LICENSE and URL above *
15 // * for the full disclaimer and the limitation of liability. *
16 // * *
17 // * This code implementation is the result of the scientific and *
18 // * technical work of the GEANT4 collaboration. *
19 // * By using, copying, modifying or distributing the software (or *
20 // * any work based on the software) you agree to acknowledge its *
21 // * use in resulting scientific publications, and indicate your *
22 // * acceptance of all terms of the Geant4 Software license. *
23 // ********************************************************************
24 //
25 // G4BorisDriver
26 //
27 // Class description:
28 //
29 // G4BorisDriver is a driver class using the second order Boris
30 // method to integrate the equation of motion.
31 //
32 //
33 // Author: Divyansh Tiwari, Google Summer of Code 2022
34 // Supervision: John Apostolakis,Renee Fatemi, Soon Yung Jun
35 // --------------------------------------------------------------------
36 #ifndef G4BORIS_DRIVER_HH
37 #define G4BORIS_DRIVER_HH
38
39 #include "G4VIntegrationDriver.hh"
40 #include "G4BorisScheme.hh"
41 #include "G4ChordFinderDelegate.hh"
42
43
44 class G4BorisDriver : public G4VIntegrationDriver,
45 public G4ChordFinderDelegate<G4BorisDriver>
46 {
47 public:
48
49 G4BorisDriver( G4double hminimum,
50 G4BorisScheme* Boris,
51 G4int numberOfComponents = 6,
52 G4bool verbosity = false);
53
54 inline ~G4BorisDriver() override = default;
55
56 inline G4BorisDriver(const G4BorisDriver&) = delete;
57 inline G4BorisDriver& operator=(const G4BorisDriver&) = delete;
58
59 // 1. Core methods that advance the integration
60
61 G4bool AccurateAdvance( G4FieldTrack& track,
62 G4double stepLen,
63 G4double epsilon,
64 G4double beginStep = 0) override;
65 // Advance integration accurately
66 // - by relative accuracy better than 'epsilon'
67
68 G4bool QuickAdvance( G4FieldTrack& y_val, // In/Out
69 const G4double dydx[],
70 G4double hstep,
71 G4double& missDist, // Out: estimated sagitta
72 G4double& dyerr ) override;
73 // Attempt one integration step, and return estimated error 'dyerr'
74
75 void OneGoodStep(G4double yCurrentState[], // In/Out: state ('y')
76 G4double& curveLength, // In/Out: 'x'
77 G4double htry, // step to attempt
78 G4double epsilon_rel, // relative accuracy
79 G4double restMass,
80 G4double charge,
81 G4double& hdid, // Out: step achieved
82 G4double& hnext); // Out: proposed next step
83 // Method to implement Accurate Advance
84
85 // 2. Methods needed to co-work with G4ChordFinder
86
87 G4double AdvanceChordLimited(G4FieldTrack& track,
88 G4double hstep,
89 G4double eps,
90 G4double chordDistance) override
91 {
92 return ChordFinderDelegate::
93 AdvanceChordLimitedImpl(track, hstep, eps, chordDistance);
94 }
95
96 void OnStartTracking() override
97 {
98 ChordFinderDelegate::ResetStepEstimate();
99 }
100
101 void OnComputeStep(const G4FieldTrack*) override {}
102
103 // 3. Does the method redo integrations when called to obtain values for
104 // internal, smaller intervals? (when needed to identify an intersection)
105
106 G4bool DoesReIntegrate() const override { return true; }
107 // It would be no if it just used interpolation to provide a result.
108
109 // 4. Relevant for calculating a new step size to achieve required accuracy
110
111 inline G4double ComputeNewStepSize(G4double errMaxNorm, // normalised error
112 G4double hstepCurrent) override; // current step size
113
114 G4double ShrinkStepSize2(G4double h, G4double error2) const;
115 G4double GrowStepSize2(G4double h, G4double error2) const;
116 // Calculate the next step size given the square of the relative error
117
118 // 5. Auxiliary Methods ...
119
120 void GetDerivatives( const G4FieldTrack& track,
121 G4double dydx[] ) const override;
122
123 void GetDerivatives( const G4FieldTrack& track,
124 G4double dydx[],
125 G4double field[] ) const override;
126
127 inline void SetVerboseLevel(G4int level) override;
128 inline G4int GetVerboseLevel() const override;
129
130 inline G4EquationOfMotion* GetEquationOfMotion() override;
131 inline const G4EquationOfMotion* GetEquationOfMotion() const;
132 void SetEquationOfMotion(G4EquationOfMotion* equation) override;
133
134 void StreamInfo( std::ostream& os ) const override;
135 // Write out the parameters / state of the driver
136
137 // 6. Not relevant for Boris and other non-RK methods
138
139 inline const G4MagIntegratorStepper* GetStepper() const override;
140 inline G4MagIntegratorStepper* GetStepper() override;
141
142 private:
143
144 inline G4int GetNumberOfVariables() const;
145
146 inline void CheckStep(const G4ThreeVector& posIn,
147 const G4ThreeVector& posOut,
148 G4double hdid) const;
149
150 private:
151
152 // INVARIANTS -- remain unchanged during tracking / integration
153 // Parameters
154 G4double fMinimumStep;
155 G4bool fVerbosity;
156
157 // State -- The core stepping algorithm
158 G4BorisScheme* boris;
159
160 // STATE -- intermediate state (to avoid creation / churn )
161 G4double yIn[G4FieldTrack::ncompSVEC],
162 yMid[G4FieldTrack::ncompSVEC],
163 yOut[G4FieldTrack::ncompSVEC],
164 yError[G4FieldTrack::ncompSVEC];
165
166 G4double yCurrent[G4FieldTrack::ncompSVEC];
167
168 // - Unused 2022.11.03:
169 // G4double derivs[2][6][G4FieldTrack::ncompSVEC];
170 // const G4int interval_sequence[2];
171
172 // INVARIANTS -- Parameters for ensuring that one call has finite number of integration steps
173 static constexpr G4int fMaxNoSteps = 300;
174 static constexpr G4double fSmallestFraction= 1e-12; // To avoid FP underflow ! ( 1.e-6 for single prec)
175
176 static constexpr G4int fIntegratorOrder= 2; // 2nd order method -- needed for error control
177 static constexpr G4double fSafetyFactor = 0.9; //
178
179 static constexpr G4double fMaxSteppingIncrease= 10.0; // Increase no more than 10x
180 static constexpr G4double fMaxSteppingDecrease= 0.1; // Reduce no more than 10x
181 static constexpr G4double fPowerShrink = -1.0 / fIntegratorOrder;
182 static constexpr G4double fPowerGrow = -1.0 / (1.0 + fIntegratorOrder);
183
184 static const G4double fErrorConstraintShrink;
185 static const G4double fErrorConstraintGrow;
186
187 using ChordFinderDelegate = G4ChordFinderDelegate<G4BorisDriver>;
188 };
189
190 #include "G4BorisDriver.icc"
191
192 #endif
193