ImpactX
Loading...
Searching...
No Matches
impactx::elements::SoftSolenoid Struct Reference

#include <SoftSol.H>

Inheritance diagram for impactx::elements::SoftSolenoid:
impactx::elements::mixin::Named impactx::elements::mixin::BeamOptic< SoftSolenoid > impactx::elements::mixin::LinearTransport< SoftSolenoid > impactx::elements::mixin::Thick impactx::elements::mixin::Alignment impactx::elements::mixin::PipeAperture amrex::simd::Vectorized< amrex::simd::native_simd_size_particlereal > amrex::simd::detail::InternalVectorized

Public Types

using PType = ImpactXParticleContainer::ParticleType
 

Public Member Functions

 SoftSolenoid (amrex::ParticleReal ds, amrex::ParticleReal bscale, std::vector< amrex::ParticleReal > cos_coef, std::vector< amrex::ParticleReal > sin_coef, int unit, amrex::ParticleReal dx=0, amrex::ParticleReal dy=0, amrex::ParticleReal rotation_degree=0, amrex::ParticleReal aperture_x=0, amrex::ParticleReal aperture_y=0, int mapsteps=1, int nslice=1, std::optional< std::string > name=std::nullopt)
 
void compute_constants (RefPart const &refpart)
 
template<typename T_Real = amrex::ParticleReal, typename T_IdCpu = uint64_t>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void operator() (T_Real &AMREX_RESTRICT x, T_Real &AMREX_RESTRICT y, T_Real &AMREX_RESTRICT t, T_Real &AMREX_RESTRICT px, T_Real &AMREX_RESTRICT py, T_Real &AMREX_RESTRICT pt, T_IdCpu &AMREX_RESTRICT idcpu, RefPart const &AMREX_RESTRICT refpart) const
 
AMREX_GPU_HOST AMREX_FORCE_INLINE void operator() (RefPart &AMREX_RESTRICT refpart) const
 
AMREX_GPU_HOST AMREX_FORCE_INLINE Map6x6 transport_map (RefPart const &AMREX_RESTRICT refpart) const
 
std::tuple< amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal > AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Sol_Bfield (amrex::ParticleReal const zeval) const
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void map1 (amrex::ParticleReal const tau, RefPart &refpart, amrex::ParticleReal &zeval) const
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void map2 (amrex::ParticleReal const tau, RefPart &refpart, amrex::ParticleReal &zeval) const
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void map3 (amrex::ParticleReal const tau, RefPart &refpart, amrex::ParticleReal &zeval) const
 
void finalize ()
 
- Public Member Functions inherited from impactx::elements::mixin::Named
AMREX_GPU_HOST void set_name (std::string const &new_name)
 
AMREX_GPU_HOST Named (std::optional< std::string > name)
 
AMREX_GPU_HOST_DEVICE ~Named ()
 
AMREX_GPU_HOST Named (Named const &other)
 
AMREX_GPU_HOST Namedoperator= (Named const &other)
 
AMREX_GPU_HOST_DEVICE Named (Named &&other) noexcept
 
AMREX_GPU_HOST_DEVICE Namedoperator= (Named &&other) noexcept
 
AMREX_FORCE_INLINE std::string name () const
 
AMREX_FORCE_INLINE bool has_name () const
 
- Public Member Functions inherited from impactx::elements::mixin::BeamOptic< SoftSolenoid >
void operator() (ImpactXParticleContainer &pc, int step, int period)
 
void operator() (ImpactXParticleContainer::iterator &pti, RefPart &AMREX_RESTRICT ref_part)
 
- Public Member Functions inherited from impactx::elements::mixin::LinearTransport< SoftSolenoid >
 LinearTransport ()=default
 
 LinearTransport (LinearTransport const &)=default
 
 LinearTransport (LinearTransport &&)=default
 
LinearTransportoperator= (LinearTransport const &)=default
 
LinearTransportoperator= (LinearTransport &&rhs)=default
 
 ~LinearTransport ()=default
 
AMREX_GPU_HOST AMREX_FORCE_INLINE void operator() (Map6x6 &AMREX_RESTRICT cm, RefPart const &AMREX_RESTRICT ref) const
 
- Public Member Functions inherited from impactx::elements::mixin::Thick
 Thick (amrex::ParticleReal ds, int nslice)
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int nslice () const
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::ParticleReal ds () const
 
- Public Member Functions inherited from impactx::elements::mixin::Alignment
 Alignment (amrex::ParticleReal dx, amrex::ParticleReal dy, amrex::ParticleReal rotation_degree)
 
 Alignment ()=default
 
 Alignment (Alignment const &)=default
 
Alignmentoperator= (Alignment const &)=default
 
 Alignment (Alignment &&)=default
 
Alignmentoperator= (Alignment &&rhs)=default
 
 ~Alignment ()=default
 
void compute_constants (RefPart const &refpart)
 
template<typename T_Real>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void shift_in (T_Real &AMREX_RESTRICT x, T_Real &AMREX_RESTRICT y, T_Real &AMREX_RESTRICT px, T_Real &AMREX_RESTRICT py) const
 
template<typename T_Real>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void shift_out (T_Real &AMREX_RESTRICT x, T_Real &AMREX_RESTRICT y, T_Real &AMREX_RESTRICT px, T_Real &AMREX_RESTRICT py) const
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::ParticleReal dx () const
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::ParticleReal dy () const
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::ParticleReal rotation () const
 
- Public Member Functions inherited from impactx::elements::mixin::PipeAperture
 PipeAperture (amrex::ParticleReal aperture_x, amrex::ParticleReal aperture_y)
 
 PipeAperture ()=default
 
 PipeAperture (PipeAperture const &)=default
 
PipeApertureoperator= (PipeAperture const &)=default
 
 PipeAperture (PipeAperture &&)=default
 
PipeApertureoperator= (PipeAperture &&rhs)=default
 
 ~PipeAperture ()=default
 
template<typename T_Real = amrex::ParticleReal, typename T_IdCpu = uint64_t>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void apply_aperture (T_Real &AMREX_RESTRICT x, T_Real &AMREX_RESTRICT y, T_IdCpu &AMREX_RESTRICT idcpu) const
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::ParticleReal aperture_x () const
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::ParticleReal aperture_y () const
 

Public Attributes

amrex::ParticleReal m_bscale
 
int m_unit
 scaling factor for solenoid Bz field
 
int m_mapsteps
 unit specification for quad strength
 
int m_id
 number of map integration steps per slice
 
int m_ncoef = 0
 unique soft solenoid id used for data lookup map
 
amrex::ParticleReal * m_cos_h_data = nullptr
 number of Fourier coefficients
 
amrex::ParticleReal * m_sin_h_data = nullptr
 non-owning pointer to host cosine coefficients
 
amrex::ParticleReal * m_cos_d_data = nullptr
 non-owning pointer to host sine coefficients
 
amrex::ParticleReal * m_sin_d_data = nullptr
 non-owning pointer to device cosine coefficients
 
- Public Attributes inherited from impactx::elements::mixin::Thick
amrex::ParticleReal m_ds
 
int m_nslice
 segment length in m
 
- Public Attributes inherited from impactx::elements::mixin::Alignment
amrex::ParticleReal m_dx = 0
 
amrex::ParticleReal m_dy = 0
 horizontal translation error [m]
 
amrex::ParticleReal m_rotation = 0
 vertical translation error [m]
 
- Public Attributes inherited from impactx::elements::mixin::PipeAperture
amrex::ParticleReal m_inv_aperture_x = 0
 
amrex::ParticleReal m_inv_aperture_y = 0
 inverse of the horizontal aperture size [1/m]
 

Static Public Attributes

static constexpr auto type = "SoftSolenoid"
 
- Static Public Attributes inherited from impactx::elements::mixin::Alignment
static constexpr amrex::ParticleReal degree2rad = ablastr::constant::math::pi / 180.0
 
- Static Public Attributes inherited from amrex::simd::Vectorized< amrex::simd::native_simd_size_particlereal >
static constexpr int simd_width
 

Member Typedef Documentation

◆ PType

Constructor & Destructor Documentation

◆ SoftSolenoid()

impactx::elements::SoftSolenoid::SoftSolenoid ( amrex::ParticleReal ds,
amrex::ParticleReal bscale,
std::vector< amrex::ParticleReal > cos_coef,
std::vector< amrex::ParticleReal > sin_coef,
int unit,
amrex::ParticleReal dx = 0,
amrex::ParticleReal dy = 0,
amrex::ParticleReal rotation_degree = 0,
amrex::ParticleReal aperture_x = 0,
amrex::ParticleReal aperture_y = 0,
int mapsteps = 1,
int nslice = 1,
std::optional< std::string > name = std::nullopt )
inline

A soft-edge solenoid

Parameters
dsSegment length in m
bscaleScaling factor for on-axis longitudinal magnetic field in 1/m (MAD-X convention) = (Bz in T) / (rigidity in T-m) OR Magnetic field Bz in T (SI units)
cos_coefcosine coefficients in Fourier expansion of on-axis magnetic field Bz
sin_coefsine coefficients in Fourier expansion of on-axis magnetic field Bz
unitUnit specification unit = 0 MADX convention (default) unit = 1 SI units
dxhorizontal translation error in m
dyvertical translation error in m
rotation_degreerotation error in the transverse plane [degrees]
aperture_xhorizontal half-aperture in m
aperture_yvertical half-aperture in m
mapstepsnumber of integration steps per slice used for map and reference particle push in applied fields
nslicenumber of slices used for the application of space charge
namea user defined and not necessarily unique name of the element

Member Function Documentation

◆ compute_constants()

void impactx::elements::SoftSolenoid::compute_constants ( RefPart const & refpart)
inline

Compute and cache the constants for the push.

In particular, used to pre-compute and cache variables that are independent of the individually tracked particle.

Parameters
refpartreference particle

◆ finalize()

void impactx::elements::SoftSolenoid::finalize ( )
inline

Close and deallocate all data and handles.

◆ map1()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void impactx::elements::SoftSolenoid::map1 ( amrex::ParticleReal const tau,
RefPart & refpart,
amrex::ParticleReal & zeval ) const
inline

This pushes the reference particle and the linear map matrix elements for a solenoid through the symplectic map associated with H_3 in the Hamiltonian splitting H = H_1 + H_2 + H_3.

Parameters
tauMap step size in m
[in,out]refpartreference particle
[in,out]zevalLongitudinal on-axis location in m

◆ map2()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void impactx::elements::SoftSolenoid::map2 ( amrex::ParticleReal const tau,
RefPart & refpart,
amrex::ParticleReal & zeval ) const
inline

This pushes the reference particle and the linear map matrix elements for a solenoid through the symplectic map associated with H_2 in the Hamiltonian splitting H = H_1 + H_2 + H_3.

Parameters
tauMap step size in m
[in,out]refpartreference particle
[in,out]zevalLongitudinal on-axis location in m

◆ map3()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void impactx::elements::SoftSolenoid::map3 ( amrex::ParticleReal const tau,
RefPart & refpart,
amrex::ParticleReal & zeval ) const
inline

This pushes the reference particle and the linear map matrix elements for a solenoid through the symplectic map associated with H_1 in the Hamiltonian splitting H = H_1 + H_2 + H_3.

Parameters
tauMap step size in m
[in,out]refpartreference particle
[in,out]zevalLongitudinal on-axis location in m

◆ operator()() [1/2]

AMREX_GPU_HOST AMREX_FORCE_INLINE void impactx::elements::SoftSolenoid::operator() ( RefPart &AMREX_RESTRICT refpart) const
inline

This pushes the reference particle.

Parameters
[in,out]refpartreference particle

◆ operator()() [2/2]

template<typename T_Real = amrex::ParticleReal, typename T_IdCpu = uint64_t>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void impactx::elements::SoftSolenoid::operator() ( T_Real &AMREX_RESTRICT x,
T_Real &AMREX_RESTRICT y,
T_Real &AMREX_RESTRICT t,
T_Real &AMREX_RESTRICT px,
T_Real &AMREX_RESTRICT py,
T_Real &AMREX_RESTRICT pt,
T_IdCpu &AMREX_RESTRICT idcpu,
RefPart const &AMREX_RESTRICT refpart ) const
inline

This is a soft-edge solenoid functor, so that a variable of this type can be used like a soft-edge solenoid function.

The

See also
compute_constants method must be called before pushing particles through this operator.
Parameters
xparticle position in x
yparticle position in y
tparticle position in t
pxparticle momentum in x
pyparticle momentum in y
ptparticle momentum in t
idcpuparticle global index
refpartreference particle

◆ Sol_Bfield()

std::tuple< amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal > AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE impactx::elements::SoftSolenoid::Sol_Bfield ( amrex::ParticleReal const zeval) const
inline

This evaluates the on-axis magnetic field Bz at a fixed location z, together with certain required integrals and derivatives. The field returned is normalized to a peak value of 1.

Parameters
zevalLongitudinal on-axis location in m

◆ transport_map()

AMREX_GPU_HOST AMREX_FORCE_INLINE Map6x6 impactx::elements::SoftSolenoid::transport_map ( RefPart const &AMREX_RESTRICT refpart) const
inline

This function returns the linear transport map.

Parameters
[in]refpartreference particle
Returns
6x6 transport matrix

Member Data Documentation

◆ m_bscale

amrex::ParticleReal impactx::elements::SoftSolenoid::m_bscale

◆ m_cos_d_data

amrex::ParticleReal* impactx::elements::SoftSolenoid::m_cos_d_data = nullptr

non-owning pointer to host sine coefficients

◆ m_cos_h_data

amrex::ParticleReal* impactx::elements::SoftSolenoid::m_cos_h_data = nullptr

number of Fourier coefficients

◆ m_id

int impactx::elements::SoftSolenoid::m_id

number of map integration steps per slice

◆ m_mapsteps

int impactx::elements::SoftSolenoid::m_mapsteps

unit specification for quad strength

◆ m_ncoef

int impactx::elements::SoftSolenoid::m_ncoef = 0

unique soft solenoid id used for data lookup map

◆ m_sin_d_data

amrex::ParticleReal* impactx::elements::SoftSolenoid::m_sin_d_data = nullptr

non-owning pointer to device cosine coefficients

◆ m_sin_h_data

amrex::ParticleReal* impactx::elements::SoftSolenoid::m_sin_h_data = nullptr

non-owning pointer to host cosine coefficients

◆ m_unit

int impactx::elements::SoftSolenoid::m_unit

scaling factor for solenoid Bz field

◆ type

auto impactx::elements::SoftSolenoid::type = "SoftSolenoid"
staticconstexpr

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