BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSFieldMagDipoleOuterOld.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 "BDSFieldMagDipole.hh"
20#include "BDSFieldMagDipoleOuterOld.hh"
21#include "BDSMagnetStrength.hh"
22
23#include "globals.hh"
24#include "G4ThreeVector.hh"
25
26#include "CLHEP/Units/PhysicalConstants.h"
27#include "CLHEP/Units/SystemOfUnits.h"
28
29#include <cmath>
30#include <vector>
31
32G4double BDSFieldMagDipoleOuterOld::transitionLengthScale = 1*CLHEP::cm;
33
34BDSFieldMagDipoleOuterOld::BDSFieldMagDipoleOuterOld(const BDSMagnetStrength* strength,
35 G4double poleTipRadiusIn,
36 G4double arbitraryScaling):
37 poleTipRadius(poleTipRadiusIn),
38 normalisation(1)
39{
40 BDSFieldMagDipole* innerField = new BDSFieldMagDipole(strength); // encapsulates logic about field parameters
41 // store copy of nominal field strength in vector form
42 localField = innerField->GetField(G4ThreeVector(0,0,0)); // doesn't vary with position
43 m = localField.unit(); // ensure unit vector
44
45 G4ThreeVector normalisationPoint = m*poleTipRadius;
46 G4ThreeVector innerFieldValue = innerField->GetField(normalisationPoint);
47 G4ThreeVector outerFieldValue = GetField(normalisationPoint);
48
49 normalisation = arbitraryScaling * innerFieldValue.mag() / outerFieldValue.mag();
50 if (std::isnan(normalisation))
51 {
52 normalisation = 0; // possible for 0 strength -> b inner = 0 / b outer = 0;
53 finiteStrength = false;
54 }
55
56 delete innerField;
57}
58
59G4ThreeVector BDSFieldMagDipoleOuterOld::GetField(const G4ThreeVector& position,
60 const G4double /*t*/) const
61{
62 if (!finiteStrength)
63 {return G4ThreeVector();} // 0,0,0
64 G4double rmag = position.mag();
65
66 if (rmag < 1) // closer than 1mm from dipole
67 {return localField;}
68
69 // calculate the field according to a magnetic dipole m at position r.
70 G4ThreeVector b = 3*position*(m.dot(position))/std::pow(rmag,5) - m/std::pow(rmag,3);
71
72 b *= normalisation;
73
74 // normalise with tanh (sigmoid-type function) between -3 and 3 in x and 0,1 in y.
75 // scaled spatially in r across the transition length scale
76 // tanh(x) goes from -1 to 1 -> change to 0,1 via +1 and /2
77 G4double weight = (std::tanh(3*(rmag - std::abs(0.5*poleTipRadius))/(transitionLengthScale)) + 1) / 2.0;
78 if (std::abs(weight) > 0.99 || std::abs(weight) < 0.01)
79 {weight = std::round(weight);} // tanh is asymptotic - snap to 1.0 when beyond 0.99
80 b = weight*b + (1-weight)*localField; // weighted sum of two components
81
82 return b;
83}
G4ThreeVector m
Dipole moment as unit vector of field direction.
G4double normalisation
Storage of the overal normalisation factor.
virtual G4ThreeVector GetField(const G4ThreeVector &position, const double t=0) const
Access the field value.
G4double poleTipRadius
Used as radial limit for returning normal field.
G4ThreeVector localField
Nominal dipole field.
A uniform dipole field.
virtual G4ThreeVector GetField(const G4ThreeVector &position, const G4double t=0) const
Access the field value.
G4bool finiteStrength
Flag to cache whether finite nor not.
Definition: BDSFieldMag.hh:83
Efficient storage of magnet strengths.