BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSMagUsualEqRhs.cc
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2023.
4
5This file is part of BDSIM.
6
7BDSIM is free software: you can redistribute it and/or modify
8it under the terms of the GNU General Public License as published
9by the Free Software Foundation version 3 of the License.
10
11BDSIM is distributed in the hope that it will be useful, but
12WITHOUT ANY WARRANTY; without even the implied warranty of
13MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14GNU General Public License for more details.
15
16You should have received a copy of the GNU General Public License
17along with BDSIM. If not, see <http://www.gnu.org/licenses/>.
18*/
19#include "BDSGlobalConstants.hh"
20#include "BDSMagUsualEqRhs.hh"
21#include "BDSStep.hh"
22#include "G4Mag_EqRhs.hh"
23
24#include "globals.hh" // geant4 types / globals
25#include "G4ClassicalRK4.hh"
26#include "G4ChargeState.hh"
27#include "G4MagneticField.hh"
28#include "G4SystemOfUnits.hh"
29#include "G4PhysicalConstants.hh"
30#include "BDSGlobalConstants.hh"
31
32BDSMagUsualEqRhs::BDSMagUsualEqRhs(G4MagneticField* MagField):
33 G4Mag_UsualEqRhs(MagField),
34 fMassCof(0.),
35 fCof_val(0.)
36{;}
37
38void BDSMagUsualEqRhs::SetChargeMomentumMass( G4ChargeState particleCharge,
39 G4double MomentumXc,
40 G4double particleMass)
41{
42 G4Mag_EqRhs::SetChargeMomentumMass(particleCharge, MomentumXc, particleMass);
43 G4double pcharge = particleCharge.GetCharge();
44
45 fCof_val = pcharge*eplus*c_light ; // B must be in Tesla
46 fMassCof = particleMass*particleMass;
47}
48
49G4double BDSMagUsualEqRhs::Beta(const G4double y[6])
50{
51 G4ThreeVector mom = G4ThreeVector(y[3], y[4], y[5]);
52 G4double beta = Beta(mom);
53 return beta;
54}
55
56G4double BDSMagUsualEqRhs::Beta(const G4ThreeVector& mom)
57{
58 G4double momentum_mag_square = mom.mag2();
59 G4double Energy = std::sqrt(momentum_mag_square + fMassCof);
60 G4double beta = std::sqrt(momentum_mag_square) / Energy;
61 return beta;
62}
63
64G4double BDSMagUsualEqRhs::TotalEnergy(const G4double y[])
65{
66 G4ThreeVector mom = G4ThreeVector(y[3], y[4], y[5]);
67 G4double Energy = TotalEnergy(mom);
68 return Energy;
69}
70
71G4double BDSMagUsualEqRhs::TotalEnergy(const G4ThreeVector& mom)
72{
73 G4double momentum_mag_square = mom.mag2();
74 G4double Energy = std::sqrt(momentum_mag_square + fMassCof);
75 return Energy;
76}
G4double fMassCof
Particle mass squared.
G4double TotalEnergy(const G4double y[])
Calculate total particle energy.
G4double Beta(const G4double y[6])
Calculate particle velocity W.R.T speed of light.
virtual void SetChargeMomentumMass(G4ChargeState particleCharge, G4double MomentumXc, G4double particleMass)
Copy of class method from G4Mag_UsualEqRhs.