BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSFieldMagDecapole.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 "BDSFieldMagDecapole.hh"
21#include "BDSMagnetStrength.hh"
22
23#include "globals.hh" // geant4 types / globals
24#include "G4ThreeVector.hh"
25
26#include "CLHEP/Units/SystemOfUnits.h"
27
28#include <cmath>
29
31 G4double const brho)
32{
33 // B'''' = d^4By/dx^4 = Brho * (1/Brho d^4By/dx^4) = Brho * k4
34 bQuadruplePrime = brho * (*strength)["k4"] / (CLHEP::m3*CLHEP::m2);
36#ifdef BDSDEBUG
37 G4cout << __METHOD_NAME__ << "B'''' = " << bQuadruplePrime << G4endl;
38#endif
39}
40
41G4ThreeVector BDSFieldMagDecapole::GetField(const G4ThreeVector &position,
42 const G4double /*t*/) const
43{
44 // B_x = 4xy(x^2-y^2) * ( B''''/4!)
45 // B_y = (x^4 - 6x^2y^2 + y^4) * ( B''''/4!)
46 // B_z = 0
47
48 //shortcuts
49 G4double x = position.x();
50 G4double y = position.y();
51
52 G4ThreeVector localField;
53 localField[0] = 4 * x * y * (std::pow(x,2)-std::pow(y,2)) * bQPNormed;
54 localField[1] = (std::pow(x,4) - 6 * std::pow(x,2) * std::pow(y,2) + std::pow(y,4)) * bQPNormed;
55 localField[2] = 0;
56
57 return localField;
58}
G4double bQuadruplePrime
B'''' - the fourth derivative of the magnetic field.
virtual G4ThreeVector GetField(const G4ThreeVector &position, const G4double t=0) const
Access the field value.
BDSFieldMagDecapole()
Private default constructor to force use of supplied constructor.
G4double bQPNormed
B'''' / 4!
Efficient storage of magnet strengths.