Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/geometry/magneticfield/include/G4IntegrationDriver.hh

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 // G4IntegrationDriver
 27 //
 28 // Class description:
 29 //
 30 // Templated driver class which controls the integration error of a 
 31 // Runge-Kutta stepper.
 32 // It's purpose is to provide a replacement of G4MagIntegratorDriver
 33 // and work for all types of steppers.  
 34 // It will serve as the driver of choice for steppers which do not 
 35 // have extra capabilities, in particular First Same As Last (FSAL)
 36 // and/or interpolation.
 37 
 38 // Author: Dmitry Sorokin, Google Summer of Code 2017
 39 // Supervision: John Apostolakis, CERN
 40 // --------------------------------------------------------------------
 41 #ifndef G4INTEGRATIONDRIVER_HH
 42 #define G4INTEGRATIONDRIVER_HH
 43 
 44 #include "G4RKIntegrationDriver.hh"
 45 #include "G4ChordFinderDelegate.hh"
 46 
 47 template <class T>
 48 class G4IntegrationDriver : public G4RKIntegrationDriver<T>,
 49                             public G4ChordFinderDelegate<G4IntegrationDriver<T>>
 50 {
 51   public:
 52 
 53     G4IntegrationDriver( G4double hminimum,
 54                          T*       stepper,
 55                          G4int    numberOfComponents = 6,
 56                          G4int    statisticsVerbosity = 0 );
 57    ~G4IntegrationDriver() override;
 58 
 59     G4IntegrationDriver(const G4IntegrationDriver &) = delete;
 60     const G4IntegrationDriver& operator =(const G4IntegrationDriver &) = delete;
 61 
 62     G4double AdvanceChordLimited(G4FieldTrack& track, 
 63                                          G4double stepMax, 
 64                                          G4double epsStep,
 65                                          G4double chordDistance) override;
 66 
 67     void OnStartTracking() override;
 68     void OnComputeStep(const G4FieldTrack* /*track*/ = nullptr) override {}
 69     G4bool DoesReIntegrate() const override { return true; }
 70 
 71     G4bool AccurateAdvance(G4FieldTrack& track,
 72                            G4double hstep,
 73                            G4double eps, // Requested y_err/hstep
 74                            G4double hinitial = 0 ) override;
 75       // Integrates ODE from current s (s=s0) to s=s0+h with accuracy eps.
 76       // On output track is replaced by value at end of interval.
 77       // The concept is similar to the odeint routine from NRC p.721-722.
 78 
 79     G4bool QuickAdvance(      G4FieldTrack& fieldTrack,
 80                         const G4double dydx[],
 81                               G4double hstep,
 82                               G4double& dchord_step,
 83                               G4double& dyerr) override;
 84       // QuickAdvance just tries one Step - it does not ensure accuracy.
 85 
 86     void SetVerboseLevel(G4int newLevel) override;
 87     G4int GetVerboseLevel() const override;
 88 
 89     void  StreamInfo( std::ostream& os ) const override;
 90      // Write out the parameters / state of the driver
 91    
 92     // Accessors
 93     //
 94     G4double GetMinimumStep() const;
 95     void SetMinimumStep(G4double newval);
 96 
 97     void OneGoodStep(      G4double  yVar[],  // InOut
 98                      const G4double  dydx[],
 99                            G4double& curveLength,
100                            G4double  htry,
101                            G4double  eps,
102                            G4double& hdid,
103                            G4double& hnext);
104       // This takes one Step that is of size htry, or as large 
105       // as possible while satisfying the accuracy criterion of:
106       //     yerr < eps * |y_end-y_start|
107 
108     G4double GetSmallestFraction() const;
109     void SetSmallestFraction(G4double val);
110 
111   protected:
112 
113     void IncrementQuickAdvanceCalls();
114 
115   private:
116 
117     void CheckStep(const G4ThreeVector& posIn, 
118                    const G4ThreeVector& posOut, 
119                          G4double hdid);
120 
121     G4double fMinimumStep;
122       // Minimum Step allowed in a Step (in absolute units)
123 
124     G4double fSmallestFraction{1e-12};
125       // Smallest fraction of (existing) curve length - in relative units
126       // below this fraction the current step will be the last
127       // Expected range: smaller than 0.1 * epsilon and bigger than 5e-13
128       // Note: this range is not enforced.
129 
130     G4int fVerboseLevel;
131       // Verbosity level for printing (debug, ..)
132       // Could be varied during tracking - to help identify issues
133 
134     G4int fNoQuickAvanceCalls{0};
135     G4int fNoAccurateAdvanceCalls{0};
136     G4int fNoAccurateAdvanceBadSteps{0};
137     G4int fNoAccurateAdvanceGoodSteps{0};
138 
139     using Base = G4RKIntegrationDriver<T>;
140     using ChordFinderDelegate = G4ChordFinderDelegate<G4IntegrationDriver<T>>;
141 };
142 
143 #include "G4IntegrationDriver.icc"
144 
145 #endif
146