Geant4 Cross Reference

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