OPAL (Object Oriented Parallel Accelerator Library)  2.2.0
OPAL
Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
SigmaGenerator< Value_type, Size_type > Class Template Reference

This class computes the matched distribution. More...

#include <SigmaGenerator.h>

Collaboration diagram for SigmaGenerator< Value_type, Size_type >:
Collaboration graph
[legend]

Public Types

typedef Value_type value_type
 Type of variables. More...
 
typedef const value_type const_value_type
 Type of constant variables. More...
 
typedef Size_type size_type
 Type for specifying sizes. More...
 
typedef
boost::numeric::ublas::matrix
< value_type
matrix_type
 Type for storing maps. More...
 
typedef std::complex< value_typecomplex_t
 
typedef
boost::numeric::ublas::matrix
< complex_t
cmatrix_type
 Type for storing complex matrices. More...
 
typedef
boost::numeric::ublas::compressed_matrix
< complex_t,
boost::numeric::ublas::row_major > 
sparse_matrix_type
 Type for storing the sparse maps. More...
 
typedef
boost::numeric::ublas::vector
< value_type
vector_type
 Type for storing vectors. More...
 
typedef std::vector< value_typecontainer_type
 Container for storing the properties for each angle. More...
 
typedef FTps< value_type, 2 *3 > Series
 Type of the truncated powere series. More...
 
typedef FVps< value_type, 2 *3 > Map
 Type of a map. More...
 
typedef std::function< Series(value_type,
value_type, value_type)> 
Hamiltonian
 Type of the Hamiltonian for the cyclotron. More...
 
typedef std::function< Series(value_type,
value_type, value_type)> 
SpaceCharge
 Type of the Hamiltonian for the space charge. More...
 

Public Member Functions

 SigmaGenerator (value_type I, value_type ex, value_type ey, value_type ez, value_type E, value_type m, const Cyclotron *cycl, size_type N, size_type Nsectors, size_type truncOrder, bool write=true)
 Constructs an object of type SigmaGenerator. More...
 
bool match (value_type accuracy, size_type maxit, size_type maxitOrbit, Cyclotron *cycl, value_type denergy, value_type rguess, bool harmonic, bool full)
 Searches for a matched distribution. More...
 
void eigsolve_m (const matrix_type &Mturn, sparse_matrix_type &R)
 
bool invertMatrix_m (const sparse_matrix_type &R, sparse_matrix_type &invR)
 
void decouple (const matrix_type &Mturn, sparse_matrix_type &R, sparse_matrix_type &invR)
 Block diagonalizes the symplex part of the one turn transfer matrix. More...
 
value_type isEigenEllipse (const matrix_type &Mturn, const matrix_type &sigma)
 Checks if the sigma-matrix is an eigenellipse and returns the L2 error. More...
 
matrix_typegetSigma ()
 Returns the converged stationary distribution. More...
 
size_type getIterations () const
 Returns the number of iterations needed for convergence (if not converged, it returns zero) More...
 
value_type getError () const
 Returns the error (if the program didn't converged, one can call this function to check the error) More...
 
std::array< value_type, 3 > getEmittances () const
 Returns the emittances (ex,ey,ez) in \( \pi\ mm\ mrad \) for which the code converged (since the whole simulation is done on normalized emittances) More...
 
const double & getInjectionRadius () const
 
const double & getInjectionMomentum () const
 

Private Member Functions

void initialize (value_type, value_type)
 
matrix_type updateInitialSigma (const matrix_type &, sparse_matrix_type &, sparse_matrix_type &)
 Computes the new initial sigma matrix. More...
 
void updateSigma (const std::vector< matrix_type > &, const std::vector< matrix_type > &)
 Computes new sigma matrices (one for each angle) More...
 
value_type L2ErrorNorm (const matrix_type &, const matrix_type &)
 Returns the L2-error norm between the old and new sigma-matrix. More...
 
value_type L1ErrorNorm (const matrix_type &, const matrix_type &)
 Returns the Lp-error norm between the old and new sigma-matrix. More...
 
std::string float2string (value_type val)
 Transforms a floating point value to a string. More...
 
void writeOrbitOutput_m (const std::pair< value_type, value_type > &tunes, const value_type &ravg, const value_type &freqError, const container_type &r_turn, const container_type &peo, const container_type &h_turn, const container_type &fidx_turn, const container_type &ds_turn)
 Called within SigmaGenerator::match(). More...
 

Private Attributes

value_type I_m
 Stores the value of the current, \( \left[I\right] = A \). More...
 
std::array< value_type, 3 > emittance_m
 Stores the desired emittances, \( \left[\varepsilon_{x}\right] = \left[\varepsilon_{y}\right] = \left[\varepsilon_{z}\right] = mm \ mrad \). More...
 
value_type wo_m
 Is the orbital frequency, \( \left[\omega_{o}\right] = \frac{1}{s} \). More...
 
value_type E_m
 Stores the user-define energy, \( \left[E\right] = MeV \). More...
 
value_type gamma_m
 Relativistic factor (which can be computed out ot the kinetic energy and rest mass (potential energy)), \( \left[\gamma\right] = 1 \). More...
 
value_type gamma2_m
 Relativistic factor squared. More...
 
value_type nh_m
 Harmonic number, \( \left[N_{h}\right] = 1 \). More...
 
value_type beta_m
 Velocity (c/v), \( \left[\beta\right] = 1 \). More...
 
value_type m_m
 Is the mass of the particles, \( \left[m\right] = \frac{MeV}{c^{2}} \). More...
 
size_type niterations_m
 Is the number of iterations needed for convergence. More...
 
bool converged_m
 Is true if converged, false otherwise. More...
 
value_type Emin_m
 Minimum energy needed in cyclotron, \( \left[E_{min}\right] = MeV \). More...
 
value_type Emax_m
 Maximum energy reached in cyclotron, \( \left[E_{max}\right] = MeV \). More...
 
size_type N_m
 Number of integration steps for closed orbit computation. More...
 
size_type nSectors_m
 Number of (symmetric) sectors the field is averaged over. More...
 
size_type nStepsPerSector_m
 Number of integration steps per sector (–> also: number of maps per sector) More...
 
size_type nSteps_m
 Number of integration steps in total. More...
 
value_type error_m
 Error of computation. More...
 
size_type truncOrder_m
 Truncation order of the power series for the Hamiltonian (cyclotron and space charge) More...
 
bool write_m
 Decides for writing output (default: true) More...
 
matrix_type sigma_m
 Stores the stationary distribution (sigma matrix) More...
 
std::vector< matrix_typesigmas_m
 Vector storing the second moment matrix for each angle. More...
 
Hamiltonian H_m
 Stores the Hamiltonian of the cyclotron. More...
 
SpaceCharge Hsc_m
 Stores the Hamiltonian for the space charge. More...
 
Series x_m
 All variables x, px, y, py, z, delta. More...
 
Series px_m
 
Series y_m
 
Series py_m
 
Series z_m
 
Series delta_m
 
double rinit_m
 
double prinit_m
 

Detailed Description

template<typename Value_type, typename Size_type>
class SigmaGenerator< Value_type, Size_type >

This class computes the matched distribution.

Definition at line 64 of file SigmaGenerator.h.

Member Typedef Documentation

template<typename Value_type, typename Size_type>
typedef boost::numeric::ublas::matrix<complex_t> SigmaGenerator< Value_type, Size_type >::cmatrix_type

Type for storing complex matrices.

Definition at line 79 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef std::complex<value_type> SigmaGenerator< Value_type, Size_type >::complex_t

Definition at line 76 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef const value_type SigmaGenerator< Value_type, Size_type >::const_value_type

Type of constant variables.

Definition at line 70 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef std::vector<value_type> SigmaGenerator< Value_type, Size_type >::container_type

Container for storing the properties for each angle.

Definition at line 88 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef std::function<Series(value_type,value_type,value_type)> SigmaGenerator< Value_type, Size_type >::Hamiltonian

Type of the Hamiltonian for the cyclotron.

Definition at line 94 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef FVps<value_type,2*3> SigmaGenerator< Value_type, Size_type >::Map

Type of a map.

Definition at line 92 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef boost::numeric::ublas::matrix<value_type> SigmaGenerator< Value_type, Size_type >::matrix_type

Type for storing maps.

Definition at line 74 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef FTps<value_type,2*3> SigmaGenerator< Value_type, Size_type >::Series

Type of the truncated powere series.

Definition at line 90 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef Size_type SigmaGenerator< Value_type, Size_type >::size_type

Type for specifying sizes.

Definition at line 72 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef std::function<Series(value_type,value_type,value_type)> SigmaGenerator< Value_type, Size_type >::SpaceCharge

Type of the Hamiltonian for the space charge.

Definition at line 96 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef boost::numeric::ublas::compressed_matrix<complex_t, boost::numeric::ublas::row_major > SigmaGenerator< Value_type, Size_type >::sparse_matrix_type

Type for storing the sparse maps.

Definition at line 84 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef Value_type SigmaGenerator< Value_type, Size_type >::value_type

Type of variables.

Definition at line 68 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
typedef boost::numeric::ublas::vector<value_type> SigmaGenerator< Value_type, Size_type >::vector_type

Type for storing vectors.

Definition at line 86 of file SigmaGenerator.h.

Constructor & Destructor Documentation

template<typename Value_type , typename Size_type >
SigmaGenerator< Value_type, Size_type >::SigmaGenerator ( value_type  I,
value_type  ex,
value_type  ey,
value_type  ez,
value_type  E,
value_type  m,
const Cyclotron cycl,
size_type  N,
size_type  Nsectors,
size_type  truncOrder,
bool  write = true 
)

Constructs an object of type SigmaGenerator.

Parameters
Ispecifies the current for which a matched distribution should be found, \( [I] = A \)
exis the emittance in x-direction (horizontal), \( \left[\varepsilon_{x}\right] = \pi\ mm\ mrad \)
eyis the emittance in y-direction (longitudinal), \( \left[\varepsilon_{y}\right] = \pi\ mm\ mrad \)
ezis the emittance in z-direction (vertical), \( \left[\varepsilon_{z}\right] = \pi\ mm\ mrad \)
Eis the energy, \( \left[E\right] = MeV \)
mis the mass of the particles \( \left[m\right] = \frac{MeV}{c^{2}} \)
cyclis the cyclotron element
Nis the number of integration steps (closed orbit computation). That's why its also the number of maps (for each integration step a map)
Nsectorsis the number of sectors that the field map is averaged over (closed orbit computation)
truncOrderis the truncation order for power series of the Hamiltonian
writeis a boolean (default: true). If true all maps of all iterations are stored, otherwise not.

Definition at line 323 of file SigmaGenerator.h.

References SigmaGenerator< Value_type, Size_type >::beta_m, Physics::c, SigmaGenerator< Value_type, Size_type >::delta_m, SigmaGenerator< Value_type, Size_type >::Emin_m, SigmaGenerator< Value_type, Size_type >::emittance_m, Physics::epsilon_0, fabs(), SigmaGenerator< Value_type, Size_type >::gamma2_m, SigmaGenerator< Value_type, Size_type >::gamma_m, SigmaGenerator< Value_type, Size_type >::H_m, SigmaGenerator< Value_type, Size_type >::Hsc_m, SigmaGenerator< Value_type, Size_type >::I_m, SigmaGenerator< Value_type, Size_type >::m_m, FTps< value_type, 2 *3 >::makeVariable(), SigmaGenerator< Value_type, Size_type >::nh_m, Physics::pi, SigmaGenerator< Value_type, Size_type >::px_m, SigmaGenerator< Value_type, Size_type >::py_m, FTps< value_type, 2 *3 >::setGlobalTruncOrder(), sqrt(), SigmaGenerator< Value_type, Size_type >::truncOrder_m, SigmaGenerator< Value_type, Size_type >::wo_m, SigmaGenerator< Value_type, Size_type >::x_m, SigmaGenerator< Value_type, Size_type >::y_m, and SigmaGenerator< Value_type, Size_type >::z_m.

Here is the call graph for this function:

Member Function Documentation

template<typename Value_type , typename Size_type >
void SigmaGenerator< Value_type, Size_type >::decouple ( const matrix_type Mturn,
sparse_matrix_type R,
sparse_matrix_type invR 
)

Block diagonalizes the symplex part of the one turn transfer matrix.

It computes the transformation matrix R and its inverse invR.

Parameters
Mturnis a 6x6 dimensional one turn transfer matrix
Ris the 6x6 dimensional transformation matrix (gets computed)
invRis the 6x6 dimensional inverse transformation (gets computed)

Definition at line 816 of file SigmaGenerator.h.

template<typename Value_type , typename Size_type >
void SigmaGenerator< Value_type, Size_type >::eigsolve_m ( const matrix_type Mturn,
sparse_matrix_type R 
)

Eigenvalue / eigenvector solver

Parameters
Mturnis a 6x6 dimensional one turn transfer matrix
Ris the 6x6 dimensional transformation matrix (gets computed)

Definition at line 681 of file SigmaGenerator.h.

References abs(), Physics::e, Attrib::Distribution::R, and sqrt().

Here is the call graph for this function:

template<typename Value_type , typename Size_type >
std::string SigmaGenerator< Value_type, Size_type >::float2string ( value_type  val)
private

Transforms a floating point value to a string.

Parameters
valis the floating point value which is transformed to a string

Definition at line 1175 of file SigmaGenerator.h.

template<typename Value_type , typename Size_type >
std::array< Value_type, 3 > SigmaGenerator< Value_type, Size_type >::getEmittances ( ) const
inline

Returns the emittances (ex,ey,ez) in \( \pi\ mm\ mrad \) for which the code converged (since the whole simulation is done on normalized emittances)

Definition at line 864 of file SigmaGenerator.h.

References Physics::pi.

template<typename Value_type , typename Size_type >
SigmaGenerator< Value_type, Size_type >::value_type SigmaGenerator< Value_type, Size_type >::getError ( ) const
inline

Returns the error (if the program didn't converged, one can call this function to check the error)

Definition at line 857 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
const double& SigmaGenerator< Value_type, Size_type >::getInjectionMomentum ( ) const
inline
template<typename Value_type, typename Size_type>
const double& SigmaGenerator< Value_type, Size_type >::getInjectionRadius ( ) const
inline
template<typename Value_type , typename Size_type >
SigmaGenerator< Value_type, Size_type >::size_type SigmaGenerator< Value_type, Size_type >::getIterations ( ) const
inline

Returns the number of iterations needed for convergence (if not converged, it returns zero)

Definition at line 850 of file SigmaGenerator.h.

template<typename Value_type , typename Size_type >
SigmaGenerator< Value_type, Size_type >::matrix_type & SigmaGenerator< Value_type, Size_type >::getSigma ( )
inline

Returns the converged stationary distribution.

Definition at line 843 of file SigmaGenerator.h.

template<typename Value_type , typename Size_type >
void SigmaGenerator< Value_type, Size_type >::initialize ( value_type  nuz,
value_type  ravg 
)
private

Initializes a first guess of the sigma matrix with the assumption of a spherical symmetric beam (ex = ey = ez). For each angle split the same initial guess is taken.

Parameters
nuzis the vertical tune
ravgis the average radius of the closed orbit

Definition at line 879 of file SigmaGenerator.h.

References Physics::alpha, Physics::c, Physics::e, endl(), Physics::epsilon_0, K, Physics::mu_0, Physics::pi, and sqrt().

Here is the call graph for this function:

template<typename Value_type , typename Size_type >
bool SigmaGenerator< Value_type, Size_type >::invertMatrix_m ( const sparse_matrix_type R,
sparse_matrix_type invR 
)
Parameters
Ris the 6x6 dimensional transformation matrix
invRis the 6x6 dimensional inverse transformation (gets computed)
Returns
true if success

Definition at line 788 of file SigmaGenerator.h.

template<typename Value_type , typename Size_type >
SigmaGenerator< Value_type, Size_type >::value_type SigmaGenerator< Value_type, Size_type >::isEigenEllipse ( const matrix_type Mturn,
const matrix_type sigma 
)

Checks if the sigma-matrix is an eigenellipse and returns the L2 error.

The idea of this function is taken from Dr. Christian Baumgarten's program.

Parameters
Mturnis the one turn transfer matrix
sigmais the sigma matrix to be tested

Definition at line 829 of file SigmaGenerator.h.

template<typename Value_type , typename Size_type >
SigmaGenerator< Value_type, Size_type >::value_type SigmaGenerator< Value_type, Size_type >::L1ErrorNorm ( const matrix_type oldS,
const matrix_type newS 
)
private

Returns the Lp-error norm between the old and new sigma-matrix.

Parameters
oldSis the old sigma matrix
newSis the new sigma matrix

Definition at line 1158 of file SigmaGenerator.h.

References abs(), and for_each().

Here is the call graph for this function:

template<typename Value_type , typename Size_type >
SigmaGenerator< Value_type, Size_type >::value_type SigmaGenerator< Value_type, Size_type >::L2ErrorNorm ( const matrix_type oldS,
const matrix_type newS 
)
private

Returns the L2-error norm between the old and new sigma-matrix.

Parameters
oldSis the old sigma matrix
newSis the new sigma matrix

Definition at line 1143 of file SigmaGenerator.h.

References sqrt().

Here is the call graph for this function:

template<typename Value_type , typename Size_type >
bool SigmaGenerator< Value_type, Size_type >::match ( value_type  accuracy,
size_type  maxit,
size_type  maxitOrbit,
Cyclotron cycl,
value_type  denergy,
value_type  rguess,
bool  harmonic,
bool  full 
)

Searches for a matched distribution.

Returns "true" if a matched distribution within the accuracy could be found, returns "false" otherwise.

Parameters
accuracyis used for computing the equilibrium orbit (to a certain accuracy)
maxitis the maximum number of iterations (the program stops if within this value no stationary distribution was found)
maxitOrbitis the maximum number of iterations for finding closed orbit
angledefines the start of the sector (one can choose any angle between 0° and 360°)
denergyenergy step size for closed orbit finder [MeV]
rguessvalue of radius for closed orbit finder
typespecifies the magnetic field format (e.g. CARBONCYCL)
harmonicis a boolean. If "true" the harmonics are used instead of the closed orbit finder.
fullmatch over full turn not just single sector

Definition at line 424 of file SigmaGenerator.h.

References Physics::e, endl(), Cyclotron::getPHIinit(), Attrib::Distribution::R, and value_type().

Here is the call graph for this function:

template<typename Value_type , typename Size_type >
SigmaGenerator< Value_type, Size_type >::matrix_type SigmaGenerator< Value_type, Size_type >::updateInitialSigma ( const matrix_type M,
sparse_matrix_type R,
sparse_matrix_type invR 
)
private

Computes the new initial sigma matrix.

Parameters
Mis the 6x6 one turn transfer matrix
Ris the transformation matrix
invRis the inverse transformation matrix

Definition at line 1041 of file SigmaGenerator.h.

References endl(), interpolation::im(), and prod().

Here is the call graph for this function:

template<typename Value_type , typename Size_type >
void SigmaGenerator< Value_type, Size_type >::updateSigma ( const std::vector< matrix_type > &  Mscs,
const std::vector< matrix_type > &  Mcycs 
)
private

Computes new sigma matrices (one for each angle)

Mscs is a vector of all space charge maps Mcycs is a vector of all cyclotron maps

Definition at line 1111 of file SigmaGenerator.h.

References endl(), and prod().

Here is the call graph for this function:

template<typename Value_type , typename Size_type >
void SigmaGenerator< Value_type, Size_type >::writeOrbitOutput_m ( const std::pair< value_type, value_type > &  tunes,
const value_type ravg,
const value_type freqError,
const container_type r_turn,
const container_type peo,
const container_type h_turn,
const container_type fidx_turn,
const container_type ds_turn 
)
private

Called within SigmaGenerator::match().

Parameters
tunes
ravgis the average radius [m]
r_turnis the radius [m]
peois the momentum
h_turnis the inverse bending radius
fidx_turnis the field index
ds_turnis the path length element

Definition at line 1183 of file SigmaGenerator.h.

References endl().

Here is the call graph for this function:

Member Data Documentation

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::beta_m
private

Velocity (c/v), \( \left[\beta\right] = 1 \).

Definition at line 202 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
bool SigmaGenerator< Value_type, Size_type >::converged_m
private

Is true if converged, false otherwise.

Definition at line 208 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
Series SigmaGenerator< Value_type, Size_type >::delta_m
private
template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::E_m
private

Stores the user-define energy, \( \left[E\right] = MeV \).

Definition at line 194 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::Emax_m
private

Maximum energy reached in cyclotron, \( \left[E_{max}\right] = MeV \).

Definition at line 212 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::Emin_m
private

Minimum energy needed in cyclotron, \( \left[E_{min}\right] = MeV \).

Definition at line 210 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
std::array<value_type,3> SigmaGenerator< Value_type, Size_type >::emittance_m
private

Stores the desired emittances, \( \left[\varepsilon_{x}\right] = \left[\varepsilon_{y}\right] = \left[\varepsilon_{z}\right] = mm \ mrad \).

Definition at line 190 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::error_m
private

Error of computation.

Definition at line 224 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::gamma2_m
private

Relativistic factor squared.

Definition at line 198 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::gamma_m
private

Relativistic factor (which can be computed out ot the kinetic energy and rest mass (potential energy)), \( \left[\gamma\right] = 1 \).

Definition at line 196 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
Hamiltonian SigmaGenerator< Value_type, Size_type >::H_m
private

Stores the Hamiltonian of the cyclotron.

Definition at line 239 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
SpaceCharge SigmaGenerator< Value_type, Size_type >::Hsc_m
private

Stores the Hamiltonian for the space charge.

Definition at line 242 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::I_m
private

Stores the value of the current, \( \left[I\right] = A \).

Definition at line 188 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::m_m
private

Is the mass of the particles, \( \left[m\right] = \frac{MeV}{c^{2}} \).

Definition at line 204 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
size_type SigmaGenerator< Value_type, Size_type >::N_m
private

Number of integration steps for closed orbit computation.

Definition at line 214 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::nh_m
private

Harmonic number, \( \left[N_{h}\right] = 1 \).

Definition at line 200 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
size_type SigmaGenerator< Value_type, Size_type >::niterations_m
private

Is the number of iterations needed for convergence.

Definition at line 206 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
size_type SigmaGenerator< Value_type, Size_type >::nSectors_m
private

Number of (symmetric) sectors the field is averaged over.

Definition at line 216 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
size_type SigmaGenerator< Value_type, Size_type >::nSteps_m
private

Number of integration steps in total.

Definition at line 221 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
size_type SigmaGenerator< Value_type, Size_type >::nStepsPerSector_m
private

Number of integration steps per sector (–> also: number of maps per sector)

Definition at line 218 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
double SigmaGenerator< Value_type, Size_type >::prinit_m
private
template<typename Value_type, typename Size_type>
Series SigmaGenerator< Value_type, Size_type >::px_m
private
template<typename Value_type, typename Size_type>
Series SigmaGenerator< Value_type, Size_type >::py_m
private
template<typename Value_type, typename Size_type>
double SigmaGenerator< Value_type, Size_type >::rinit_m
private
template<typename Value_type, typename Size_type>
matrix_type SigmaGenerator< Value_type, Size_type >::sigma_m
private

Stores the stationary distribution (sigma matrix)

Definition at line 233 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
std::vector<matrix_type> SigmaGenerator< Value_type, Size_type >::sigmas_m
private

Vector storing the second moment matrix for each angle.

Definition at line 236 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
size_type SigmaGenerator< Value_type, Size_type >::truncOrder_m
private

Truncation order of the power series for the Hamiltonian (cyclotron and space charge)

Definition at line 227 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
value_type SigmaGenerator< Value_type, Size_type >::wo_m
private

Is the orbital frequency, \( \left[\omega_{o}\right] = \frac{1}{s} \).

Definition at line 192 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
bool SigmaGenerator< Value_type, Size_type >::write_m
private

Decides for writing output (default: true)

Definition at line 230 of file SigmaGenerator.h.

template<typename Value_type, typename Size_type>
Series SigmaGenerator< Value_type, Size_type >::x_m
private

All variables x, px, y, py, z, delta.

Definition at line 245 of file SigmaGenerator.h.

Referenced by SigmaGenerator< Value_type, Size_type >::SigmaGenerator().

template<typename Value_type, typename Size_type>
Series SigmaGenerator< Value_type, Size_type >::y_m
private
template<typename Value_type, typename Size_type>
Series SigmaGenerator< Value_type, Size_type >::z_m
private

The documentation for this class was generated from the following file: