BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSFieldMagUndulator.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 "BDSDebug.hh"
20#include "BDSFieldMagUndulator.hh"
21#include "BDSMagnetStrength.hh"
22#include "BDSUtilities.hh"
23
24#include "globals.hh" // geant4 types / globals
25#include "G4ThreeVector.hh"
26
27#include "CLHEP/Units/PhysicalConstants.h"
28#include "CLHEP/Units/SystemOfUnits.h"
29
30#include <cmath>
31
33 G4double beamPipeRadiusIn):
34 beamPipeRadius(beamPipeRadiusIn)
35{
36 wavenumber = (CLHEP::twopi)/(*strength)["length"];
37 B = (*strength)["field"];// / CLHEP::tesla;
38 finiteStrength = BDS::IsFinite(B);
39 limit = 3 * B; // limit for validity of this
40#ifdef BDSDEBUG
41 G4cout << __METHOD_NAME__ << "B = " << B << G4endl;
42 G4cout << __METHOD_NAME__ << "B_limit = " << limit << G4endl;
43#endif
44}
45
46G4ThreeVector BDSFieldMagUndulator::GetField(const G4ThreeVector& position,
47 const G4double /*t*/) const
48{
49 G4ThreeVector field;
50
51 if (std::abs(position.x()) > beamPipeRadius || std::abs(position.y()) > beamPipeRadius)
52 {return G4ThreeVector(0,0,0);}
53
54 G4double yFactor = std::cosh(position.y() * wavenumber);
55 if (yFactor > 5)
56 {yFactor = 5;}
57 G4double zFactor = std::sinh(position.y() * wavenumber);
58 if (zFactor > 5)
59 {zFactor = 5;}
60
61 field[0] = 0;
62 field[1] = B * std::cos(position.z() * wavenumber) * yFactor;
63 field[2] = -B * std::sin(position.z() * wavenumber) * zFactor;
64
65 // hyperbolic functions can tend to infinity if we happen to query a
66 // very large point limits -> a) the max magnitude and b) no nans
67 if (field[1] > limit || std::isnan(field[1]))
68 {field[1] = limit;}
69 if (field[2] > limit || std::isnan(field[2]))
70 {field[2] = limit;}
71
72 return field;
73}
G4double wavenumber
The undulator wavenumber.
BDSFieldMagUndulator()
Private default constructor to force use of supplied constructor.
G4double beamPipeRadius
Cache of beam pipe radius to know maximum valid extent of field.
virtual G4ThreeVector GetField(const G4ThreeVector &position, const G4double t=0) const
Access the field value.
G4double B
The peak field.
Efficient storage of magnet strengths.
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())