diff --git a/include/Particle.hpp b/include/Particle.hpp index 682859c..3f5d182 100644 --- a/include/Particle.hpp +++ b/include/Particle.hpp @@ -3,7 +3,7 @@ * @author Cory Alexander Balaton (coryab) * @author Janita Ovidie Sandtrøen Willumsen (janitaws) * - * @version 0.1 + * @version 1.0 * * @brief A class that holds the properties of a particle. * @@ -22,10 +22,10 @@ class Particle { private: - vec_3d r_vec; ///< position - vec_3d v_vec; ///< velocity - double q; ///< Charge - double m; ///< Mass + vec3 r_vec; ///< position + vec3 v_vec; ///< velocity + double q; ///< Charge + double m; ///< Mass public: /** @brief Initialize the particle. @@ -38,8 +38,7 @@ public: * @param r_vec The initial position of the particle * @param v_vec The initial velocity of the particle * */ - Particle(vec_3d r_vec, vec_3d v_vec, double q = CA_CHARGE, - double m = CA_MASS); + Particle(vec3 r_vec, vec3 v_vec, double q = CA_CHARGE, double m = CA_MASS); /** @brief Make private attributes available for PenningTrap. * */ diff --git a/include/PenningTrap.hpp b/include/PenningTrap.hpp index 425afeb..0581e50 100644 --- a/include/PenningTrap.hpp +++ b/include/PenningTrap.hpp @@ -3,7 +3,7 @@ * @author Cory Alexander Balaton (coryab) * @author Janita Ovidie Sandtrøen Willumsen (janitaws) * - * @version 0.1 + * @version 1.0 * * @brief A class for simulating a Penning trap. * @@ -20,7 +20,7 @@ #include "typedefs.hpp" #include "utils.hpp" -#pragma omp declare reduction(+ : vec_3d : omp_out += omp_in) \ +#pragma omp declare reduction(+ : vec3 : omp_out += omp_in) \ initializer(omp_priv = omp_orig) /** @brief A class that simulates a Penning trap. @@ -28,14 +28,15 @@ * This class simulates a Penning trap. It can take in a number of particles * and simulate how they would behave inside a Penning trap. * */ -class PenningTrap { +class PenningTrap +{ private: - double B_0; ///< Magnetic field strength - double V_0; ///< Applied potential + double B_0; ///< Magnetic field strength + double V_0; ///< Applied potential std::function perturbation; ///< Time-dependent perturbation - double d; ///< Characteristic dimension - double t; ///< Current time - std::vector particles; ///< The particles in the Penning trap + double d; ///< Characteristic dimension + double t; ///< Current time + std::vector particles; ///< The particles in the Penning trap sim_arr k_v; ///< A 2D vector containing all \f$k_{i,j}\f$ where \f$j\f$ is ///< the index of a particle sim_arr k_r; ///< A 2D vector containing all \f$k_{i,j}\f$ where \f$j\f$ is @@ -49,9 +50,9 @@ private: * @param j Index j for \f$k_{v,i,j}\f$ * @param dt the step length (delta time) * - * @return vec_3d + * @return vec3 * */ - vec_3d v_func(uint i, uint j, double dt); + vec3 v_func(uint i, uint j, double dt); /** @brief Helper for evolve_RK4 when calculating \f$k_{r,i,j}\f$ values * @@ -61,82 +62,25 @@ private: * @param j Index j for \f$k_{r,i,j}\f$ * @param dt The step length (delta time) * - * @return vec_3d + * @return vec3 * */ - vec_3d r_func(uint i, uint j, double dt); - -public: - /** @brief Constructor for the PenningTrap class - * - * @param B_0 The magnetic field strength - * @param V_0 The time dependent applied potential - * @param d The characteristic dimension - * @param t The starting time - * */ - PenningTrap( - double B_0 = T, - double V_0 = (25. * V) / 1000., - double d = 500., double t = 0.); - - /** @brief Constructor for the PenningTrap class - * - * @param i The number of particles to generate - * @param B_0 The magnetic field strength - * @param V_0 The time dependent applied potential - * @param d The characteristic dimension - * @param t The starting time - * */ - PenningTrap( - uint i, double B_0 = T, - double V_0 = (25. * V) / 1000., - double d = 500., double t = 0.); - - /** @brief Constructor for the PenningTrap class - * - * @param particles The starting particles - * @param B_0 The magnetic field strength - * @param V_0 The time dependent applied potential - * @param d The characteristic dimension - * @param t The starting time - * */ - PenningTrap( - std::vector particles, double B_0 = T, - double V_0 = (25. * V) / 1000., - double d = 500., double t = 0.); - - void set_pertubation(double f, double omega_V); - - /** @brief Give all particles new positions and velocities, and change t - * and V_0. - * - * @param V_0 The tiome dependent applied potential - * @param t The starting time - * */ - void reinitialize( - double f, double omega_V, - double t = 0.); - - /** @brief Add a particle to the system - * - * @param particle The particle to add to the Penning trap - * */ - void add_particle(Particle particle); + vec3 r_func(uint i, uint j, double dt); /** @brief Calculate E at point r * * @param r The position where we want to calculate the E field * - * @return vec_3d + * @return vec3 * */ - vec_3d external_E_field(vec_3d r); + vec3 external_E_field(vec3 r); /** @brief Calculate B at point r * * @param r The position where we want to calculate the B field * - * @return vec_3d + * @return vec3 * */ - vec_3d external_B_field(vec_3d r); + vec3 external_B_field(vec3 r); /** @brief Calculate the force between 2 particles. * @@ -146,9 +90,9 @@ public: * @param i The index of particle p_i * @param j The index of particle p_j * - * @return vec_3d + * @return vec3 * */ - vec_3d force_on_particle(uint i, uint j); + vec3 force_on_particle(uint i, uint j); /** @brief Calculate the total external force on a particle. * @@ -157,26 +101,80 @@ public: * * @param i The index of particle p_i * - * @return vec_3d + * @return vec3 * */ - vec_3d total_force_external(uint i); + vec3 total_force_external(uint i); /** @brief Calculate the total force on a particle p_i from other * particles. * * @param i The index of particle p_i * - * @return vec_3d + * @return vec3 * */ - vec_3d total_force_particles(uint i); + vec3 total_force_particles(uint i); /** @brief calculate the total force on a particle p_i. * * @param i The index of particle p_i * - * @return vec_3d + * @return vec3 * */ - vec_3d total_force(uint i); + vec3 total_force(uint i); + +public: + /** @brief Constructor for the PenningTrap class + * + * @param B_0 The magnetic field strength + * @param V_0 The time dependent applied potential + * @param d The characteristic dimension + * @param t The starting time + * */ + PenningTrap(double B_0 = T, double V_0 = (25. * V) / 1000., double d = 500., + double t = 0.); + + /** @brief Constructor for the PenningTrap class + * + * @param i The number of particles to generate + * @param B_0 The magnetic field strength + * @param V_0 The time dependent applied potential + * @param d The characteristic dimension + * @param t The starting time + * */ + PenningTrap(uint i, double B_0 = T, double V_0 = (25. * V) / 1000., + double d = 500., double t = 0.); + + /** @brief Constructor for the PenningTrap class + * + * @param particles The starting particles + * @param B_0 The magnetic field strength + * @param V_0 The time dependent applied potential + * @param d The characteristic dimension + * @param t The starting time + * */ + PenningTrap(std::vector particles, double B_0 = T, + double V_0 = (25. * V) / 1000., double d = 500., double t = 0.); + + /** @brief Time dependent perturbation to V_0 + * + * @param f The amplitude of the perturbation + * @parma omega_V the angular frequency of the perturbation + * */ + void set_pertubation(double f, double omega_V); + + /** @brief Give all particles new positions and velocities, and change t + * and V_0. + * + * @param V_0 The tiome dependent applied potential + * @param t The starting time + * */ + void reinitialize(double f, double omega_V, double t = 0.); + + /** @brief Add a particle to the system + * + * @param particle The particle to add to the Penning trap + * */ + void add_particle(Particle particle); /** @brief Go forward one timestep using the RK4 method * @@ -202,8 +200,7 @@ public: * * @return simulation_t * */ - simulation_t simulate(double time, uint steps, - std::string method = "rk4", + simulation_t simulate(double time, uint steps, std::string method = "rk4", bool particle_interaction = true); /** @brief Simulate and write the displacement of all particles to files. @@ -214,8 +211,7 @@ public: * @param method The method to use when moving forward a timestep * @param particle_interaction Turn particle interactions on/off * */ - void write_simulation_to_dir(std::string path, double time, - uint steps, + void write_simulation_to_dir(std::string path, double time, uint steps, std::string method = "rk4", bool particle_interaction = true); @@ -232,6 +228,8 @@ public: double fraction_of_particles_left(double time, uint steps, std::string method = "rk4", bool particle_interaction = true); + + friend class PenningTrapTest; }; #endif diff --git a/include/constants.hpp b/include/constants.hpp index 0546a4f..85fc7d2 100644 --- a/include/constants.hpp +++ b/include/constants.hpp @@ -3,7 +3,7 @@ * @author Cory Alexander Balaton (coryab) * @author Janita Ovidie Sandtrøen Willumsen (janitaws) * - * @version 0.1 + * @version 1.0 * * @brief Library of constants * diff --git a/include/typedefs.hpp b/include/typedefs.hpp index 64234e8..688e942 100644 --- a/include/typedefs.hpp +++ b/include/typedefs.hpp @@ -18,23 +18,23 @@ #include #include +/** @brief Typedef for a fixed 3d arma vector. + * */ +typedef arma::vec::fixed<3> vec3; + /** @brief Typedef for the column of the result vector from simulating * particles. * */ -typedef std::vector> sim_cols; +typedef std::vector sim_cols; /** @brief Typedef for the row of the result vector from simulating particles. * */ -typedef std::vector> sim_rows; +typedef std::vector sim_rows; /** @brief Typedef for the result of the simulate method. * */ typedef std::vector sim_arr; -/** @brief Typedef for a fixed 3d arma vector. - * */ -typedef arma::vec::fixed<3> vec_3d; - /** @brief Typedef for PenningTrap::simulation return value. * */ typedef struct simulation { diff --git a/include/utils.hpp b/include/utils.hpp index e2e28c2..8816e8f 100644 --- a/include/utils.hpp +++ b/include/utils.hpp @@ -31,7 +31,7 @@ * messages by adding the -DDBG flag when compiling. * */ #ifdef DBG -#define DEBUG(msg) \ +#define DEBUG(msg) \ std::cout << __FILE__ << " " << __LINE__ << ": " << msg << std::endl #else #define DEBUG(msg) @@ -43,7 +43,7 @@ * This macro calls the m_assert function which is a more informative * assertion function than the regular assert function from cassert. * */ -#define ASSERT(expr, msg) \ +#define ASSERT(expr, msg) \ m_assert(expr, #expr, __METHOD_NAME__, __FILE__, __LINE__, msg) /** @def __METHOD_NAME__