OPAL (Object Oriented Parallel Accelerator Library) 2022.1
OPAL
MapGenerator.h
Go to the documentation of this file.
1
10#ifndef MAPGENERATOR_H
11#define MAPGENERATOR_H
12
13#include <cmath>
14#include <stdexcept>
15#include <vector>
16
17#include <boost/numeric/ublas/matrix.hpp>
18
20
21#include "FixedAlgebra/FTps.h"
23
24#include "Physics/Physics.h"
26
28template<typename Value_type, typename Size_type, typename Series_type,
29 typename Map_type, typename Hamiltonian_type, typename Space_charge_type>
31{
32 public:
34 typedef Value_type value_type;
36 typedef Size_type size_type;
38 typedef Series_type series_type;
40 typedef Map_type map_type;
42 typedef Hamiltonian_type hamiltonian_type;
44 typedef Space_charge_type space_charge_type;
45
47 typedef boost::numeric::ublas::matrix<value_type> matrix_type;
49 typedef std::vector<value_type> vector_type;
50
52
56
58
64
66
72 void combine(std::vector<matrix_type>&, std::vector<matrix_type>&);
73
78 matrix_type combine(std::vector<matrix_type>& maps);
79
82
88 std::pair<value_type, value_type> computeTunes(const matrix_type& map);
89
90 private:
95// /// Stores the one turn transfer map
96// map_type M_m;
97};
98
99// -----------------------------------------------------------------------------------------------------------------------
100// PUBLIC MEMBER FUNCTIONS
101// -----------------------------------------------------------------------------------------------------------------------
102
103
104template<typename Value_type,
105 typename Size_type,
106 typename Series_type,
107 typename Map_type,
108 typename Hamiltonian_type,
109 typename Space_charge_type>
110MapGenerator<Value_type,
111 Size_type,
112 Series_type,
113 Map_type,
114 Hamiltonian_type,
115 Space_charge_type>::MapGenerator(size_type nMaps)
116: nMaps_m(nMaps)
117{}
118
119template<typename Value_type,
120 typename Size_type,
121 typename Series_type,
122 typename Map_type,
123 typename Hamiltonian_type,
124 typename Space_charge_type>
125typename MapGenerator<Value_type,
126 Size_type,
127 Series_type,
128 Map_type,
129 Hamiltonian_type,
130 Space_charge_type>::matrix_type
131MapGenerator<Value_type,
132 Size_type,
133 Series_type,
134 Map_type,
135 Hamiltonian_type,
136 Space_charge_type>::generateMap/*cyc*/(const series_type& H,
137 value_type ds,
138 size_type order) {
139
140 // expand
141 map_type M = ExpMap(-H * ds, order);
142
143 // get linear part
144 FMatrix<value_type,6,6> matrix = M.linearTerms();
145
146 matrix_type map = boost::numeric::ublas::zero_matrix<value_type>(6);
147
148 for (size_type i = 0; i < 6; ++i) {
149 for (size_type j = 0; j < 6; ++j) {
150 map(i,j) = matrix[i][j];
151 }
152 }
153
154 return map;
155}
156
157template<typename Value_type,
158 typename Size_type,
159 typename Series_type,
160 typename Map_type,
161 typename Hamiltonian_type,
162 typename Space_charge_type>
163void MapGenerator<Value_type,
164 Size_type,
165 Series_type,
166 Map_type,
167 Hamiltonian_type,
168 Space_charge_type>::combine(std::vector<matrix_type>& Mscs,
169 std::vector<matrix_type>& Mcycs) {
170
171 if (nMaps_m != Mscs.size() || nMaps_m != Mcycs.size())
172 throw OpalException("MapGenerator::combine()", "Wrong vector dimensions.");
173
174 Mturn_m = boost::numeric::ublas::identity_matrix<value_type>(6);
175
176 for (size_type i = 0; i < nMaps_m; ++i)
177 Mturn_m = matt_boost::gemmm<matrix_type>(Mscs[i],Mcycs[i],Mturn_m);
178}
179
180template<typename Value_type,
181 typename Size_type,
182 typename Series_type,
183 typename Map_type,
184 typename Hamiltonian_type,
185 typename Space_charge_type>
186typename MapGenerator<Value_type,
187 Size_type,
188 Series_type,
189 Map_type,
190 Hamiltonian_type,
191 Space_charge_type>::matrix_type
192MapGenerator<Value_type,
193 Size_type,
194 Series_type,
195 Map_type,
196 Hamiltonian_type,
197 Space_charge_type>::combine(std::vector<matrix_type>& maps)
198{
199 matrix_type map = boost::numeric::ublas::identity_matrix<value_type>(6);
200
201 for (std::size_t i = 0; i < maps.size(); ++i) {
202 matrix_type tmp = prod(maps[i], map);
203 map = tmp;
204 }
205
206 return map;
207}
208
209
210
211template<typename Value_type,
212 typename Size_type,
213 typename Series_type,
214 typename Map_type,
215 typename Hamiltonian_type,
216 typename Space_charge_type>
217typename MapGenerator<Value_type,
218 Size_type,
219 Series_type,
220 Map_type,
221 Hamiltonian_type,
222 Space_charge_type>::matrix_type
223MapGenerator<Value_type,
224 Size_type,
225 Series_type,
226 Map_type,
227 Hamiltonian_type,
228 Space_charge_type>::getMap() {
229 return Mturn_m;
230}
231
232
233template<typename Value_type,
234 typename Size_type,
235 typename Series_type,
236 typename Map_type,
237 typename Hamiltonian_type,
238 typename Space_charge_type>
239
240std::pair<Value_type,
241 Value_type> MapGenerator<Value_type,
242 Size_type,
243 Series_type,
244 Map_type,
245 Hamiltonian_type,
246 Space_charge_type>::computeTunes(const matrix_type& map)
247{
248 /*
249 * M = [ cos(mu) + alpha * sin(mu) beta * sin(mu)
250 * -gamma * sin(mu) cos(mu) - alpha * sin(mu)]
251 *
252 * i = 0, 2, 4
253 * --> cos(mu) = 0.5 * [ M(i, i) + M(i + 1, i + 1) ]
254 */
255
256 value_type arg = 0.0;
257
258 // horizontal phase advance [rad]
259 arg = 0.5 * ( map(0, 0) + map(1, 1) );
260
261 if ( std::abs( arg ) > 1 )
262 throw OpalException("MapGenerator::computeTunes()",
263 "Horizontal phase advance: Acos argument " + std::to_string(arg) + " out of range.");
264
265 value_type mux = std::acos( arg );
266
267 // vertical phase advance [rad]
268 arg = 0.5 * ( map(2, 2) + map(3, 3) );
269
270 if ( std::abs( arg ) > 1 )
271 throw OpalException("MapGenerator::computeTunes()",
272 "Vertical phase advance: Acos argument " + std::to_string(arg) + " out of range.");
273
274 value_type muz = std::acos( arg );
275
276
277
278 return std::make_pair<value_type,
280 muz * Physics::u_two_pi );
281}
282
283#endif
FVps< T, N > ExpMap(const FTps< T, N > &H, int trunc=FTps< T, N >::EXACT)
Build the exponential series.
Definition: FTps.hpp:1994
T * value_type(const SliceIterator< T > &)
T::PETE_Expr_t::PETE_Return_t prod(const PETE_Expr< T > &expr)
Definition: PETE.h:1121
PETE_TUTree< FnArcCos, typename T::PETE_Expr_t > acos(const PETE_Expr< T > &l)
Definition: PETE.h:725
PETE_TUTree< FnAbs, typename T::PETE_Expr_t > abs(const PETE_Expr< T > &l)
arg(a))
constexpr double u_two_pi
The value of.
Definition: Physics.h:36
A templated representation for matrices.
Definition: FMatrix.h:39
This class generates the matrices for the one turn matrix of a cyclotron.
Definition: MapGenerator.h:31
Series_type series_type
Type of truncated power series.
Definition: MapGenerator.h:38
boost::numeric::ublas::matrix< value_type > matrix_type
Type for specifying matrices.
Definition: MapGenerator.h:47
Hamiltonian_type hamiltonian_type
Type of the Hamiltonian.
Definition: MapGenerator.h:42
Value_type value_type
Type of variables.
Definition: MapGenerator.h:34
MapGenerator(size_type)
Initialize.
Definition: MapGenerator.h:115
std::pair< value_type, value_type > computeTunes(const matrix_type &map)
Definition: MapGenerator.h:246
matrix_type Mturn_m
One-turn matrix.
Definition: MapGenerator.h:94
Size_type size_type
Type for specifying sizes.
Definition: MapGenerator.h:36
matrix_type generateMap(const series_type &, value_type, size_type)
Generates a map based on the Hamiltonian for a given angle.
Definition: MapGenerator.h:136
matrix_type getMap()
Returns the one turn map.
Definition: MapGenerator.h:228
Space_charge_type space_charge_type
Type of the Hamiltonian representing the space charge.
Definition: MapGenerator.h:44
std::vector< value_type > vector_type
Type for specifying vectors.
Definition: MapGenerator.h:49
void combine(std::vector< matrix_type > &, std::vector< matrix_type > &)
Combine given maps.
Definition: MapGenerator.h:168
size_type nMaps_m
Number of maps.
Definition: MapGenerator.h:92
Map_type map_type
Type of a map.
Definition: MapGenerator.h:40
The base class for all OPAL exceptions.
Definition: OpalException.h:28