BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSFieldMagSextupole.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 "BDSFieldMagSextupole.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/SystemOfUnits.h"
28
29#include <cmath>
30
32 G4double const brho)
33{
34 // B'' = d^2By/dx^2 = Brho * (1/Brho d^2By/dx^2) = Brho * k2
35 bDoublePrime = brho * (*strength)["k2"] / CLHEP::m3;
38#ifdef BDSDEBUG
39 G4cout << __METHOD_NAME__ << "B'' = " << bDoublePrime << G4endl;
40#endif
41}
42
43G4ThreeVector BDSFieldMagSextupole::GetField(const G4ThreeVector& position,
44 const G4double /*t*/) const
45{
46 // B_x = 2*x*y * (B''/2!)
47 // B_y = (x^2 - y^2) * (B''/2!)
48 // B_z = 0
49
50 G4ThreeVector localField;
51 localField[0] = position.x() * position.y() * bDoublePrime;
52 localField[1] = (std::pow(position.x(),2) - std::pow(position.y(),2)) * halfBDoublePrime;
53 localField[2] = 0;
54
55 return localField;
56}
BDSFieldMagSextupole()
Private default constructor to avoid usage.
G4double bDoublePrime
B'' - second derivative of the magnetic field.
virtual G4ThreeVector GetField(const G4ThreeVector &position, const G4double t=0) const
Access the field value.
G4double halfBDoublePrime
bDoublePrime / 2! cached
G4bool finiteStrength
Flag to cache whether finite nor not.
Definition: BDSFieldMag.hh:83
Efficient storage of magnet strengths.
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())