OPAL (Object Oriented Parallel Accelerator Library)  2021.1.99
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 
19 #include "matrix_vector_operation.h"
20 
21 #include "FixedAlgebra/FTps.h"
22 #include "FixedAlgebra/FMatrix.h"
23 
24 #include "Physics/Physics.h"
25 #include "Utilities/OpalException.h"
26 
28 template<typename Value_type, typename Size_type, typename Series_type,
29  typename Map_type, typename Hamiltonian_type, typename Space_charge_type>
30 class MapGenerator
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 
55  MapGenerator(size_type);
56 
58 
63  matrix_type generateMap(const series_type&, value_type, size_type);
64 
66 
72  void combine(std::vector<matrix_type>&, std::vector<matrix_type>&);
73 
78  matrix_type combine(std::vector<matrix_type>& maps);
79 
81  matrix_type getMap();
82 
88  std::pair<value_type, value_type> computeTunes(const matrix_type& map);
89 
90  private:
92  size_type nMaps_m;
94  matrix_type Mturn_m;
95 // /// Stores the one turn transfer map
96 // map_type M_m;
97 };
98 
99 // -----------------------------------------------------------------------------------------------------------------------
100 // PUBLIC MEMBER FUNCTIONS
101 // -----------------------------------------------------------------------------------------------------------------------
102 
103 
104 template<typename Value_type,
105  typename Size_type,
106  typename Series_type,
107  typename Map_type,
108  typename Hamiltonian_type,
109  typename Space_charge_type>
110 MapGenerator<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 
119 template<typename Value_type,
120  typename Size_type,
121  typename Series_type,
122  typename Map_type,
123  typename Hamiltonian_type,
124  typename Space_charge_type>
125 typename MapGenerator<Value_type,
126  Size_type,
127  Series_type,
128  Map_type,
129  Hamiltonian_type,
130  Space_charge_type>::matrix_type
131 MapGenerator<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 
157 template<typename Value_type,
158  typename Size_type,
159  typename Series_type,
160  typename Map_type,
161  typename Hamiltonian_type,
162  typename Space_charge_type>
163 void 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 
180 template<typename Value_type,
181  typename Size_type,
182  typename Series_type,
183  typename Map_type,
184  typename Hamiltonian_type,
185  typename Space_charge_type>
186 typename MapGenerator<Value_type,
187  Size_type,
188  Series_type,
189  Map_type,
190  Hamiltonian_type,
191  Space_charge_type>::matrix_type
192 MapGenerator<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 
211 template<typename Value_type,
212  typename Size_type,
213  typename Series_type,
214  typename Map_type,
215  typename Hamiltonian_type,
216  typename Space_charge_type>
217 typename MapGenerator<Value_type,
218  Size_type,
219  Series_type,
220  Map_type,
221  Hamiltonian_type,
222  Space_charge_type>::matrix_type
223 MapGenerator<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 
233 template<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 
240 std::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,
279  value_type>( mux * Physics::u_two_pi,
280  muz * Physics::u_two_pi );
281 }
282 
283 #endif
value_type
T * value_type(const SliceIterator< T > &)
Definition: SliceIterator.h:191
ExpMap
FVps< T, N > ExpMap(const FTps< T, N > &H, int trunc=FTps< T, N >::EXACT)
Build the exponential series.
Definition: FTps.hpp:1994
abs
PETE_TUTree< FnAbs, typename T::PETE_Expr_t > abs(const PETE_Expr< T > &l)
Definition: IpplExpressions.h:62
arg
arg(a))
prod
T::PETE_Expr_t::PETE_Return_t prod(const PETE_Expr< T > &expr)
Definition: PETE.h:1121
acos
PETE_TUTree< FnArcCos, typename T::PETE_Expr_t > acos(const PETE_Expr< T > &l)
Definition: PETE.h:725
Physics::u_two_pi
constexpr double u_two_pi
The value of.
Definition: Physics.h:36
FMatrix
A templated representation for matrices.
Definition: FMatrix.h:39
MapGenerator
This class generates the matrices for the one turn matrix of a cyclotron.
Definition: MapGenerator.h:31
MapGenerator::series_type
Series_type series_type
Type of truncated power series.
Definition: MapGenerator.h:38
MapGenerator::matrix_type
boost::numeric::ublas::matrix< value_type > matrix_type
Type for specifying matrices.
Definition: MapGenerator.h:47
MapGenerator::hamiltonian_type
Hamiltonian_type hamiltonian_type
Type of the Hamiltonian.
Definition: MapGenerator.h:42
MapGenerator::value_type
Value_type value_type
Type of variables.
Definition: MapGenerator.h:34
MapGenerator::MapGenerator
MapGenerator(size_type)
Initialize.
Definition: MapGenerator.h:115
MapGenerator::computeTunes
std::pair< value_type, value_type > computeTunes(const matrix_type &map)
Definition: MapGenerator.h:246
MapGenerator::Mturn_m
matrix_type Mturn_m
One-turn matrix.
Definition: MapGenerator.h:94
MapGenerator::size_type
Size_type size_type
Type for specifying sizes.
Definition: MapGenerator.h:36
MapGenerator::generateMap
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
MapGenerator::getMap
matrix_type getMap()
Returns the one turn map.
Definition: MapGenerator.h:228
MapGenerator::space_charge_type
Space_charge_type space_charge_type
Type of the Hamiltonian representing the space charge.
Definition: MapGenerator.h:44
MapGenerator::vector_type
std::vector< value_type > vector_type
Type for specifying vectors.
Definition: MapGenerator.h:49
MapGenerator::combine
void combine(std::vector< matrix_type > &, std::vector< matrix_type > &)
Combine given maps.
Definition: MapGenerator.h:168
MapGenerator::nMaps_m
size_type nMaps_m
Number of maps.
Definition: MapGenerator.h:92
MapGenerator::map_type
Map_type map_type
Type of a map.
Definition: MapGenerator.h:40
OpalException
The base class for all OPAL exceptions.
Definition: OpalException.h:28
Generated on Wed Aug 25 2021 16:41:03 for OPAL (Object Oriented Parallel Accelerator Library) by doxygen 1.9.1