include/Geant4/G4TsitourasRK45.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 // G4TsitourasRK45
0027 //
0028 // Class description:
0029 //
0030 // Tsitouras - 5(4) RK stepper (non-FSAL version)
0031 // Implements RK tableau from 'Table 1' of:
0032 //    C. Tsitouras, "Runge-Kutta pairs of order 5(4) satisfying only
0033 //    the first column simplifying assumption"
0034 //    Computers & Mathematics with Applications,
0035 //    vol. 62, no. 2, pp. 770-775, 2011.
0036
0037 // Author: Somnath Banerjee, Google Summer of Code 2015, 11.06.2015
0038 // Supervision: John Apostolakis, CERN
0039 // --------------------------------------------------------------------
0040 #ifndef TSITOURAS_RK45_HH
0041 #define TSITOURAS_RK45_HH
0042
0043 #include "G4MagIntegratorStepper.hh"
0044
0045 class G4TsitourasRK45 : public G4MagIntegratorStepper
0046 {
0047   public:
0048
0049     G4TsitourasRK45(G4EquationOfMotion* EqRhs,
0050                     G4int numberOfVariables = 6,
0051                     G4bool primary =  true);
0052    ~G4TsitourasRK45() override;
0053
0054     G4TsitourasRK45(const G4TsitourasRK45&) = delete;
0055     G4TsitourasRK45& operator=(const G4TsitourasRK45&) = delete;
0056
0057     void Stepper( const G4double y[],
0058                   const G4double dydx[],
0059                         G4double h,
0060                         G4double yout[],
0061                         G4double yerr[] ) override;
0062
0063     void SetupInterpolation( /* const G4double yInput[],
0064                               const G4double dydx[],
0065                               const G4double Step */  );
0066
0067     void Interpolate( 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     void interpolate( const G4double yInput[],
0075                       const G4double dydx[],
0076                             G4double yOut[],
0077                             G4double Step,
0078                             G4double tau);
0079
0080     G4double DistChord() const override;
0081     inline G4int IntegratorOrder() const override { return 4; }
0082
0083   private :
0084
0085     G4double *ak2, *ak3, *ak4, *ak5, *ak6, *ak7, *ak8, *yTemp, *yIn;
0086
0087     G4double fLastStepLength = 0.0;
0088     G4double *fLastInitialVector, *fLastFinalVector,
0089              *fLastDyDx, *fMidVector, *fMidError;
0090       // For DistChord calculations
0091
0092     G4TsitourasRK45* fAuxStepper = nullptr;
0093 };
0094
0095 #endif