Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/geometry/magneticfield/src/G4FieldUtils.cc

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 ]

Diff markup

Differences between /geometry/magneticfield/src/G4FieldUtils.cc (Version 11.3.0) and /geometry/magneticfield/src/G4FieldUtils.cc (Version 11.1.3)


  1 // *******************************************      1 // ********************************************************************
  2 // * License and Disclaimer                         2 // * License and Disclaimer                                           *
  3 // *                                                3 // *                                                                  *
  4 // * The  Geant4 software  is  copyright of th      4 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  5 // * the Geant4 Collaboration.  It is provided      5 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  6 // * conditions of the Geant4 Software License      6 // * conditions of the Geant4 Software License,  included in the file *
  7 // * LICENSE and available at  http://cern.ch/      7 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  8 // * include a list of copyright holders.           8 // * include a list of copyright holders.                             *
  9 // *                                                9 // *                                                                  *
 10 // * Neither the authors of this software syst     10 // * Neither the authors of this software system, nor their employing *
 11 // * institutes,nor the agencies providing fin     11 // * institutes,nor the agencies providing financial support for this *
 12 // * work  make  any representation or  warran     12 // * work  make  any representation or  warranty, express or implied, *
 13 // * regarding  this  software system or assum     13 // * regarding  this  software system or assume any liability for its *
 14 // * use.  Please see the license in the file      14 // * use.  Please see the license in the file  LICENSE  and URL above *
 15 // * for the full disclaimer and the limitatio     15 // * for the full disclaimer and the limitation of liability.         *
 16 // *                                               16 // *                                                                  *
 17 // * This  code  implementation is the result      17 // * This  code  implementation is the result of  the  scientific and *
 18 // * technical work of the GEANT4 collaboratio     18 // * technical work of the GEANT4 collaboration.                      *
 19 // * By using,  copying,  modifying or  distri     19 // * By using,  copying,  modifying or  distributing the software (or *
 20 // * any work based  on the software)  you  ag     20 // * any work based  on the software)  you  agree  to acknowledge its *
 21 // * use  in  resulting  scientific  publicati     21 // * use  in  resulting  scientific  publications,  and indicate your *
 22 // * acceptance of all terms of the Geant4 Sof     22 // * acceptance of all terms of the Geant4 Software license.          *
 23 // *******************************************     23 // ********************************************************************
 24 //                                                 24 //
 25 // Helper namespace field_utils implementation     25 // Helper namespace field_utils implementation
 26 //                                                 26 //
 27 // Author: Dmitry Sorokin, Google Summer of Co     27 // Author: Dmitry Sorokin, Google Summer of Code 2017
 28 // Supervision: John Apostolakis, CERN             28 // Supervision: John Apostolakis, CERN
 29 // -------------------------------------------     29 // --------------------------------------------------------------------
 30                                                    30 
 31 #include "G4FieldUtils.hh"                         31 #include "G4FieldUtils.hh"
 32                                                    32 
 33 #include "G4SystemOfUnits.hh"                      33 #include "G4SystemOfUnits.hh"
 34 #include "G4PhysicalConstants.hh"                  34 #include "G4PhysicalConstants.hh"
 35                                                    35 
 36 namespace field_utils {                            36 namespace field_utils {
 37                                                    37 
 38 G4double absoluteError(const G4double y[],         38 G4double absoluteError(const G4double y[], 
 39                        const G4double yError[]     39                        const G4double yError[],
 40                              G4double hstep)       40                              G4double hstep)
 41 {                                                  41 {
 42    const G4double momentum2 = getValue2(y, Val     42    const G4double momentum2 = getValue2(y, Value3D::Momentum);
 43    const G4double invMomentum2 = 1.0 / momentu     43    const G4double invMomentum2 = 1.0 / momentum2; 
 44    const G4double positionError2 = getValue2(y     44    const G4double positionError2 = getValue2(yError, Value3D::Position);
 45    const G4double momentumError2 = getValue2(y     45    const G4double momentumError2 = getValue2(yError, Value3D::Momentum);
 46    const G4double relativeMomentumError2 = mom     46    const G4double relativeMomentumError2 = momentumError2 * invMomentum2;
 47                                                    47 
 48    return std::max(std::sqrt(positionError2),      48    return std::max(std::sqrt(positionError2),
 49                    std::sqrt(relativeMomentumE     49                    std::sqrt(relativeMomentumError2) * hstep);
 50 }                                                  50 }
 51                                                    51 
 52 G4double relativeError2(const G4double y[],        52 G4double relativeError2(const G4double y[],
 53                         const G4double yerr[],     53                         const G4double yerr[],
 54                               G4double h,          54                               G4double h,
 55                               G4double eps_rel     55                               G4double eps_rel_max)
 56 {                                                  56 {
 57    G4double errmax_sq;                             57    G4double errmax_sq;
 58                                                    58 
 59    G4double inv_eps_vel_sq = 1.0 / (eps_rel_ma     59    G4double inv_eps_vel_sq = 1.0 / (eps_rel_max * eps_rel_max);
 60    G4double errvel_sq = 0.0;    // square of m     60    G4double errvel_sq = 0.0;    // square of momentum vector difference
 61                                                    61 
 62    G4double eps_pos = eps_rel_max * h;             62    G4double eps_pos = eps_rel_max * h; 
 63    G4double inv_eps_pos_sq = 1.0 / (eps_pos *      63    G4double inv_eps_pos_sq = 1.0 / (eps_pos * eps_pos); 
 64                                                    64 
 65    // Evaluate accuracy                            65    // Evaluate accuracy
 66    //                                              66    //
 67    G4double errpos_sq = getValue2(yerr, Value3     67    G4double errpos_sq = getValue2(yerr, Value3D::Position);
 68    errpos_sq *= inv_eps_pos_sq; // Scale relat     68    errpos_sq *= inv_eps_pos_sq; // Scale relative to required tolerance
 69                                                    69 
 70    // Accuracy for momentum                        70    // Accuracy for momentum
 71    //                                              71    //
 72    G4double magvel_sq = getValue2(y, Value3D::     72    G4double magvel_sq = getValue2(y, Value3D::Momentum);
 73    G4double sumerr_sq = getValue2(yerr, Value3     73    G4double sumerr_sq = getValue2(yerr, Value3D::Momentum); 
 74    if (magvel_sq > 0.0)                            74    if (magvel_sq > 0.0)
 75    {                                               75    {
 76       errvel_sq = sumerr_sq / magvel_sq;           76       errvel_sq = sumerr_sq / magvel_sq; 
 77    }                                               77    }
 78    else                                            78    else
 79    {                                               79    {
 80       G4Exception("field_utils::relativeError"     80       G4Exception("field_utils::relativeError","Field001",
 81                   JustWarning, "found case of      81                   JustWarning, "found case of zero momentum");
 82       errvel_sq = sumerr_sq;                       82       errvel_sq = sumerr_sq; 
 83    }                                               83    }
 84    errvel_sq *= inv_eps_vel_sq;                    84    errvel_sq *= inv_eps_vel_sq;
 85    errmax_sq = std::max(errpos_sq, errvel_sq);     85    errmax_sq = std::max(errpos_sq, errvel_sq);
 86                                                    86 
 87    return errmax_sq;                               87    return errmax_sq;
 88 }                                                  88 }
 89                                                    89 
 90 G4double relativeError(const G4double y[],         90 G4double relativeError(const G4double y[],
 91                        const G4double yError[]     91                        const G4double yError[],
 92                        const G4double h,           92                        const G4double h,
 93                        const G4double errorTol     93                        const G4double errorTolerance)
 94 {                                                  94 {
 95    return std::sqrt(relativeError2(y, yError,      95    return std::sqrt(relativeError2(y, yError, h, errorTolerance));
 96 }                                                  96 }
 97                                                    97 
 98 void copy(G4double dst[], const G4double src[]     98 void copy(G4double dst[], const G4double src[], std::size_t size)
 99 {                                                  99 {
100    std::memcpy(dst, src, sizeof(G4double) * si    100    std::memcpy(dst, src, sizeof(G4double) * size);
101 }                                                 101 }
102                                                   102 
103                                                   103 
104 G4double inverseCurvatureRadius(G4double parti    104 G4double inverseCurvatureRadius(G4double particleCharge,
105                                 G4double momen    105                                 G4double momentum,
106                                 G4double BFiel    106                                 G4double BField)
107 {                                                 107 {
108    return -c_light * particleCharge * BField /    108    return -c_light * particleCharge * BField / momentum;
109 }                                                 109 }
110                                                   110 
111 } // field_utils                                  111 } // field_utils
112                                                   112 
113                                                   113