Geant4  10.03
G4MagIntegratorDriver.hh
Go to the documentation of this file.
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: G4MagIntegratorDriver.hh 69699 2013-05-13 08:50:30Z gcosmo $
28 //
29 //
30 // class G4MagInt_Driver
31 //
32 // Class description:
33 //
34 // Provides a driver that talks to the Integrator Stepper, and insures that
35 // the error is within acceptable bounds.
36 
37 // History:
38 // - Created. J.Apostolakis.
39 // --------------------------------------------------------------------
40 
41 #ifndef G4MagInt_Driver_Def
42 #define G4MagInt_Driver_Def
43 
44 #include "G4Types.hh"
45 #include "G4FieldTrack.hh"
47 
49 {
50  public: // with description
51 
53  G4double hstep,
54  G4double eps, // Requested y_err/hstep
55  G4double hinitial=0.0); // Suggested 1st interval
56  // Above drivers for integrator (Runge-Kutta) with stepsize control.
57  // Integrates ODE starting values y_current
58  // from current s (s=s0) to s=s0+h with accuracy eps.
59  // On output ystart is replaced by value at end of interval.
60  // The concept is similar to the odeint routine from NRC p.721-722.
61 
62  G4bool QuickAdvance( G4FieldTrack& y_val, // INOUT
63  const G4double dydx[],
64  G4double hstep, // IN
65  G4double& dchord_step,
66  G4double& dyerr ) ;
67  // QuickAdvance just tries one Step - it does not ensure accuracy.
68 
69  G4bool QuickAdvance( G4FieldTrack& y_posvel, // INOUT
70  const G4double dydx[],
71  G4double hstep, // IN
72  G4double& dchord_step,
73  G4double& dyerr_pos_sq,
74  G4double& dyerr_mom_rel_sq ) ;
75  // New QuickAdvance that also just tries one Step
76  // (so also does not ensure accuracy)
77  // but does return the errors in position and
78  // momentum (normalised: Delta_Integration(p^2)/(p^2) )
79 
80  G4MagInt_Driver( G4double hminimum,
81  G4MagIntegratorStepper *pItsStepper,
82  G4int numberOfComponents=6,
83  G4int statisticsVerbosity=1);
85  // Constructor, destructor.
86 
87  inline G4double GetHmin() const;
88  inline G4double Hmin() const; // Obsolete
89  inline G4double GetSafety() const;
90  inline G4double GetPshrnk() const;
91  inline G4double GetPgrow() const;
92  inline G4double GetErrcon() const;
93  inline void GetDerivatives( const G4FieldTrack &y_curr, // const, INput
94  G4double dydx[] ); // OUTput
95  // Accessors.
96 
97  inline void RenewStepperAndAdjust(G4MagIntegratorStepper *pItsStepper);
98  // Sets a new stepper pItsStepper for this driver. Then it calls
99  // ReSetParameters to reset its parameters accordingly.
100 
101  inline void ReSetParameters(G4double new_safety= 0.9 );
102  // i) sets the exponents (pgrow & pshrnk),
103  // using the current Stepper's order,
104  // ii) sets the safety
105  // ii) calculates "errcon" according to the above values.
106 
107  inline void SetSafety(G4double valS);
108  inline void SetPshrnk(G4double valPs);
109  inline void SetPgrow (G4double valPg);
110  inline void SetErrcon(G4double valEc);
111  // When setting safety or pgrow, errcon will be set to a
112  // compatible value.
113 
114  inline G4double ComputeAndSetErrcon();
115 
116  inline const G4MagIntegratorStepper* GetStepper() const;
118 
119  void OneGoodStep( G4double ystart[], // Like old RKF45step()
120  const G4double dydx[],
121  G4double& x,
122  G4double htry,
123  G4double eps, // memb variables ?
124  G4double& hdid,
125  G4double& hnext ) ;
126  // This takes one Step that is as large as possible while
127  // satisfying the accuracy criterion of:
128  // yerr < eps * |y_end-y_start|
129 
130  G4double ComputeNewStepSize( G4double errMaxNorm, // normalised error
131  G4double hstepCurrent); // current step size
132  // Taking the last step's normalised error, calculate
133  // a step size for the next step.
134  // Do not limit the next step's size within a factor of the
135  // current one.
136 
138  G4double errMaxNorm, // normalised error
139  G4double hstepCurrent); // current step size
140  // Taking the last step's normalised error, calculate
141  // a step size for the next step.
142  // Limit the next step's size within a range around the current one.
143 
144  inline G4int GetMaxNoSteps() const;
145  inline void SetMaxNoSteps( G4int val);
146  // Modify and Get the Maximum number of Steps that can be
147  // taken for the integration of a single segment -
148  // (ie a single call to AccurateAdvance).
149 
150  public: // without description
151 
152  inline void SetHmin(G4double newval);
153  inline void SetVerboseLevel(G4int newLevel);
154  inline G4double GetVerboseLevel() const;
155 
156  inline G4double GetSmallestFraction() const;
157  void SetSmallestFraction( G4double val );
158 
159  protected: // without description
160  void WarnSmallStepSize( G4double hnext, G4double hstep,
161  G4double h, G4double xDone,
162  G4int noSteps);
163  void WarnTooManyStep( G4double x1start, G4double x2end, G4double xCurrent);
164  void WarnEndPointTooFar (G4double endPointDist,
165  G4double hStepSize ,
166  G4double epsilonRelative,
167  G4int debugFlag);
168  // Issue warnings for undesirable situations
169 
170  void PrintStatus( const G4double* StartArr,
171  G4double xstart,
172  const G4double* CurrentArr,
173  G4double xcurrent,
174  G4double requestStep,
175  G4int subStepNo );
176  void PrintStatus( const G4FieldTrack& StartFT,
177  const G4FieldTrack& CurrentFT,
178  G4double requestStep,
179  G4int subStepNo );
180  void PrintStat_Aux( const G4FieldTrack& aFieldTrack,
181  G4double requestStep,
182  G4double actualStep,
183  G4int subStepNo,
184  G4double subStepSize,
185  G4double dotVelocities );
186  // Verbose output for debugging
187 
188  void PrintStatisticsReport() ;
189  // Report on the number of steps, maximum errors etc.
190 
191 #ifdef QUICK_ADV_TWO
192  G4bool QuickAdvance( G4double yarrin[], // In
193  const G4double dydx[],
194  G4double hstep,
195  G4double yarrout[], // Out
196  G4double& dchord_step, // Out
197  G4double& dyerr ); // in length
198 #endif
199 
200  private:
201 
204  // Private copy constructor and assignment operator.
205 
206  private:
207 
208  // ---------------------------------------------------------------
209  // INVARIANTS
210 
212  // Minimum Step allowed in a Step (in absolute units)
213  G4double fSmallestFraction; // Expected range 1e-12 to 5e-15;
214  // Smallest fraction of (existing) curve length - in relative units
215  // below this fraction the current step will be the last
216 
217  const G4int fNoIntegrationVariables; // Number of Variables in integration
218  const G4int fMinNoVars; // Minimum number for FieldTrack
219  const G4int fNoVars; // Full number of variable
220 
222  static const G4int fMaxStepBase;
223 
225  G4double pshrnk; // exponent for shrinking
226  G4double pgrow; // exponent for growth
228  // Parameters used to grow and shrink trial stepsize.
229 
232  // Maximum stepsize increase/decrease factors.
233 
235 
236  // ---------------------------------------------------------------
237  // DEPENDENT Objects
239 
240  // ---------------------------------------------------------------
241  // STATE
242 
247  // Step Statistics
248 
249  G4int fVerboseLevel; // Verbosity level for printing (debug, ..)
250  // Could be varied during tracking - to help identify issues
251 
252 };
253 
254 #include "G4MagIntegratorDriver.icc"
255 
256 #endif /* G4MagInt_Driver_Def */
void SetPgrow(G4double valPg)
G4double GetSmallestFraction() const
G4MagIntegratorStepper * pIntStepper
void SetVerboseLevel(G4int newLevel)
const G4MagIntegratorStepper * GetStepper() const
void RenewStepperAndAdjust(G4MagIntegratorStepper *pItsStepper)
void WarnTooManyStep(G4double x1start, G4double x2end, G4double xCurrent)
void WarnEndPointTooFar(G4double endPointDist, G4double hStepSize, G4double epsilonRelative, G4int debugFlag)
void GetDerivatives(const G4FieldTrack &y_curr, G4double dydx[])
static const G4double eps
G4MagInt_Driver(G4double hminimum, G4MagIntegratorStepper *pItsStepper, G4int numberOfComponents=6, G4int statisticsVerbosity=1)
G4double Hmin() const
int G4int
Definition: G4Types.hh:78
void OneGoodStep(G4double ystart[], const G4double dydx[], G4double &x, G4double htry, G4double eps, G4double &hdid, G4double &hnext)
G4MagInt_Driver & operator=(const G4MagInt_Driver &)
G4double GetSafety() const
static const G4double max_stepping_increase
void SetHmin(G4double newval)
void SetMaxNoSteps(G4int val)
void SetPshrnk(G4double valPs)
G4int GetMaxNoSteps() const
void PrintStatus(const G4double *StartArr, G4double xstart, const G4double *CurrentArr, G4double xcurrent, G4double requestStep, G4int subStepNo)
G4double ComputeNewStepSize_WithinLimits(G4double errMaxNorm, G4double hstepCurrent)
bool G4bool
Definition: G4Types.hh:79
G4double GetHmin() const
G4double GetPgrow() const
void SetSafety(G4double valS)
G4double GetVerboseLevel() const
void PrintStat_Aux(const G4FieldTrack &aFieldTrack, G4double requestStep, G4double actualStep, G4int subStepNo, G4double subStepSize, G4double dotVelocities)
G4double GetPshrnk() const
G4double GetErrcon() const
G4bool QuickAdvance(G4FieldTrack &y_val, const G4double dydx[], G4double hstep, G4double &dchord_step, G4double &dyerr)
void SetErrcon(G4double valEc)
void WarnSmallStepSize(G4double hnext, G4double hstep, G4double h, G4double xDone, G4int noSteps)
void SetSmallestFraction(G4double val)
void ReSetParameters(G4double new_safety=0.9)
G4double ComputeNewStepSize(G4double errMaxNorm, G4double hstepCurrent)
double G4double
Definition: G4Types.hh:76
static const G4int fMaxStepBase
G4bool AccurateAdvance(G4FieldTrack &y_current, G4double hstep, G4double eps, G4double hinitial=0.0)
const G4int fNoIntegrationVariables
G4double ComputeAndSetErrcon()
static const G4double max_stepping_decrease