include/Geant4/G4DormandPrinceRK56.hh

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 // G4DormandPrinceRK56
0027 //
0028 // Class description:
0029 //
0030 // Dormand-Prince RK 6(5) non-FSAL method
0031
0032 // Created: Somnath Banerjee, Google Summer of Code 2015, 26 June 2015
0033 // Supervision: John Apostolakis, CERN
0034 // --------------------------------------------------------------------
0035 #ifndef G4DORMAND_PRINCE_RK56_HH
0036 #define G4DORMAND_PRINCE_RK56_HH
0037
0038 #include "G4MagIntegratorStepper.hh"
0039
0040 class G4DormandPrinceRK56 : public G4MagIntegratorStepper
0041 {
0042   public:
0043
0044     G4DormandPrinceRK56( G4EquationOfMotion* EqRhs,
0045                          G4int numberOfVariables = 6,
0046                          G4bool primary = true ) ;
0047
0048     ~G4DormandPrinceRK56() override ;
0049
0050     G4DormandPrinceRK56(const G4DormandPrinceRK56&) = delete;
0051     G4DormandPrinceRK56& operator=(const G4DormandPrinceRK56&) = delete;
0052
0053     void Stepper( const G4double y[],
0054                   const G4double dydx[],
0055                         G4double h,
0056                         G4double yout[],
0057                         G4double yerr[] ) override ;
0058
0059     G4double  DistChord()   const override;
0060     G4int IntegratorOrder() const override { return 5; }
0061
0062     void SetupInterpolate_low( const G4double yInput[],
0063                                const G4double dydx[],
0064                                const G4double Step );
0065       // For preparing the Interpolant and calculating the extra stages
0066
0067     void Interpolate_low( const G4double yInput[],
0068                           const G4double dydx[],
0069                           const G4double Step,
0070                                 G4double yOut[],
0071                                 G4double tau );
0072       // For calculating the output at the tau fraction of Step
0073
0074     inline void SetupInterpolation()
0075     {
0076       SetupInterpolate( fLastInitialVector, fLastDyDx, fLastStepLength);
0077     }
0078
0079     inline void SetupInterpolate( const G4double yInput[],
0080                                   const G4double dydx[],
0081                                   const G4double Step )
0082     {
0083       SetupInterpolate_low( yInput, dydx, Step);
0084     }
0085
0086     inline void Interpolate( const G4double yInput[],
0087                              const G4double dydx[],
0088                              const G4double Step,
0089                                    G4double yOut[],
0090                                    G4double tau )
0091     {
0092       Interpolate_low( yInput, dydx, Step, yOut, tau);
0093     }
0094       // For calculating the output at the tau fraction of Step
0095
0096     inline void Interpolate( G4double tau, G4double yOut[])
0097     {
0098       Interpolate( fLastInitialVector, fLastDyDx, fLastStepLength, yOut, tau );
0099     }
0100
0101     void SetupInterpolate_high( const G4double yInput[],
0102                                 const G4double dydx[],
0103                                 const G4double Step );
0104
0105     void Interpolate_high( const G4double yInput[],
0106                            const G4double dydx[],
0107                            const G4double Step,
0108                                  G4double yOut[],
0109                                  G4double tau );
0110       // For calculating the output at the tau fraction of Step
0111
0112   private:
0113
0114     G4double *ak2, *ak3, *ak4, *ak5, *ak6, *ak7, *ak8, *ak9;
0115       // For storing intermediate 'k' values in stepper
0116     G4double *ak10_low, *ak10, *ak11, * ak12;
0117       // For the additional stages of Interpolant
0118     G4double *yTemp, *yIn;
0119
0120     G4double fLastStepLength = -1.0;
0121     G4double *fLastInitialVector, *fLastFinalVector,
0122              *fLastDyDx, *fMidVector, *fMidError;
0123       // For DistChord calculations
0124
0125     G4DormandPrinceRK56* fAuxStepper = nullptr;
0126 };
0127
0128 #endif /* G4DormandPrinceRK56 */