Coryab/code #9

Merged
coryab merged 23 commits from coryab/code into develop 2023-10-22 15:46:55 +00:00
2 changed files with 69 additions and 35 deletions
Showing only changes of commit 76dd30df23 - Show all commits

View File

@ -26,9 +26,15 @@
#define MASS 40. // unit: amu
// Particles used for testing
Particle p1(CHARGE, MASS, vec_3d{20., 0., 20.}, vec_3d{0., 25., 0.});
Particle p2(CHARGE, MASS, vec_3d{25., 25., 0.}, vec_3d{0., 40., 5.});
Particle p1(CHARGE, MASS, vec_3d{20., 0., 20.}, vec_3d{0., 25., 0.}); ///< Particle 1
Particle p2(CHARGE, MASS, vec_3d{25., 25., 0.}, vec_3d{0., 40., 5.}); ///< Particle 2
/** @brief The analytical solution for particle p1
*
* @param t Time
*
* @return vec_3d
* */
vec_3d analytical_solution_particle_1(double t)
{
double w_0 = T / MASS;
@ -44,6 +50,8 @@ vec_3d analytical_solution_particle_1(double t)
return res;
}
/** @brief Simulate a single particle over the period of 50 \f$ \mu s \f$.
* */
void simulate_single_particle()
{
// Initialize trap with particle 1
@ -56,6 +64,9 @@ void simulate_single_particle()
"rk4", false);
}
/** @brief Simulate 2 particles over the period of 50 \f$ \mu s \f$ with and
* without particle interactions.
* */
void simulate_two_particles()
{
// Initialize traps with particles
@ -71,6 +82,9 @@ void simulate_two_particles()
"output/simulate_2_particles/with_interaction", time, N);
}
/** @brief Simulate a single particle over 50 \f$ \mu s \f$ using different
* amount of steps and different methods.
* */
void simulate_single_particle_with_different_steps()
{
double time = 50.; // microseconds
@ -112,16 +126,23 @@ void simulate_single_particle_with_different_steps()
}
}
/** @brief Simulate 100 particles over 50 \f$ \mu s \f$.
* */
void simulate_100_particles()
{
PenningTrap trap((unsigned)100);
double time = 50.; // microseconds
trap.write_simulation_to_dir("output/simulate_100_particles", time, N, "rk4", false);
trap.write_simulation_to_dir("output/simulate_100_particles", time, N);
}
// Wide sweep
/** @brief Simulate 100 particles over 500 \f$ \mu s \f$ using a time
* dependent potential.
*
* @details The simulation sweeps over different frequencies in [0.2, 2.5] MHz.
*
* */
void simulate_100_particles_with_time_potential()
{
double time = 500.;

View File

@ -17,32 +17,36 @@
#include "PenningTrap.hpp"
#include "utils.hpp"
/** @brief Test class for the Penning trap.
* */
class PenningTrapTest {
public:
/** @brief Test that the external E field gives correct values.
* */
static void test_external_E_field()
{
PenningTrap trap;
// Vector containing inputs and expected results
std::vector<std::pair<arma::vec, arma::vec>> tests;
std::vector<std::pair<vec_3d, vec_3d>> tests;
tests.push_back(
std::make_pair(arma::vec{0., 0., 0.}, arma::vec{0., 0., 0.}));
std::make_pair(vec_3d{0., 0., 0.}, vec_3d{0., 0., 0.}));
tests.push_back(std::make_pair(arma::vec{10., 0., 0.},
arma::vec{96.4852558, 0., 0.}));
tests.push_back(std::make_pair(vec_3d{10., 0., 0.},
vec_3d{96.4852558, 0., 0.}));
tests.push_back(std::make_pair(arma::vec{10., 0., 0.},
arma::vec{96.4852558, 0., 0.}));
tests.push_back(std::make_pair(vec_3d{10., 0., 0.},
vec_3d{96.4852558, 0., 0.}));
tests.push_back(std::make_pair(arma::vec{0., 10., 0.},
arma::vec{0., 96.4852558, 0.}));
tests.push_back(std::make_pair(vec_3d{0., 10., 0.},
vec_3d{0., 96.4852558, 0.}));
tests.push_back(std::make_pair(arma::vec{0., 0., 10.},
arma::vec{0., 0., -192.9705116}));
tests.push_back(std::make_pair(vec_3d{0., 0., 10.},
vec_3d{0., 0., -192.9705116}));
arma::vec result;
arma::vec v;
vec_3d result;
vec_3d v;
std::stringstream msg;
for (size_t i = 0; i < tests.size(); i++) {
v = tests.at(i).first;
@ -56,75 +60,84 @@ public:
}
}
/** @brief Test that the external B field gives correct values.
* */
static void test_external_B_field()
{
// No point in testing at different points since it's not dependent
// on position.
PenningTrap trap;
arma::vec expected{0., 0., T};
arma::vec result = trap.external_B_field(arma::vec{0., 0., 0.});
vec_3d expected{0., 0., T};
vec_3d result = trap.external_B_field(vec_3d{0., 0., 0.});
ASSERT(close_to(expected, result),
"Testing the external B field at (0,0,0)");
}
/** @brief Test that the force between particles gives expected results.
* */
static void test_force_on_particle()
{
PenningTrap trap;
arma::vec v{0., 0., 0.};
vec_3d v{0., 0., 0.};
// Add particles to test
trap.add_particle(Particle(1., 40., arma::vec{0., 0., 0.}, v));
trap.add_particle(Particle(1., 40., arma::vec{1., 0., 0.}, v));
trap.add_particle(Particle(1., 40., arma::vec{0., 3., 4.}, v));
trap.add_particle(Particle(1., 40., vec_3d{0., 0., 0.}, v));
trap.add_particle(Particle(1., 40., vec_3d{1., 0., 0.}, v));
trap.add_particle(Particle(1., 40., vec_3d{0., 3., 4.}, v));
// Test p0 and p1
arma::vec expected{-1., 0., 0.};
arma::vec result = trap.force_on_particle(0, 1);
vec_3d expected{-1., 0., 0.};
vec_3d result = trap.force_on_particle(0, 1);
ASSERT(close_to(expected, result),
"Testing the force on a particle at (0,0,0) from a "
"particle at (1,0,0).");
// Test p0 and p2
expected = arma::vec{0, -.024, -.032};
expected = vec_3d{0, -.024, -.032};
result = trap.force_on_particle(0, 2);
ASSERT(close_to(expected, result),
"Testing the force on a particle at (0,0,0) from a "
"particle at (0,3,4).");
}
/** @brief Test that the total external force returns expected results
* */
static void test_total_force_external()
{
PenningTrap trap;
trap.add_particle(
Particle(1., 40., arma::vec{1., 2., 3.}, arma::vec{3., 4., 5.}));
Particle(1., 40., vec_3d{1., 2., 3.}, vec_3d{3., 4., 5.}));
arma::vec expected{395.58954878, -270.15871624, -57.89115348};
arma::vec result = trap.total_force_external(0);
vec_3d expected{395.58954878, -270.15871624, -57.89115348};
vec_3d result = trap.total_force_external(0);
ASSERT(close_to(expected, result),
"Testing the total external force on a particle at "
"(1,2,3) with velocity (3,4,5)");
}
/** @brief Test that the total force of all particles on a single particle
* returns expected results.
* */
static void test_total_force_particles()
{
PenningTrap trap;
trap.add_particle(
Particle(1., 40., arma::vec{0., 0., 0.}, arma::vec{0., 0., 0.}));
Particle(1., 40., vec_3d{0., 0., 0.}, vec_3d{0., 0., 0.}));
arma::vec expected{0., 0., 0.};
arma::vec result = trap.total_force_particles(0);
vec_3d expected{0., 0., 0.};
vec_3d result = trap.total_force_particles(0);
ASSERT(close_to(expected, result),
"Testing the total force of all particles on particle 0 "
"with only a single particle");
trap.add_particle(
Particle(1., 40., arma::vec{1., 0., 0.}, arma::vec{0., 0., 0.}));
Particle(1., 40., vec_3d{1., 0., 0.}, vec_3d{0., 0., 0.}));
trap.add_particle(
Particle(1., 40., arma::vec{0., 1., 0.}, arma::vec{0., 0., 0.}));
Particle(1., 40., vec_3d{0., 1., 0.}, vec_3d{0., 0., 0.}));
trap.add_particle(
Particle(1., 40., arma::vec{0., 0., 1.}, arma::vec{0., 0., 0.}));
Particle(1., 40., vec_3d{0., 0., 1.}, vec_3d{0., 0., 0.}));
expected = arma::vec(3, arma::fill::value(-3473.383325));
expected = vec_3d().fill(-3473.383325);
result = trap.total_force_particles(0);
ASSERT(close_to(expected, result),
"Testing the total force of all particles on particle 0 "