Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4IntegrationDriver.hh
이 파일의 문서화 페이지로 가기
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4IntegrationDriver.hh 106739 2017-10-20 14:45:37Z dmsoroki $
28 //
29 //
30 // class G4IntegrationDriver
31 //
32 // Class description:
33 //
34 // Driver class which controls the integration error of a
35 // Runge-Kutta stepper
36 
37 // History:
38 // - Created. D.Sorokin
39 // --------------------------------------------------------------------
40 
41 #ifndef G4IntegrationDriver_HH
42 #define G4IntegrationDriver_HH
43 
44 #include "G4Types.hh"
45 #include "G4FieldTrack.hh"
46 #include "G4VIntegrationDriver.hh"
47 
48 template <class T>
50 public:
51  G4IntegrationDriver( G4double hminimum,
52  T* stepper,
53  G4int numberOfComponents = 6,
54  G4int statisticsVerbosity = 1);
55 
56  virtual ~G4IntegrationDriver() override;
57 
58  G4IntegrationDriver(const G4IntegrationDriver &) = delete;
59  const G4IntegrationDriver& operator =(const G4IntegrationDriver &) = delete;
60 
61  // Integrates ODE from current s (s=s0) to s=s0+h with accuracy eps.
62  // On output track is replaced by value at end of interval.
63  // The concept is similar to the odeint routine from NRC p.721-722.
64  virtual G4bool AccurateAdvance(
66  G4double hstep,
67  G4double eps, // Requested y_err/hstep
68  G4double hinitial = 0) override; // Suggested 1st interval
69 
70  // QuickAdvance just tries one Step - it does not ensure accuracy.
71  virtual G4bool QuickAdvance(
72  G4FieldTrack& fieldTrack,
73  const G4double dydx[],
74  G4double hstep,
75  G4double& dchord_step,
76  G4double& dyerr) override;
77 
78  virtual void GetDerivatives(
79  const G4FieldTrack &track,
80  G4double dydx[]) const override;
81 
82  // Taking the last step's normalised error, calculate
83  // a step size for the next step.
84  // Do not limit the next step's size within a factor of the
85  // current one.
87  G4double errMaxNorm, // normalised error
88  G4double hstepCurrent) override; // current step size
89 
90  virtual void SetVerboseLevel(G4int newLevel) override;
91  virtual G4int GetVerboseLevel() const override;
92 
93  virtual G4EquationOfMotion* GetEquationOfMotion() override;
94  virtual void SetEquationOfMotion(G4EquationOfMotion* equation) override;
95 
96  virtual const T* GetStepper() const override;
97  virtual T* GetStepper() override;
98 
99  // Accessors.
100  G4double GetMinimumStep() const;
101  G4double GetSafety() const;
102  G4double GetPshrnk() const;
103  G4double GetPgrow() const;
104 
105  virtual void RenewStepperAndAdjust(G4MagIntegratorStepper *pItsStepper) override;
106 
107  // Sets a new stepper pItsStepper for this driver. Then it calls
108  // ReSetParameters to reset its parameters accordingly.
109  inline void RenewStepperAndAdjustStrict(T *pItsStepper);
110 
111  // i) sets the exponents (pgrow & pshrnk),
112  // using the current Stepper's order,
113  // ii) sets the safety
114  void ReSetParameters(G4double safety = 0.9);
115 
116  void SetMinimumStep(G4double newval);
117  void SetSafety(G4double valS);
118 
119  // This takes one Step that is of size htry, or as large
120  // as possible while satisfying the accuracy criterion of:
121  // yerr < eps * |y_end-y_start|
122  void OneGoodStep( G4double yVar[], // InOut
123  const G4double dydx[],
124  G4double& curveLength,
125  G4double htry,
126  G4double eps,
127  G4double& hdid,
128  G4double& hnext);
129 
130  // Modify and Get the Maximum number of Steps that can be
131  // taken for the integration of a single segment -
132  // (ie a single call to AccurateAdvance).
133  G4int GetMaxNoSteps() const;
134  void SetMaxNoSteps( G4int val);
135 
137  void SetSmallestFraction(G4double val);
138 
139 private:
141  G4double GrowStepSize(G4double h, G4double error) const;
142  void UpdateErrorConstraints();
143 
144  void CheckStep(
145  const G4ThreeVector& posIn, const G4ThreeVector& posOut, G4double hdid);
146 
147  // Minimum Step allowed in a Step (in absolute units)
149 
150  // Smallest fraction of (existing) curve length - in relative units
151  // below this fraction the current step will be the last
153  // Expected range: smaller than 0.1 * epsilon and bigger than 5e-13
154  // ( Note: this range is not enforced. )
156 
157  // The (default) maximum number of steps is Base
158  // divided by the order of Stepper
160 
161  // Parameters used to grow and shrink trial stepsize.
163  G4double pshrnk; // exponent for shrinking
164  G4double pgrow; // exponent for growth
165 
166  // muximum error values for shrinking / growing (optimisation).
169 
171 
172  // Step Statistics
174 
175  G4int fVerboseLevel; // Verbosity level for printing (debug, ..)
176  // Could be varied during tracking - to help identify issues
177 
179 };
180 
181 #include "G4IntegrationDriver.icc"
182 
183 #endif
void SetSafety(G4double valS)
void RenewStepperAndAdjustStrict(T *pItsStepper)
void SetMinimumStep(G4double newval)
G4double GetPgrow() const
void OneGoodStep(G4double yVar[], const G4double dydx[], G4double &curveLength, G4double htry, G4double eps, G4double &hdid, G4double &hnext)
G4int GetMaxNoSteps() const
virtual void SetVerboseLevel(G4int newLevel) override
G4double GetSmallestFraction() const
G4double GetMinimumStep() const
virtual G4EquationOfMotion * GetEquationOfMotion() override
void SetSmallestFraction(G4double val)
virtual G4bool AccurateAdvance(G4FieldTrack &track, G4double hstep, G4double eps, G4double hinitial=0) override
virtual G4bool QuickAdvance(G4FieldTrack &fieldTrack, const G4double dydx[], G4double hstep, G4double &dchord_step, G4double &dyerr) override
virtual G4int GetVerboseLevel() const override
void SetMaxNoSteps(G4int val)
void UpdateErrorConstraints()
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
void CheckStep(const G4ThreeVector &posIn, const G4ThreeVector &posOut, G4double hdid)
G4IntegrationDriver(G4double hminimum, T *stepper, G4int numberOfComponents=6, G4int statisticsVerbosity=1)
virtual ~G4IntegrationDriver() override
void ReSetParameters(G4double safety=0.9)
G4double ShrinkStepSize(G4double h, G4double error) const
int G4int
Definition: G4Types.hh:78
virtual const T * GetStepper() const override
virtual G4double ComputeNewStepSize(G4double errMaxNorm, G4double hstepCurrent) override
G4double GetPshrnk() const
G4double GrowStepSize(G4double h, G4double error) const
virtual void RenewStepperAndAdjust(G4MagIntegratorStepper *pItsStepper) override
static PROLOG_HANDLER error
Definition: xmlrole.cc:112
static const G4double eps
virtual void GetDerivatives(const G4FieldTrack &track, G4double dydx[]) const override
const G4IntegrationDriver & operator=(const G4IntegrationDriver &)=delete
virtual void SetEquationOfMotion(G4EquationOfMotion *equation) override
G4double GetSafety() const