ImpactX
beamoptic.H
Go to the documentation of this file.
1 /* Copyright 2022-2023 The Regents of the University of California, through Lawrence
2  * Berkeley National Laboratory (subject to receipt of any required
3  * approvals from the U.S. Dept. of Energy). All rights reserved.
4  *
5  * This file is part of ImpactX.
6  *
7  * Authors: Axel Huebl
8  * License: BSD-3-Clause-LBNL
9  */
10 #ifndef IMPACTX_ELEMENTS_MIXIN_BEAMOPTIC_H
11 #define IMPACTX_ELEMENTS_MIXIN_BEAMOPTIC_H
12 
14 #include "particles/PushAll.H"
15 
16 #include <AMReX_Extension.H> // for AMREX_RESTRICT
17 #include <AMReX_REAL.H>
18 
19 #include <type_traits>
20 
21 
22 namespace impactx::elements
23 {
24 namespace detail
25 {
39  template <typename T_Element>
41  {
43 
53  PushSingleParticle (T_Element element,
54  PType* AMREX_RESTRICT aos_ptr,
55  amrex::ParticleReal* AMREX_RESTRICT part_px,
56  amrex::ParticleReal* AMREX_RESTRICT part_py,
57  amrex::ParticleReal* AMREX_RESTRICT part_pt,
58  RefPart ref_part)
59  : m_element(std::move(element)), m_aos_ptr(aos_ptr),
60  m_part_px(part_px), m_part_py(part_py), m_part_pt(part_pt),
61  m_ref_part(ref_part)
62  {
63  }
64 
65  PushSingleParticle () = delete;
68  ~PushSingleParticle () = default;
69 
75  void
76  operator() (long i) const
77  {
78  // access AoS data such as positions and cpu/id
80 
81  // access SoA Real data
82  amrex::ParticleReal & AMREX_RESTRICT px = m_part_px[i];
83  amrex::ParticleReal & AMREX_RESTRICT py = m_part_py[i];
84  amrex::ParticleReal & AMREX_RESTRICT pt = m_part_pt[i];
85 
86  // push through element
87  m_element(p, px, py, pt, m_ref_part);
88 
89  }
90 
91  private:
92  T_Element const m_element;
94  amrex::ParticleReal* const AMREX_RESTRICT m_part_px;
95  amrex::ParticleReal* const AMREX_RESTRICT m_part_py;
96  amrex::ParticleReal* const AMREX_RESTRICT m_part_pt;
98  };
99 
102  template< typename T_Element >
105  RefPart & AMREX_RESTRICT ref_part,
106  T_Element & element
107  ) {
108  const int np = pti.numParticles();
109 
110  // preparing access to particle data: AoS
112  auto& aos = pti.GetArrayOfStructs();
113  PType* AMREX_RESTRICT aos_ptr = aos().dataPtr();
114 
115  // preparing access to particle data: SoA of Reals
116  auto& soa_real = pti.GetStructOfArrays().GetRealData();
117  amrex::ParticleReal* const AMREX_RESTRICT part_px = soa_real[RealSoA::px].dataPtr();
118  amrex::ParticleReal* const AMREX_RESTRICT part_py = soa_real[RealSoA::py].dataPtr();
119  amrex::ParticleReal* const AMREX_RESTRICT part_pt = soa_real[RealSoA::pt].dataPtr();
120 
121  detail::PushSingleParticle<T_Element> const pushSingleParticle(
122  element, aos_ptr, part_px, part_py, part_pt, ref_part);
123  // loop over beam particles in the box
124  amrex::ParallelFor(np, pushSingleParticle);
125  }
126 } // namespace detail
127 
133  template<typename T_Element>
134  struct BeamOptic
135  {
139  int step
140  )
141  {
142  static_assert(
143  std::is_base_of_v<BeamOptic, T_Element>,
144  "BeamOptic can only be used as a mixin class!"
145  );
146 
147  T_Element& element = *static_cast<T_Element*>(this);
148  push_all(pc, element, step);
149  }
150 
160  RefPart & AMREX_RESTRICT ref_part
161  )
162  {
163  static_assert(
164  std::is_base_of_v<BeamOptic, T_Element>,
165  "BeamOptic can only be used as a mixin class!"
166  );
167 
168  T_Element& element = *static_cast<T_Element*>(this);
169  detail::push_all_particles<T_Element>(pti, ref_part, element);
170  }
171  };
172 
173 } // namespace impactx::elements
174 
175 #endif // IMPACTX_ELEMENTS_MIXIN_BEAMOPTIC_H
#define AMREX_FORCE_INLINE
#define AMREX_RESTRICT
#define AMREX_GPU_DEVICE
Definition: ImpactXParticleContainer.H:145
Definition: ImpactXParticleContainer.H:114
std::enable_if_t< std::is_integral< T >::value > ParallelFor(TypeList< CTOs... >, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
i
void push_all_particles(ImpactXParticleContainer::iterator &pti, RefPart &AMREX_RESTRICT ref_part, T_Element &element)
Definition: beamoptic.H:103
Definition: alignment.H:23
void push_all(ImpactXParticleContainer &pc, T_Element &element, [[maybe_unused]] int step, [[maybe_unused]] bool omp_parallel=true)
Definition: PushAll.H:32
@ pt
energy deviation, scaled by speed of light * the magnitude of the reference momentum [unitless] (at f...
Definition: ImpactXParticleContainer.H:77
@ px
momentum in x, scaled by the magnitude of the reference momentum [unitless] (at fixed s or t)
Definition: ImpactXParticleContainer.H:75
@ py
momentum in y, scaled by the magnitude of the reference momentum [unitless] (at fixed s or t)
Definition: ImpactXParticleContainer.H:76
Definition: ReferenceParticle.H:30
Definition: beamoptic.H:135
void operator()(ImpactXParticleContainer &pc, int step)
Definition: beamoptic.H:137
amrex::ParticleReal *const AMREX_RESTRICT m_part_px
Definition: beamoptic.H:94
amrex::ParticleReal *const AMREX_RESTRICT m_part_py
Definition: beamoptic.H:95
T_Element const m_element
Definition: beamoptic.H:92
PType *const AMREX_RESTRICT m_aos_ptr
Definition: beamoptic.H:93
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void operator()(long i) const
Definition: beamoptic.H:76
PushSingleParticle(PushSingleParticle const &)=default
amrex::ParticleReal *const AMREX_RESTRICT m_part_pt
Definition: beamoptic.H:96
ImpactXParticleContainer::ParticleType PType
Definition: beamoptic.H:42
RefPart const m_ref_part
Definition: beamoptic.H:97
PushSingleParticle(PushSingleParticle &&)=default
PushSingleParticle(T_Element element, PType *AMREX_RESTRICT aos_ptr, amrex::ParticleReal *AMREX_RESTRICT part_px, amrex::ParticleReal *AMREX_RESTRICT part_py, amrex::ParticleReal *AMREX_RESTRICT part_pt, RefPart ref_part)
Definition: beamoptic.H:53