BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSFieldEMRFCavity.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 "BDSCavityInfo.hh"
20#include "BDSDebug.hh"
21#include "BDSException.hh"
22#include "BDSFieldEMRFCavity.hh"
23#include "BDSMagnetStrength.hh"
24#include "BDSUtilities.hh"
25
26#include "CLHEP/Units/PhysicalConstants.h"
27#include "globals.hh"
28#include "G4ThreeVector.hh"
29
30#include "TMath.h"
31
32#include <cmath>
33#include <utility>
34
35const G4double BDSFieldEMRFCavity::j0FirstZero = 2.404825557695772768622;
36
37const G4double BDSFieldEMRFCavity::Z0 = CLHEP::mu0 * CLHEP::c_light;
38
40 G4double brho):
41 BDSFieldEMRFCavity((*strength)["efield"],
42 (*strength)["frequency"],
43 (*strength)["phase"],
44 (*strength)["equatorradius"])
45{
46 eFieldMax *= BDS::Sign(brho);
47}
48
49BDSFieldEMRFCavity::BDSFieldEMRFCavity(G4double eFieldAmplitude,
50 G4double frequencyIn,
51 G4double phaseOffset,
52 G4double cavityRadiusIn):
53 eFieldMax(eFieldAmplitude),
54 phase(phaseOffset),
55 cavityRadius(cavityRadiusIn),
56 wavelength(CLHEP::c_light / frequencyIn),
57 normalisedCavityRadius(j0FirstZero/cavityRadius),
58 angularFrequency(CLHEP::twopi * frequencyIn)
59{
60 // this would cause NANs to be propagated into tracking which is really bad
61 if (!BDS::IsFinite(cavityRadiusIn) || std::isnan(normalisedCavityRadius) || std::isinf(normalisedCavityRadius))
62 {throw BDSException(__METHOD_NAME__, "no cavity radius supplied - required for pill box model");}
63}
64
65std::pair<G4ThreeVector, G4ThreeVector> BDSFieldEMRFCavity::GetField(const G4ThreeVector& position,
66 const G4double t) const
67{
68 // Converting from Local Cartesian to Local Cylindrical
69 G4double phi = std::atan2(position.y(),position.x());
70 G4double r = std::hypot(position.x(),position.y());
71
72 G4double rNormalised = normalisedCavityRadius * r;
73
74 // In case a point outside the cavity is queried, ensure the bessel will return 0
75 if (rNormalised > j0FirstZero)
76 {rNormalised = j0FirstZero - 1e-6;}
77
78 // Source for calculating the TM010 mode: Gerigk, Frank.
79 // "Cavity types." arXiv preprint arXiv:1111.4897 (2011).
80
81 G4double J0r = TMath::BesselJ0(rNormalised);
82 G4double J1r = TMath::BesselJ1(rNormalised);
83
84 // Calculating free-space impedance and scale factor for Bphi:
85 G4double hMax = -eFieldMax/Z0;
86 G4double Bmax = hMax * CLHEP::mu0;
87
88 // Calculating field components.
89 G4double zFactor = std::cos(CLHEP::twopi*position.z() / wavelength);
90 G4double Ez = eFieldMax * J0r * std::cos(angularFrequency*t + phase) * zFactor;
91 G4double Bphi = Bmax * J1r * std::sin(angularFrequency*t + phase) * zFactor;
92
93 // Converting Bphi into cartesian coordinates:
94 G4TwoVector bxby(0,Bphi); // this is equivalent to a pi/2 rotation of (1,0)
95 bxby.rotate(phi);
96 G4double Bx = bxby.x();
97 G4double By = bxby.y();
98
99 // Local B and E fields:
100 G4ThreeVector LocalB = G4ThreeVector(Bx, By, 0);
101 G4ThreeVector LocalE = G4ThreeVector(0, 0, Ez);
102
103 auto result = std::make_pair(LocalB, LocalE);
104 return result;
105}
General exception with possible name of object and message.
Definition: BDSException.hh:35
Pill box cavity electromagnetic field.
virtual std::pair< G4ThreeVector, G4ThreeVector > GetField(const G4ThreeVector &position, const G4double t) const
Accessor to get B and E field.
const G4double angularFrequency
Angular frequency calculated from frequency - cached to avoid repeated calculation.
G4double phase
Phase offset of the oscillator.
BDSFieldEMRFCavity()
Private constructor to force use of provided one.
static const G4double j0FirstZero
X coordinate of first 0 point for bessel J0.
static const G4double Z0
Impedance of free space.
G4double eFieldMax
Maximum field in V/m.
const G4double normalisedCavityRadius
Pre-calculated normalised calculated radius w.r.t. bessel first 0.
Efficient storage of magnet strengths.
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())