Navigation

Main Page
Download
Support
Installation
Tutorial
Examples
Reference Manual
   Version 1.0.6dev
      Class Index
      File List
   Version 1.0.6
   Version 1.0.5new_solver
   Version 1.0.5dev
   Version 1.0.5b
   Version 1.0.4dev
   Version 1.0.4
Publications


Hosted by Get Ion Beam Simulator at SourceForge.net. Fast, secure and Free Open Source software downloads
ParticleDataBase3D Class Reference

Particle database class for three dimensions. More...

#include <particledatabase.hpp>

Inheritance diagram for ParticleDataBase3D:
ParticleDataBase

Public Member Functions

 ParticleDataBase3D (const Geometry &geom)
 Constructor. More...
 
 ParticleDataBase3D (const ParticleDataBase3D &pdb)
 Copy constructor. More...
 
 ParticleDataBase3D (std::istream &s, const Geometry &geom)
 Constructor for loading particle statistics from a file. More...
 
 ~ParticleDataBase3D ()
 Destructor. More...
 
const ParticleDataBase3Doperator= (const ParticleDataBase3D &pdb)
 Assignment. More...
 
virtual Particle3Dparticle (uint32_t i)
 Returns a reference to particle i. More...
 
virtual const Particle3Dparticle (uint32_t i) const
 Returns a const reference to particle i. More...
 
virtual const ParticleP3Dtrajectory_point (uint32_t i, uint32_t j) const
 Gets the particle i trajectory point j as particle point. More...
 
void trajectories_at_plane (std::vector< Particle3D > &tdata, coordinate_axis_e axis, double val) const
 Gets trajectory data at plane axis = val as particles in vector tdata. More...
 
void add_particle (double IQ, double q, double m, const ParticleP3D &x)
 Add one particle. More...
 
void add_particle (const Particle3D &p)
 Add one particle. More...
 
void add_cylindrical_beam_with_total_energy (uint32_t N, double J, double q, double m, double Etot, const ScalarField &epot, double Tp, double Tt, Vec3D c, Vec3D dir1, Vec3D dir2, double r)
 Add a cylindrical beam with energies. More...
 
void add_cylindrical_beam_with_energy (uint32_t N, double J, double q, double m, double E, double Tp, double Tt, Vec3D c, Vec3D dir1, Vec3D dir2, double r)
 Add a cylindrical beam with energies. More...
 
void add_cylindrical_beam_with_velocity (uint32_t N, double J, double q, double m, double v, double dvp, double dvt, Vec3D c, Vec3D dir1, Vec3D dir2, double r)
 Add a cylindrical beam with velocities. More...
 
void add_rectangular_beam_with_energy (uint32_t N, double J, double q, double m, double E, double Tp, double Tt, Vec3D c, Vec3D dir1, Vec3D dir2, double size1, double size2)
 Add a rectangular beam with energies. More...
 
void add_rectangular_beam_with_velocity (uint32_t N, double J, double q, double m, double v, double dvp, double dvt, Vec3D c, Vec3D dir1, Vec3D dir2, double size1, double size2)
 Add a rectangular beam with velocities. More...
 
void add_3d_KV_beam_with_emittance (uint32_t N, double I, double q, double m, double E0, double a1, double b1, double e1, double a2, double b2, double e2, Vec3D c, Vec3D dir1, Vec3D dir2)
 Add a 3d beam with defined KV emittance. More...
 
void add_3d_waterbag_beam_with_emittance (uint32_t N, double I, double q, double m, double E0, double a1, double b1, double e1, double a2, double b2, double e2, Vec3D c, Vec3D dir1, Vec3D dir2)
 Add a 3d beam with defined waterbag emittance. More...
 
void add_3d_gaussian_beam_with_emittance (uint32_t N, double I, double q, double m, double E0, double a1, double b1, double e1, double a2, double b2, double e2, Vec3D c, Vec3D dir1, Vec3D dir2)
 Add a 3d beam with defined gaussian emittance. More...
 
void trajectories_at_free_plane (TrajectoryDiagnosticData &tdata, const Vec3D &c, const Vec3D &o, const Vec3D &p, const std::vector< trajectory_diagnostic_e > &diagnostics) const
 Get trajectory diagnostic on a plane. More...
 
void export_path_manager_data (const std::string &filename, double ref_E, double ref_q, double ref_m, const Vec3D &c, const Vec3D &o, const Vec3D &p) const
 Export particle data as Path Manager data. More...
 
virtual void save (const std::string &filename) const
 Saves data to a new file filename. More...
 
virtual void save (std::ostream &s) const
 Saves data to stream. More...
 
virtual void debug_print (std::ostream &os) const
 Print debugging information to os. More...
 
virtual const ParticlePBasetrajectory_point (uint32_t i, uint32_t j) const=0
 Gets the particle i trajectory point j as particle point. More...
 
void trajectory_point (double &t, Vec3D &loc, Vec3D &vel, uint32_t i, uint32_t j) const
 Gets the particle i trajectory point j into vel, loc and t. More...
 
void trajectories_at_plane (TrajectoryDiagnosticData &tdata, coordinate_axis_e axis, double val, const std::vector< trajectory_diagnostic_e > &diagnostics) const
 Gets trajectory diagnostic diagnostics at plane axis = val in trajectory diagnostic data object tdata. More...
 
- Public Member Functions inherited from ParticleDataBase
virtual ~ParticleDataBase ()
 Virtual destructor. More...
 
void set_thread_count (uint32_t threadcount)
 Set the number of threads used for calculation. More...
 
void set_accuracy (double epsabs, double epsrel)
 Set the accuracy requirement for calculation. More...
 
void set_bfield_suppression (const CallbackFunctorD_V *functor)
 Set magnetic field suppression location depedent callback functor. More...
 
void set_trajectory_handler_callback (TrajectoryHandlerCallback *thand_cb)
 Set trajectory handler callback. More...
 
void set_trajectory_end_callback (TrajectoryEndCallback *tend_cb)
 Set trajectory end callback. More...
 
void set_trajectory_surface_collision_callback (TrajectorySurfaceCollisionCallback *tsur_cb)
 Set trajectory surface collision callback. More...
 
void set_relativistic (bool enable)
 Set relativistic particle iteration. More...
 
void set_surface_collision (bool surface_collision)
 Set surface collision model to be used. More...
 
void set_polyint (bool polyint)
 Set the interpolation type to polynomial(true) or linear(false). More...
 
bool get_polyint (void) const
 Get current interpolation type. More...
 
void set_trajectory_interpolation (trajectory_interpolation_e intrp)
 Set trajectory interpolation type. More...
 
trajectory_interpolation_e get_trajectory_interpolation (void) const
 Get trajectory interpolation type. More...
 
void set_scharge_deposition (scharge_deposition_e type)
 Set space charge deposition type. More...
 
scharge_deposition_e get_scharge_deposition (void) const
 Get space charge deposition type. More...
 
void set_max_steps (uint32_t maxsteps)
 Set maximum number of steps to iterate. More...
 
void set_max_time (double maxt)
 Set maximum lifetime of particle in simulation. More...
 
void set_save_all_points (bool save_points)
 Set trajectory saving. More...
 
void set_save_trajectories (uint32_t div)
 Set trajectory saving. More...
 
uint32_t get_save_trajectories (void) const
 Get trajectory saving. More...
 
void set_mirror (const bool mirror[6])
 Set particle mirroring on boundaries. More...
 
void get_mirror (bool mirror[6]) const
 Get particle mirroring on boundaries. More...
 
int get_iteration_number (void) const
 Get the number of iteration rounds done. More...
 
double get_rhosum (void) const
 Return sum of defined beam space charge density. More...
 
void set_rhosum (double rhosum)
 Set sum of defined beam space charge density. More...
 
const ParticleStatisticsget_statistics (void) const
 Get particle iterator statistics. More...
 
geom_mode_e geom_mode () const
 Return geometry mode. More...
 
size_t size (void) const
 Returns particle count. More...
 
double traj_length (uint32_t i) const
 Returns the length of trajectory i. More...
 
size_t traj_size (uint32_t i) const
 Returns number of trajectory points for particle i. More...
 
void trajectory_point (double &t, Vec3D &loc, Vec3D &vel, uint32_t i, uint32_t j) const
 Gets the particle i trajectory point j into vel, loc and t. More...
 
void trajectories_at_plane (TrajectoryDiagnosticData &tdata, coordinate_axis_e axis, double val, const std::vector< trajectory_diagnostic_e > &diagnostics) const
 Gets trajectory diagnostic diagnostics at plane axis = val in trajectory diagnostic data object tdata. More...
 
void build_trajectory_density_field (MeshScalarField &tdens) const
 Build trajectory density field. More...
 
void clear (void)
 Clears the particle database of all particles. More...
 
void clear_trajectories (void)
 Clears the particle trajectory database. More...
 
void clear_trajectory (size_t a)
 Clears particle a in the particle trajectory database. More...
 
void reset_trajectories (void)
 Clears the particle trajectory database and set particles to initial coordinates. More...
 
void reset_trajectory (size_t a)
 Clears the particle trajectory and set the particle to initial coordinates. More...
 
void reserve (size_t size)
 Reserve memory for size particles. More...
 
void iterate_trajectories (MeshScalarField &scharge, const VectorField &efield, const VectorField &bfield)
 Iterate particles through the geometry. More...
 
void step_particles (MeshScalarField &scharge, const VectorField &efield, const VectorField &bfield, double dt)
 Step particles forward by time step dt. More...
 

Additional Inherited Members

- Protected Member Functions inherited from ParticleDataBase
 ParticleDataBase ()
 Constructor. More...
 
 ParticleDataBase (const ParticleDataBase &pdb)
 Copy constructor. More...
 
const ParticleDataBaseoperator= (const ParticleDataBase &pdb)
 Assignment. More...
 
void set_implementation_pointer (class ParticleDataBaseImp *imp)
 Set particle database implementation pointer. More...
 

Detailed Description

Particle database class for three dimensions.

Particle database holds the definitions of particles and possibly the trajectories of particles it the particle iterator has saved them. ParticleDataBase3D provides a variety of convenience functions for defining particle beams.

Particles are always stored in the database in the order they are defined. When reading back the simulation results, the order can be used to identify the particles.

Constructor & Destructor Documentation

◆ ParticleDataBase3D() [1/3]

ParticleDataBase3D::ParticleDataBase3D ( const Geometry geom)

Constructor.

◆ ParticleDataBase3D() [2/3]

ParticleDataBase3D::ParticleDataBase3D ( const ParticleDataBase3D pdb)

Copy constructor.

◆ ParticleDataBase3D() [3/3]

ParticleDataBase3D::ParticleDataBase3D ( std::istream &  s,
const Geometry geom 
)

Constructor for loading particle statistics from a file.

◆ ~ParticleDataBase3D()

ParticleDataBase3D::~ParticleDataBase3D ( )

Destructor.

Member Function Documentation

◆ add_3d_gaussian_beam_with_emittance()

void ParticleDataBase3D::add_3d_gaussian_beam_with_emittance ( uint32_t  N,
double  I,
double  q,
double  m,
double  E0,
double  a1,
double  b1,
double  e1,
double  a2,
double  b2,
double  e2,
Vec3D  c,
Vec3D  dir1,
Vec3D  dir2 
)

Add a 3d beam with defined gaussian emittance.

Adds a beam consisting of N particles and a total current of I (A). The particles are defined to have equal currents, charge q (in multiples of e) and mass m (u). The starting energy of the beam is E0 (eV) and the starting location c (center point). The beam propagates to the direction dir3 = dir1 x dir2, where dir1 and dir2 are the emittance axes. The beam is made to match Twiss parameters $ \alpha_1 $ (a1), $ \beta_1 $ (b1), $ \epsilon_{1,\mathrm{rms}} $ (e1) in the first projectional direction dir1 and $ \alpha_2 $ (a2), $ \beta_2 $ (b2), $ \epsilon_{2,\mathrm{rms}} $ (e2) in the second projectional direction dir2.

The beam spread in the transverse phase space is made according to Gaussian distribution.

◆ add_3d_KV_beam_with_emittance()

void ParticleDataBase3D::add_3d_KV_beam_with_emittance ( uint32_t  N,
double  I,
double  q,
double  m,
double  E0,
double  a1,
double  b1,
double  e1,
double  a2,
double  b2,
double  e2,
Vec3D  c,
Vec3D  dir1,
Vec3D  dir2 
)

Add a 3d beam with defined KV emittance.

Adds a beam consisting of N particles and a total current of I (A). The particles are defined to have equal currents, charge q (in multiples of e) and mass m (u). The starting energy of the beam is E0 (eV) and the starting location c (center point). The beam propagates to the direction dir3 = dir1 x dir2, where dir1 and dir2 are the emittance axes. The beam is made to match Twiss parameters $ \alpha_1 $ (a1), $ \beta_1 $ (b1), $ \epsilon_{1,\mathrm{rms}} $ (e1) in the first projectional direction dir1 and $ \alpha_2 $ (a2), $ \beta_2 $ (b2), $ \epsilon_{2,\mathrm{rms}} $ (e2) in the second projectional direction dir2.

The beam spread in the transverse phase space is made according to KV/hard-edged (Kapchinsky-Vladimirsky) distribution.

◆ add_3d_waterbag_beam_with_emittance()

void ParticleDataBase3D::add_3d_waterbag_beam_with_emittance ( uint32_t  N,
double  I,
double  q,
double  m,
double  E0,
double  a1,
double  b1,
double  e1,
double  a2,
double  b2,
double  e2,
Vec3D  c,
Vec3D  dir1,
Vec3D  dir2 
)

Add a 3d beam with defined waterbag emittance.

Adds a beam consisting of N particles and a total current of I (A). The particles are defined to have equal currents, charge q (in multiples of e) and mass m (u). The starting energy of the beam is E0 (eV) and the starting location c (center point). The beam propagates to the direction dir3 = dir1 x dir2, where dir1 and dir2 are the emittance axes. The beam is made to match Twiss parameters $ \alpha_1 $ (a1), $ \beta_1 $ (b1), $ \epsilon_{1,\mathrm{rms}} $ (e1) in the first projectional direction dir1 and $ \alpha_2 $ (a2), $ \beta_2 $ (b2), $ \epsilon_{2,\mathrm{rms}} $ (e2) in the second projectional direction dir2.

The beam spread in the transverse phase space is made according to waterbag distribution.

◆ add_cylindrical_beam_with_energy()

void ParticleDataBase3D::add_cylindrical_beam_with_energy ( uint32_t  N,
double  J,
double  q,
double  m,
double  E,
double  Tp,
double  Tt,
Vec3D  c,
Vec3D  dir1,
Vec3D  dir2,
double  r 
)

Add a cylindrical beam with energies.

Adds a beam consisting of N particles. The beam current density is J (A/m^2), charge of beam particles is q (in multiples of e), mass is m (u). The beam starting surface is a disc of radius r centered at c. The normal direction of the disc is dir3 = dir1 x dir2. The first tangent direction is dir1 and the second is dir1 x dir3. If you want beam to go to positive x-direction, dir1 could be (0,1,0) and dir2 (0,0,1) for example. The beam energy E (eV) is defined in the normal direction. Temperatures are defined in normal (parallel) direction and transverse direction as Tp (eV) and Tt (eV), respectively.

The particle speeds of the beam in direction i are sampled from a gaussian distribution with standard deviation dv_i = sqrt(T_i*e/m), where T_i is the beam temperature in direction (eV), e is electron charge (C) and m is the mass of the ion (kg).

Space charge J/v is added to the rhosum variable.

◆ add_cylindrical_beam_with_total_energy()

void ParticleDataBase3D::add_cylindrical_beam_with_total_energy ( uint32_t  N,
double  J,
double  q,
double  m,
double  Etot,
const ScalarField epot,
double  Tp,
double  Tt,
Vec3D  c,
Vec3D  dir1,
Vec3D  dir2,
double  r 
)

Add a cylindrical beam with energies.

Adds a beam consisting of N particles. The beam current density is J (A/m^2), charge of beam particles is q (in multiples of e), mass is m (u). The beam starting surface is a disc of radius r centered at c. The normal direction of the disc is dir3 = dir1 x dir2. The first tangent direction is dir1 and the second is dir1 x dir3. If you want beam to go to positive x-direction, dir1 could be (0,1,0) and dir2 (0,0,1) for example. The beam total energy Etot (eV) is defined in the normal direction. The total energy is calculated as Etot = 0.5mv^2 + qU. Temperatures are defined in normal (parallel) direction and transverse direction as Tp (eV) and Tt (eV), respectively.

The particle speeds of the beam in direction i are sampled from a gaussian distribution with standard deviation dv_i = sqrt(T_i*e/m), where T_i is the beam temperature in direction (eV), e is electron charge (C) and m is the mass of the ion (kg).

◆ add_cylindrical_beam_with_velocity()

void ParticleDataBase3D::add_cylindrical_beam_with_velocity ( uint32_t  N,
double  J,
double  q,
double  m,
double  v,
double  dvp,
double  dvt,
Vec3D  c,
Vec3D  dir1,
Vec3D  dir2,
double  r 
)

Add a cylindrical beam with velocities.

Adds a beam consisting of N particles. The beam current density is J (A/m^2), charge of beam particles is q (in multiples of e), mass is m (u). The beam starting surface is a disc of radius r centered at c. The normal direction of the disc is dir3 = dir1 x dir2. The first tangent direction is dir1 and the second is dir1 x dir3. The beam velocity v (m/s) is defined in the normal direction. The velocity distribution standard deviation are defined in normal (parallel) direction and transverse directions as dvp (m/s) and dvt (m/s), respectively

The particle velocities of the beam in direction i are sampled from a gaussian distribution with standard deviation dv_i, which is related to beam temperature by dv_i = sqrt(T_i*e/m), where T_i is the beam temperature in direction (eV), e is electron charge (C) and m is the mass of the ion (kg).

Space charge J/v is added to the rhosum variable.

◆ add_particle() [1/2]

void ParticleDataBase3D::add_particle ( const Particle3D p)

Add one particle.

Adds one particle to database.

◆ add_particle() [2/2]

void ParticleDataBase3D::add_particle ( double  IQ,
double  q,
double  m,
const ParticleP3D x 
)

Add one particle.

Adds one particle to particle database. Particle properties are: IQ is the current (A) in time-independent or charge (C) in time-dependent simulations carried by the particle cloud that the simulated particle represents, q is the charge state of the microscopic particle (in multiples of e), m is the mass of the microscopic particle (u) and x contains the time, position (m) and velocity (m/s) of the particle.

◆ add_rectangular_beam_with_energy()

void ParticleDataBase3D::add_rectangular_beam_with_energy ( uint32_t  N,
double  J,
double  q,
double  m,
double  E,
double  Tp,
double  Tt,
Vec3D  c,
Vec3D  dir1,
Vec3D  dir2,
double  size1,
double  size2 
)

Add a rectangular beam with energies.

Very much like add_cylindrical_beam_with_energy(), but the beam shape is rectangular. The beam size is defined with two variables, sizex and sizey. These sizes are half the range, i.e. the beam ranges in the first direction from center point - size1 to center point + size1.

◆ add_rectangular_beam_with_velocity()

void ParticleDataBase3D::add_rectangular_beam_with_velocity ( uint32_t  N,
double  J,
double  q,
double  m,
double  v,
double  dvp,
double  dvt,
Vec3D  c,
Vec3D  dir1,
Vec3D  dir2,
double  size1,
double  size2 
)

Add a rectangular beam with velocities.

Very much like add_cylindrical_beam_with_velocity(), but the beam shape is rectangular. The beam size is defined with two variables, sizex and sizey. These sizes are half the range, i.e. the beam ranges in the first direction from center point - size1 to center point + size1.

◆ debug_print()

void ParticleDataBase3D::debug_print ( std::ostream &  os) const
virtual

Print debugging information to os.

Implements ParticleDataBase.

◆ export_path_manager_data()

void ParticleDataBase3D::export_path_manager_data ( const std::string &  filename,
double  ref_E,
double  ref_q,
double  ref_m,
const Vec3D c,
const Vec3D o,
const Vec3D p 
) const

Export particle data as Path Manager data.

Makes trajectory diagnostics on a plane defined in three dimensional space by center point vector c and two basis vectors o and p. The vectors o and p are normalized and p is adjusted to be orthogonal to o. This is done by calculating q = cross(o,p) and p = cross(q,o).

The particle properties on the diagnostics plane are written to file filename in the format required by the Path Manager program. Reference particle for the output has energy ref_E (in eV), charge state ref_q (in electron charges) and mass ref_m (in atomic mass units).

The file contains a header with information about the reference particle. After the header, the file has a line for each particle, which contain the following columns: x, x', y, y', phase, (p_z-p_ref)/p_ref, 0-flag, q, m. The particles are mirrored according the particle database mirroring settings. The phase is randomized for each particle between 0 and 2 pi corresponding to continuous beam.

◆ operator=()

const ParticleDataBase3D & ParticleDataBase3D::operator= ( const ParticleDataBase3D pdb)

Assignment.

◆ particle() [1/2]

Particle3D & ParticleDataBase3D::particle ( uint32_t  i)
virtual

Returns a reference to particle i.

Implements ParticleDataBase.

◆ particle() [2/2]

const Particle3D & ParticleDataBase3D::particle ( uint32_t  i) const
virtual

Returns a const reference to particle i.

Implements ParticleDataBase.

◆ save() [1/2]

void ParticleDataBase3D::save ( const std::string &  filename) const
virtual

Saves data to a new file filename.

Implements ParticleDataBase.

◆ save() [2/2]

void ParticleDataBase3D::save ( std::ostream &  s) const
virtual

Saves data to stream.

Implements ParticleDataBase.

◆ trajectories_at_free_plane()

void ParticleDataBase3D::trajectories_at_free_plane ( TrajectoryDiagnosticData tdata,
const Vec3D c,
const Vec3D o,
const Vec3D p,
const std::vector< trajectory_diagnostic_e > &  diagnostics 
) const

Get trajectory diagnostic on a plane.

Builds trajectory diagnostics on a plane defined in three dimensional space by center point vector c and two basis vectors o and p. The vectors o and p are normalized and p is adjusted to be orthogonal to o. This is done by calculating q = cross(o,p) and p = cross(q,o).

The diagnostic gathered is defined by diagnostics. Data is returned in tdata.

Valid diagnostic types are DIAG_T, DIAG_X, DIAG_VX, DIAG_Y, DIAG_VY, DIAG_Z, DIAG_VZ, DIAG_O, DIAG_VO, DIAG_P, DIAG_VP, DIAG_Q, DIAG_VQ, DIAG_OP, DIAG_PP, DIAG_CURR, DIAG_EK and DIAG_QM, DIAG_MASS and DIAG_CHARGE.

◆ trajectories_at_plane() [1/2]

void ParticleDataBase3D::trajectories_at_plane ( std::vector< Particle3D > &  tdata,
coordinate_axis_e  axis,
double  val 
) const

Gets trajectory data at plane axis = val as particles in vector tdata.

◆ trajectories_at_plane() [2/2]

void ParticleDataBase::trajectories_at_plane

Gets trajectory diagnostic diagnostics at plane axis = val in trajectory diagnostic data object tdata.

◆ trajectory_point() [1/3]

void ParticleDataBase::trajectory_point

Gets the particle i trajectory point j into vel, loc and t.

◆ trajectory_point() [2/3]

const ParticleP3D & ParticleDataBase3D::trajectory_point ( uint32_t  i,
uint32_t  j 
) const
virtual

Gets the particle i trajectory point j as particle point.

Implements ParticleDataBase.

◆ trajectory_point() [3/3]

virtual const ParticlePBase& ParticleDataBase::trajectory_point

Gets the particle i trajectory point j as particle point.


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


Reference manual for Ion Beam Simulator 1.0.6dev
Generated by Doxygen 1.9.1 on Thu Sep 11 2025 09:37:24.