|
|
|||
File indexing completed on 2026-10-04 09:07:21
0001 // 0002 // ******************************************************************** 0003 // * License and Disclaimer * 0004 // * * 0005 // * The Geant4 software is copyright of the Copyright Holders of * 0006 // * the Geant4 Collaboration. It is provided under the terms and * 0007 // * conditions of the Geant4 Software License, included in the file * 0008 // * LICENSE and available at http://cern.ch/geant4/license . These * 0009 // * include a list of copyright holders. * 0010 // * * 0011 // * Neither the authors of this software system, nor their employing * 0012 // * institutes,nor the agencies providing financial support for this * 0013 // * work make any representation or warranty, express or implied, * 0014 // * regarding this software system or assume any liability for its * 0015 // * use. Please see the license in the file LICENSE and URL above * 0016 // * for the full disclaimer and the limitation of liability. * 0017 // * * 0018 // * This code implementation is the result of the scientific and * 0019 // * technical work of the GEANT4 collaboration. * 0020 // * By using, copying, modifying or distributing the software (or * 0021 // * any work based on the software) you agree to acknowledge its * 0022 // * use in resulting scientific publications, and indicate your * 0023 // * acceptance of all terms of the Geant4 Software license. * 0024 // ******************************************************************** 0025 // 0026 // G4MagInt_Driver 0027 // 0028 // Class description: 0029 // 0030 // Provides a driver that talks to the Integrator Stepper, and insures that 0031 // the error is within acceptable bounds. 0032 0033 // Author: Vladimir Grichine (CERN), 07.10.1996 - Created 0034 // W.Wander (MIT), 28.01.1998 - Added ability for low order integrators 0035 // -------------------------------------------------------------------- 0036 #ifndef G4MAGINT_DRIVER_HH 0037 #define G4MAGINT_DRIVER_HH 0038 0039 #include "G4VIntegrationDriver.hh" 0040 #include "G4MagIntegratorStepper.hh" 0041 #include "G4ChordFinderDelegate.hh" 0042 0043 /** 0044 * @brief G4MagInt_Driver provides a driver that talks to the Integrator 0045 * Stepper and insures that the error is within acceptable bounds. 0046 */ 0047 0048 class G4MagInt_Driver : public G4VIntegrationDriver, 0049 public G4ChordFinderDelegate<G4MagInt_Driver> 0050 { 0051 public: 0052 0053 /** 0054 * Constructor for G4MagInt_Driver. 0055 * @param[in] hminimum The minumum allowed step. 0056 * @param[in] pItsStepper Pointer to the integrator stepper. 0057 * @param[in] numberOfComponents The number of integration variables. 0058 * @param[in] statisticsVerbosity Flag for verbosity. 0059 */ 0060 G4MagInt_Driver(G4double hminimum, 0061 G4MagIntegratorStepper* pItsStepper, 0062 G4int numberOfComponents = 6, 0063 G4int statisticsVerbosity = 0); 0064 0065 /** 0066 * Destructor. Provides statistics if verbosity level is greater than 1. 0067 */ 0068 ~G4MagInt_Driver() override; 0069 0070 /** 0071 * Copy constructor and assignment operator not allowed. 0072 */ 0073 G4MagInt_Driver(const G4MagInt_Driver&) = delete; 0074 G4MagInt_Driver& operator=(const G4MagInt_Driver&) = delete; 0075 0076 /** 0077 * Computes the step to take, based on chord limits. 0078 * @param[in,out] track The current track in field. 0079 * @param[in] stepMax Proposed maximum step length. 0080 * @param[in] epsStep Requested accuracy, y_err/hstep. 0081 * @param[in] chordDistance Maximum sagitta distance. 0082 * @returns The length of step taken. 0083 */ 0084 inline G4double AdvanceChordLimited(G4FieldTrack& track, 0085 G4double stepMax, 0086 G4double epsStep, 0087 G4double chordDistance) override; 0088 0089 /** 0090 * Dispatch interface method for initialisation/reset of driver. 0091 */ 0092 inline void OnStartTracking() override; 0093 0094 /** 0095 * Dispatch interface method for computing step. Does nothing here. 0096 */ 0097 inline void OnComputeStep(const G4FieldTrack* = nullptr) override {} 0098 0099 /** 0100 * The driver implements re-integration, so returns true. 0101 */ 0102 G4bool DoesReIntegrate() const override { return true; } 0103 0104 /** 0105 * Advances integration accurately by relative accuracy better than 'eps'. 0106 * @param[in,out] y_current The current track in field. 0107 * @param[in] hstep Proposed step length. 0108 * @param[in] eps Requested accuracy, y_err/hstep. 0109 * @param[in] hinitial Initial minimum integration step. 0110 * @returns true if integration succeeds. 0111 */ 0112 G4bool AccurateAdvance(G4FieldTrack& y_current, 0113 G4double hstep, 0114 G4double eps, // Requested y_err/hstep 0115 G4double hinitial = 0.0) override; 0116 0117 /** 0118 * Attempts one integration step, and returns estimated error 'dyerr'. 0119 * It does not ensure accuracy. 0120 * @param[in,out] y_val The current track in field. 0121 * @param[in] dydx dydx array. 0122 * @param[in] hstep Proposed step length. 0123 * @param[out] dchord_step Estimated sagitta distance. 0124 * @param[out] dyerr Estimated error. 0125 * @returns true if integration succeeds. 0126 */ 0127 G4bool QuickAdvance(G4FieldTrack& y_val, // In/Out 0128 const G4double dydx[], 0129 G4double hstep, 0130 G4double& dchord_step, 0131 G4double& dyerr) override; 0132 0133 /** 0134 * Writes out to stream the parameters/state of the driver. 0135 */ 0136 void StreamInfo( std::ostream& os ) const override; 0137 0138 /** 0139 * Attempts one integration step, and returns estimated error 'dyerr'. 0140 * It does not ensure accuracy. 0141 * @param[in,out] y_posvel The current track in field. 0142 * @param[in] dydx dydx array. 0143 * @param[in] hstep Proposed step length. 0144 * @param[out] dchord_step Estimated sagitta distance. 0145 * @param[out] dyerr_pos_sq Estimated error in position. 0146 * @param[out] dyerr_mom_rel_sq Estimated error in momentum 0147 * (normalised: Delta_Integration(p^2)/(p^2)). 0148 * @returns true if integration succeeds. 0149 */ 0150 G4bool QuickAdvance(G4FieldTrack& y_posvel, // In/Out 0151 const G4double dydx[], 0152 G4double hstep, // In 0153 G4double& dchord_step, 0154 G4double& dyerr_pos_sq, 0155 G4double& dyerr_mom_rel_sq ); 0156 0157 /** 0158 * Accessors. 0159 */ 0160 inline G4double GetHmin() const; 0161 inline G4double Hmin() const; // Obsolete 0162 inline G4double GetSafety() const; 0163 inline G4double GetPshrnk() const; 0164 inline G4double GetPgrow() const; 0165 inline G4double GetErrcon() const; 0166 void GetDerivatives(const G4FieldTrack& y_curr, // INput 0167 G4double dydx[]) const override; // OUTput 0168 void GetDerivatives(const G4FieldTrack& track, 0169 G4double dydx[], 0170 G4double field[]) const override; 0171 0172 /** 0173 * Getter and setter for the equation of motion. 0174 */ 0175 G4EquationOfMotion* GetEquationOfMotion() override; 0176 void SetEquationOfMotion(G4EquationOfMotion* equation) override; 0177 0178 /** 0179 * Sets a new stepper 'pItsStepper' for this driver. Then it calls 0180 * ResetParameters() to update its parameters accordingly. 0181 */ 0182 void RenewStepperAndAdjust(G4MagIntegratorStepper* pItsStepper) override; 0183 0184 /** 0185 * Resets the qarameters according to the new provided safety value. 0186 * i) sets the exponents (pgrow & pshrnk), using the current order; 0187 * ii) sets the safety and calculates "errcon" according to the above values. 0188 */ 0189 inline void ReSetParameters(G4double new_safety = 0.9); 0190 0191 /** 0192 * Modifiers. When setting safety or pgrow, errcon will be set 0193 * to a compatible value. 0194 */ 0195 inline void SetSafety(G4double valS); 0196 inline void SetPshrnk(G4double valPs); 0197 inline void SetPgrow (G4double valPg); 0198 inline void SetErrcon(G4double valEc); 0199 inline G4double ComputeAndSetErrcon(); 0200 0201 /** 0202 * Accessors for the integrator stepper. 0203 */ 0204 const G4MagIntegratorStepper* GetStepper() const override; 0205 G4MagIntegratorStepper* GetStepper() override; 0206 0207 /** 0208 * Takes one Step that is as large as possible while satisfying the 0209 * accuracy criterion of: yerr < eps * |y_end-y_start|. 0210 * @param[in,out] ystart The current track state, y. 0211 * @param[in] dydx The derivatives array. 0212 * @param[in,out] x Step start, x. 0213 * @param[in] htry Step to attempt. 0214 * @param[in] eps The relative accuracy. 0215 * @param[out] hdid Step achieved. 0216 * @param[out] hnext Proposed next step. 0217 * @returns true if integration succeeds. 0218 */ 0219 void OneGoodStep(G4double ystart[], // Like old RKF45step() 0220 const G4double dydx[], 0221 G4double& x, 0222 G4double htry, 0223 G4double eps, 0224 G4double& hdid, 0225 G4double& hnext ) ; 0226 0227 /** 0228 * Takes the last step's normalised error and calculates a step size 0229 * for the next step. Does it limit the next step's size within a factor 0230 * of the current? 0231 * -- DOES NOT limit for very bad steps 0232 * -- DOES limit for very good (x5). 0233 */ 0234 G4double ComputeNewStepSize(G4double errMaxNorm, // normalised 0235 G4double hstepCurrent) override; 0236 0237 /** 0238 * Taking the last step's normalised error, calculates a step size for 0239 * the next step. Does not limit the next step's size within a factor of 0240 * the current one when *reducing* the size, i.e. for badly failing steps. 0241 */ 0242 G4double ComputeNewStepSize_WithoutReductionLimit(G4double errMaxNorm, 0243 G4double hstepCurrent); 0244 0245 /** 0246 * Taking the last step's normalised error, calculates a step size for 0247 * the next step. Limits the next step's size within a range around the 0248 * current one. 0249 */ 0250 G4double ComputeNewStepSize_WithinLimits(G4double errMaxNorm, // normalised 0251 G4double hstepCurrent); 0252 0253 /** 0254 * Modifier and accessor for the maximum number of steps that can be taken 0255 * for the integration of a single segment, i.e. a single call to 0256 * AccurateAdvance(). 0257 */ 0258 inline G4int GetMaxNoSteps() const; 0259 inline void SetMaxNoSteps(G4int val); 0260 0261 /** 0262 * More modifiers and accessors. 0263 */ 0264 inline void SetHmin(G4double newval); 0265 void SetVerboseLevel(G4int newLevel) override; 0266 G4int GetVerboseLevel() const override; 0267 inline G4double GetSmallestFraction() const; 0268 void SetSmallestFraction( G4double val ); 0269 0270 protected: 0271 0272 /** 0273 * Loggers, issuing warnings for undesirable situations. 0274 */ 0275 void WarnSmallStepSize(G4double hnext, G4double hstep, 0276 G4double h, G4double xDone, 0277 G4int noSteps); 0278 void WarnTooManyStep(G4double x1start, G4double x2end, G4double xCurrent); 0279 void WarnEndPointTooFar(G4double endPointDist, 0280 G4double hStepSize , 0281 G4double epsilonRelative, 0282 G4int debugFlag); 0283 0284 /** 0285 * Loggers for verbosity printouts. 0286 */ 0287 void PrintStatus(const G4double* StartArr, 0288 G4double xstart, 0289 const G4double* CurrentArr, 0290 G4double xcurrent, 0291 G4double requestStep, 0292 G4int subStepNo); 0293 void PrintStatus(const G4FieldTrack& StartFT, 0294 const G4FieldTrack& CurrentFT, 0295 G4double requestStep, 0296 G4int subStepNo); 0297 void PrintStat_Aux(const G4FieldTrack& aFieldTrack, 0298 G4double requestStep, 0299 G4double actualStep, 0300 G4int subStepNo, 0301 G4double subStepSize, 0302 G4double dotVelocities); 0303 /** 0304 * Reports on the number of steps, maximum errors etc. 0305 */ 0306 void PrintStatisticsReport(); 0307 0308 #ifdef QUICK_ADV_TWO 0309 G4bool QuickAdvance( G4double yarrin[], // In 0310 const G4double dydx[], 0311 G4double hstep, 0312 G4double yarrout[], // Out 0313 G4double& dchord_step, // Out 0314 G4double& dyerr ); // in length 0315 #endif 0316 0317 private: 0318 0319 // --------------------------------------------------------------- 0320 // INVARIANTS 0321 0322 /** Minimum Step allowed in a Step (in absolute units). */ 0323 G4double fMinimumStep = 0.0; 0324 0325 /** Smallest fraction of (existing) curve length, in relative units. 0326 Below this fraction the current step will be the last. */ 0327 G4double fSmallestFraction = 1.0e-12; // Expected range 1e-12 to 5e-15 0328 0329 /** Variables in integration. */ 0330 const G4int fNoIntegrationVariables = 0; 0331 0332 /** Minimum number for FieldTrack. */ 0333 const G4int fMinNoVars = 12; 0334 0335 /** Full number of variable. */ 0336 const G4int fNoVars = 0; 0337 0338 /** Default maximum number of steps is Base divided by the order of Stepper. */ 0339 G4int fMaxNoSteps; 0340 G4int fMaxStepBase = 250; // was 5000 0341 0342 /** Parameters used to grow and shrink trial stepsize. */ 0343 G4double safety; 0344 G4double pshrnk; // exponent for shrinking 0345 G4double pgrow; // exponent for growth 0346 G4double errcon; 0347 0348 G4int fStatisticsVerboseLevel = 0; 0349 0350 // --------------------------------------------------------------- 0351 // DEPENDENT Objects 0352 0353 G4MagIntegratorStepper* pIntStepper = nullptr; 0354 0355 // --------------------------------------------------------------- 0356 // STATE 0357 0358 /** Step Statistics. */ 0359 unsigned long fNoTotalSteps=0, fNoBadSteps=0; 0360 unsigned long fNoSmallSteps=0, fNoInitialSmallSteps=0, fNoCalls=0; 0361 G4double fDyerr_max=0.0, fDyerr_mx2=0.0; 0362 G4double fDyerrPos_smTot=0.0, fDyerrPos_lgTot=0.0, fDyerrVel_lgTot=0.0; 0363 G4double fSumH_sm=0.0, fSumH_lg=0.0; 0364 0365 /** Could be varied during tracking - to help identify issues. */ 0366 G4int fVerboseLevel = 0; // Verbosity level for printing (debug, ..) 0367 0368 using ChordFinderDelegate = G4ChordFinderDelegate<G4MagInt_Driver>; 0369 }; 0370 0371 #include "G4MagIntegratorDriver.icc" 0372 0373 #endif
| [ Source navigation ] | [ Diff markup ] | [ Identifier search ] | [ general search ] |
|
This page was automatically generated by the 2.3.7 LXR engine. The LXR team |
|