BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
BDSIntegratorSet.cc
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2022.
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 "BDSIntegratorSet.hh"
21#include "BDSIntegratorSetType.hh"
22#include "BDSIntegratorType.hh"
23
24#include "globals.hh" // geant4 types / globals
25#include "G4Version.hh"
26
27BDSIntegratorSet::BDSIntegratorSet(BDSIntegratorType solenoidIn,
28 BDSIntegratorType dipoleIn,
29 BDSIntegratorType dipoleQuadrupoleIn,
30 BDSIntegratorType quadrupoleIn,
31 BDSIntegratorType sextupoleIn,
32 BDSIntegratorType octupoleIn,
33 BDSIntegratorType decapoleIn,
34 BDSIntegratorType multipoleThickIn,
35 BDSIntegratorType muonSpoilerIn,
36 BDSIntegratorType rfcavityIn,
38 BDSIntegratorType generalIn,
39 BDSIntegratorType skewQuadrupoleIn,
40 BDSIntegratorType skewSextupoleIn,
41 BDSIntegratorType skewOctupoleIn,
42 BDSIntegratorType skewDecapoleIn,
43 BDSIntegratorType dipoleFringeIn,
44 BDSIntegratorType multipoleThinIn,
45 BDSIntegratorType multipoleOuterIn,
46 BDSIntegratorType rmatrixThinIn,
47 BDSIntegratorType parallelTransporterIn,
48 BDSIntegratorType undulatorIn,
49 BDSIntegratorType cavityFringeIn):
50 solenoid(solenoidIn),
51 dipole(dipoleIn),
52 dipoleQuadrupole(dipoleQuadrupoleIn),
53 quadrupole(quadrupoleIn),
54 sextupole(sextupoleIn),
55 octupole(octupoleIn),
56 decapole(decapoleIn),
57 multipoleThick(multipoleThickIn),
58 muonSpoiler(muonSpoilerIn),
59 rfcavity(rfcavityIn),
60 rf(rfIn),
61 general(generalIn),
62 skewQuadrupole(skewQuadrupoleIn),
63 skewSextupole(skewSextupoleIn),
64 skewOctupole(skewOctupoleIn),
65 skewDecapole(skewDecapoleIn),
66 dipoleFringe(dipoleFringeIn),
67 multipoleThin(multipoleThinIn),
68 multipoleOuter(multipoleOuterIn),
69 rmatrixThin(rmatrixThinIn),
70 parallelTransporter(parallelTransporterIn),
71 undulator(undulatorIn),
72 cavityFringe(cavityFringeIn)
73{
74 isMatrix = false; //default
75 // use dipolematrix integrator to check if matrix style as it is the
76 // only integrator exclusive to matrix style integrator sets
77 if (dipoleIn == BDSIntegratorType::dipolematrix)
78 {isMatrix = true;}
79}
80
81namespace BDS
82{
83 const BDSIntegratorSet* integratorsBDSIMOne =
84 new BDSIntegratorSet(BDSIntegratorType::solenoid, // solenoid
85 BDSIntegratorType::dipolerodrigues, // dipole
86 BDSIntegratorType::dipolematrix, // dipole quadrupole
87 BDSIntegratorType::quadrupole, // quadrupole
88 BDSIntegratorType::sextupole, // sextupole
89 BDSIntegratorType::octupole, // octupole
90 BDSIntegratorType::decapole, // decapole
91 BDSIntegratorType::g4classicalrk4, // thick multipole
92 BDSIntegratorType::g4classicalrk4, // muon spoiler
93 BDSIntegratorType::g4classicalrk4, // rfcavity
94 BDSIntegratorType::g4classicalrk4, // rf
95 BDSIntegratorType::g4classicalrk4, // general
96 BDSIntegratorType::g4classicalrk4, // skew quadrupole
97 BDSIntegratorType::g4classicalrk4, // skew sextupole
98 BDSIntegratorType::g4classicalrk4, // skew octupole
99 BDSIntegratorType::g4classicalrk4, // skew decapole
100 BDSIntegratorType::dipolefringe, // dipole fringe field
101 BDSIntegratorType::multipolethin, // thin multipole
102 BDSIntegratorType::g4rk4minimumstep, // multipole outer
103 BDSIntegratorType::rmatrixthin, // thin rmatrix
104 BDSIntegratorType::paralleltransport,// parallel transport
105 BDSIntegratorType::g4classicalrk4, // undulator
106 BDSIntegratorType::cavityfringe); // cavity fringe
107
108 const BDSIntegratorSet* integratorsBDSIMTwo =
109 new BDSIntegratorSet(BDSIntegratorType::solenoid, // solenoid
110 BDSIntegratorType::dipolerodrigues2, // dipole
111 BDSIntegratorType::dipolematrix, // dipole quadrupole
112 BDSIntegratorType::quadrupole, // quadrupole
113 BDSIntegratorType::euler, // sextupole
114 BDSIntegratorType::euler, // octupole
115 BDSIntegratorType::euler, // decapole
116 BDSIntegratorType::g4classicalrk4, // (thick) multipole
117 BDSIntegratorType::g4classicalrk4, // muonspoiler
118 BDSIntegratorType::g4classicalrk4, // rfcavity
119 BDSIntegratorType::g4classicalrk4, // rf
120 BDSIntegratorType::g4classicalrk4, // general
121 BDSIntegratorType::g4classicalrk4, // skew quadrupole
122 BDSIntegratorType::g4classicalrk4, // skew sextupole
123 BDSIntegratorType::g4classicalrk4, // skew octupole
124 BDSIntegratorType::g4classicalrk4, // skew decapole
125 BDSIntegratorType::dipolefringe, // dipole fringe field
126 BDSIntegratorType::multipolethin, // thin multipole
127 BDSIntegratorType::g4rk4minimumstep, // multipole outer
128 BDSIntegratorType::rmatrixthin, // thin rmatrix
129 BDSIntegratorType::paralleltransport,// parallel transport
130 BDSIntegratorType::g4classicalrk4, // undulator
131 BDSIntegratorType::cavityfringe); // cavity fringe
134 new BDSIntegratorSet(BDSIntegratorType::solenoid, // solenoid
135 BDSIntegratorType::dipolematrix, // dipole
136 BDSIntegratorType::dipolematrix, // dipole quadrupole
137 BDSIntegratorType::quadrupole, // quadrupole
138 BDSIntegratorType::euler, // sextupole
139 BDSIntegratorType::euler, // octupole
140 BDSIntegratorType::euler, // decapole
141 BDSIntegratorType::g4classicalrk4, // thick multipole
142 BDSIntegratorType::g4classicalrk4, // muon spoiler
143 BDSIntegratorType::g4classicalrk4, // rfcavity
144 BDSIntegratorType::g4classicalrk4, // rf
145 BDSIntegratorType::g4classicalrk4, // general
146 BDSIntegratorType::g4classicalrk4, // skew quadrupole
147 BDSIntegratorType::g4classicalrk4, // skew sextupole
148 BDSIntegratorType::g4classicalrk4, // skew octupole
149 BDSIntegratorType::g4classicalrk4, // skew decapole
150 BDSIntegratorType::dipolefringe, // dipole fringe field
151 BDSIntegratorType::multipolethin, // thin multipole
152 BDSIntegratorType::g4rk4minimumstep, // multipole outer
153 BDSIntegratorType::rmatrixthin, // thin rmatrix
154 BDSIntegratorType::paralleltransport, // parallel transport
155 BDSIntegratorType::g4classicalrk4, // undulator
156 BDSIntegratorType::cavityfringe); // cavity fringe
157 const BDSIntegratorSet* integratorsBDSIMMatrixFringeScaling =
158 new BDSIntegratorSet(BDSIntegratorType::solenoid, // solenoid
159 BDSIntegratorType::dipolematrix, // dipole
160 BDSIntegratorType::dipolematrix, // dipole quadrupole
161 BDSIntegratorType::quadrupole, // quadrupole
162 BDSIntegratorType::euler, // sextupole
163 BDSIntegratorType::euler, // octupole
164 BDSIntegratorType::euler, // decapole
165 BDSIntegratorType::g4classicalrk4, // thick multipole
166 BDSIntegratorType::g4classicalrk4, // muon spoiler
167 BDSIntegratorType::g4classicalrk4, // rfcavity
168 BDSIntegratorType::g4classicalrk4, // rf
169 BDSIntegratorType::g4classicalrk4, // general
170 BDSIntegratorType::g4classicalrk4, // skew quadrupole
171 BDSIntegratorType::g4classicalrk4, // skew sextupole
172 BDSIntegratorType::g4classicalrk4, // skew octupole
173 BDSIntegratorType::g4classicalrk4, // skew decapole
174 BDSIntegratorType::dipolefringescaling, // dipole fringe field
175 BDSIntegratorType::multipolethin, // thin multipole
176 BDSIntegratorType::g4rk4minimumstep, // multipole outer - nystrom doesn't work in g4.10.5
177 BDSIntegratorType::rmatrixthin, // thin rmatrix
178 BDSIntegratorType::paralleltransport, // parallel transport
179 BDSIntegratorType::g4classicalrk4, // undulator
180 BDSIntegratorType::cavityfringe); // cavity fringe
181 const BDSIntegratorSet* integratorsGeant4 =
182 new BDSIntegratorSet(BDSIntegratorType::g4classicalrk4, // solenoid
183 BDSIntegratorType::g4classicalrk4, // dipole
184 BDSIntegratorType::g4classicalrk4, // dipole quadrupole
185 BDSIntegratorType::g4classicalrk4, // quadrupole
186 BDSIntegratorType::g4classicalrk4, // sextupole
187 BDSIntegratorType::g4classicalrk4, // octupole
188 BDSIntegratorType::g4classicalrk4, // decapole
189 BDSIntegratorType::g4classicalrk4, // thick multipole
190 BDSIntegratorType::g4classicalrk4, // muon spoiler
191 BDSIntegratorType::g4classicalrk4, // rfcavity
192 BDSIntegratorType::g4classicalrk4, // rf
193 BDSIntegratorType::g4classicalrk4, // general
194 BDSIntegratorType::g4classicalrk4, // skew quadrupole
195 BDSIntegratorType::g4classicalrk4, // skew sextupole
196 BDSIntegratorType::g4classicalrk4, // skew octupole
197 BDSIntegratorType::g4classicalrk4, // skew decapole
198 BDSIntegratorType::dipolefringe, // dipole fringe field
199 BDSIntegratorType::multipolethin, // thin multipole
200 BDSIntegratorType::g4rk4minimumstep, // multipole outer
201 BDSIntegratorType::rmatrixthin, // thin rmatrix
202 BDSIntegratorType::paralleltransport, // parallel transport
203 BDSIntegratorType::g4classicalrk4, // undulator
204 BDSIntegratorType::cavityfringe); // cavity fringe
205#if G4VERSION_NUMBER > 1039
206 const BDSIntegratorSet* integratorsGeant4DP =
207 new BDSIntegratorSet(BDSIntegratorType::g4dormandprince745, // solenoid
208 BDSIntegratorType::g4dormandprince745, // dipole
209 BDSIntegratorType::g4dormandprince745, // dipole quadrupole
210 BDSIntegratorType::g4dormandprince745, // quadrupole
211 BDSIntegratorType::g4dormandprince745, // sextupole
212 BDSIntegratorType::g4dormandprince745, // octupole
213 BDSIntegratorType::g4dormandprince745, // decapole
214 BDSIntegratorType::g4dormandprince745, // thick multipole
215 BDSIntegratorType::g4dormandprince745, // muon spoiler
216 BDSIntegratorType::g4dormandprince745, // rfcavity
217 BDSIntegratorType::g4dormandprince745, // rf
218 BDSIntegratorType::g4dormandprince745, // general
219 BDSIntegratorType::g4dormandprince745, // skew quadrupole
220 BDSIntegratorType::g4dormandprince745, // skew sextupole
221 BDSIntegratorType::g4dormandprince745, // skew octupole
222 BDSIntegratorType::g4dormandprince745, // skew decapole
223 BDSIntegratorType::dipolefringe, // dipole fringe field
224 BDSIntegratorType::multipolethin, // thin multipole
225 BDSIntegratorType::g4dormandprince745, // multipole outer
226 BDSIntegratorType::rmatrixthin, // thin rmatrix
227 BDSIntegratorType::paralleltransport, // parallel transport
228 BDSIntegratorType::g4dormandprince745, // undulator
229 BDSIntegratorType::cavityfringe); // cavity fringe
230#endif
231}
232
234{
236}
237
239{
240 switch (set.underlying())
241 {
242 case BDSIntegratorSetType::geant4:
243 {return BDS::integratorsGeant4; break;}
244#if G4VERSION_NUMBER > 1039
245 case BDSIntegratorSetType::geant4dp:
246 {return BDS::integratorsGeant4DP; break;}
247#endif
248 case BDSIntegratorSetType::bdsimone:
249 {return BDS::integratorsBDSIMOne; break;}
250 case BDSIntegratorSetType::bdsimtwo:
251 {return BDS::integratorsBDSIMTwo;}
252 case BDSIntegratorSetType::bdsimmatrix:
254 case BDSIntegratorSetType::bdsimmatrixfringescaling:
255 {return BDS::integratorsBDSIMMatrixFringeScaling;}
256 default:
257 {return BDS::integratorsBDSIMOne; break;}
258 }
259}
260
262{
263 switch (field.underlying())
264 {
265 case BDSFieldType::none:
266 {return general; break;}// shouldn't really happen, but for completeness.
267 case BDSFieldType::bmap1d:
268 case BDSFieldType::bmap2d:
269 case BDSFieldType::bmap3d:
270 case BDSFieldType::bmap4d:
271 case BDSFieldType::ebmap1d:
272 case BDSFieldType::ebmap2d:
273 case BDSFieldType::ebmap3d:
274 case BDSFieldType::ebmap4d:
275 case BDSFieldType::emap1d:
276 case BDSFieldType::emap2d:
277 case BDSFieldType::emap3d:
278 case BDSFieldType::emap4d:
279 case BDSFieldType::mokka:
280 {
281 G4cout << __METHOD_NAME__ << "WARNING - this is overriding the specified field maps integrator" << G4endl;
282 return general;
283 break;
284 }
285 case BDSFieldType::solenoid:
286 {return solenoid; break;}
287 case BDSFieldType::dipole:
288 {return dipole; break;}
289 case BDSFieldType::dipolequadrupole:
290 {return dipoleQuadrupole; break;}
291 case BDSFieldType::quadrupole:
292 {return quadrupole; break;}
293 case BDSFieldType::sextupole:
294 {return sextupole; break;}
295 case BDSFieldType::octupole:
296 {return octupole; break;}
297 case BDSFieldType::decapole:
298 {return decapole; break;}
299 case BDSFieldType::multipole:
300 {return multipoleThick; break;}
301 case BDSFieldType::muonspoiler:
302 {return muonSpoiler; break;}
303 case BDSFieldType::skewquadrupole:
304 {return skewQuadrupole; break;}
305 case BDSFieldType::skewsextupole:
306 {return skewSextupole; break;}
307 case BDSFieldType::skewoctupole:
308 {return skewOctupole; break;}
309 case BDSFieldType::skewdecapole:
310 {return skewDecapole; break;}
311 case BDSFieldType::rfcavity:
312 {return rfcavity; break;}
313 case BDSFieldType::rf:
314 {return rf; break;}
315 case BDSFieldType::rmatrix:
316 {return rmatrixThin; break;}
317 case BDSFieldType::paralleltransporter:
318 {return parallelTransporter; break;}
319 case BDSFieldType::cavityfringe:
320 {return cavityFringe; break;}
321 case BDSFieldType::undulator:
322 {return undulator; break;}
323 case BDSFieldType::dipole3d:
324 {return general; break;}
325 case BDSFieldType::multipoleouterdipole:
326 case BDSFieldType::multipoleouterquadrupole:
327 case BDSFieldType::multipoleoutersextupole:
328 case BDSFieldType::multipoleouteroctupole:
329 case BDSFieldType::multipoleouterdecapole:
330 case BDSFieldType::skewmultipoleouterquadrupole:
331 case BDSFieldType::skewmultipoleoutersextupole:
332 case BDSFieldType::skewmultipoleouteroctupole:
333 case BDSFieldType::skewmultipoleouterdecapole:
334 case BDSFieldType::multipoleouterdipole3d:
335 case BDSFieldType::multipoleouterdipolelhc:
336 case BDSFieldType::multipoleouterquadrupolelhc:
337 case BDSFieldType::multipoleoutersextupolelhc:
338 {return multipoleOuter; break;}
339 default:
340 {return general; break;}
341 }
342}
Which integrator to use for each type of magnet / field object.
BDSIntegratorType Integrator(const BDSFieldType field) const
Get appropriate integrator based on the field type.
type underlying() const
return underlying value (can be used in switch statement)
Return either G4Tubs or G4CutTubs depending on flat face.
const BDSIntegratorSet * IntegratorSet(G4String set)
Return the appropriate set of integrators to use for each magnet type.
BDSIntegratorSetType DetermineIntegratorSetType(G4String integratorSet)
Function that gives corresponding enum value for string (case-insensitive)
const BDSIntegratorSet * integratorsBDSIMMatrix
Mad-x style tracking.