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