Use new constructors and methods
This commit is contained in:
parent
8df35f90f5
commit
13c662c148
133
src/main.cpp
133
src/main.cpp
@ -12,143 +12,84 @@
|
||||
|
||||
#include <fstream>
|
||||
#include <omp.h>
|
||||
#include <string>
|
||||
#include <sys/stat.h>
|
||||
#include <vector>
|
||||
|
||||
#include "PenningTrap.hpp"
|
||||
#include "utils.hpp"
|
||||
|
||||
#define PARTICLES 100
|
||||
#define N 10000
|
||||
#define CHARGE 1.
|
||||
#define MASS 40. // unit: amu
|
||||
|
||||
Particle p1(1., 40., arma::vec{20., 0., 20.}, arma::vec{0., 25., 0.});
|
||||
Particle p2(1., 40., arma::vec{25., 25., 0.}, arma::vec{0., 40., 5.});
|
||||
Particle p1(CHARGE, MASS, arma::vec{20., 0., 20.}, arma::vec{0., 25., 0.});
|
||||
Particle p2(CHARGE, MASS, arma::vec{25., 25., 0.}, arma::vec{0., 40., 5.});
|
||||
|
||||
void simulate_single_particle()
|
||||
{
|
||||
PenningTrap trap;
|
||||
|
||||
trap.add_particle(p1);
|
||||
PenningTrap trap(std::vector<Particle>{p1});
|
||||
|
||||
double time = 50.; // microseconds
|
||||
double dt = time / (double)N;
|
||||
|
||||
auto res = new arma::vec::fixed<3>[N];
|
||||
|
||||
for (int i=0; i<N;i++) {
|
||||
res[i] = trap.get_particle(0);
|
||||
trap.evolve_RK4(dt);
|
||||
}
|
||||
|
||||
mkdir("output", 0777);
|
||||
mkdir("output/simulate_single_particle", 0777);
|
||||
|
||||
std::ofstream ofile;
|
||||
|
||||
ofile.open("output/simulate_single_particle/out.txt");
|
||||
for (int i=0; i<N;i++) {
|
||||
arma::vec p = res[i];
|
||||
ofile << p(0) << "," << p(1) << "," << p(3) << "\n";
|
||||
}
|
||||
ofile.close();
|
||||
delete [] res;
|
||||
trap.write_simulation_to_dir("output/simulate_single_particle", time, N);
|
||||
}
|
||||
|
||||
void simulate_two_particles()
|
||||
{
|
||||
PenningTrap trap;
|
||||
|
||||
trap.add_particle(p1);
|
||||
trap.add_particle(p2);
|
||||
PenningTrap trap_no_interaction(std::vector<Particle>{p1, p2});
|
||||
PenningTrap trap_with_interaction(std::vector<Particle>{p1, p2});
|
||||
|
||||
double time = 50.; // microseconds
|
||||
double dt = time / (double)N;
|
||||
|
||||
auto res = new arma::vec::fixed<3>[2][N];
|
||||
|
||||
|
||||
int counter = 0;
|
||||
|
||||
// Get the path of all particles
|
||||
for (int j = 0; j < N; j++) {
|
||||
#pragma omp parallel for
|
||||
for (int i = 0; i < 2; i++) {
|
||||
res[i][j] = trap.get_particle(i);
|
||||
}
|
||||
trap.evolve_RK4(dt);
|
||||
trap_no_interaction.write_simulation_to_dir(
|
||||
"output/simulate_2_particles/no_interaction", time, N, "rk4", false);
|
||||
trap_with_interaction.write_simulation_to_dir(
|
||||
"output/simulate_2_particles/with_interaction", time, N);
|
||||
}
|
||||
|
||||
std::cout << counter << std::endl;
|
||||
void simulate_single_particle_with_different_steps()
|
||||
{
|
||||
|
||||
arma::vec::fixed<3> *cur_row;
|
||||
arma::vec::fixed<3> cur_elem;
|
||||
double time = 50; // microseconds
|
||||
|
||||
mkdir("output", 0777);
|
||||
mkdir("output/simulate_2_particles", 0777);
|
||||
|
||||
std::ofstream ofile;
|
||||
|
||||
// Write particle paths to file
|
||||
#pragma omp parallel for private(cur_row, cur_elem, ofile)
|
||||
for (int i = 0; i < 2; i++) {
|
||||
cur_row = res[i];
|
||||
ofile.open("output/simulate_100_particles/p" + std::to_string(i) +
|
||||
".txt");
|
||||
for (int j = 0; j < N; j++) {
|
||||
cur_elem = cur_row[j];
|
||||
ofile << cur_elem(0) << "," << cur_elem(1) << "," << cur_elem(2)
|
||||
<< "\n";
|
||||
}
|
||||
ofile.close();
|
||||
for (int i = 0; i < 4; i++) {
|
||||
int steps = 4000 * (i + 1);
|
||||
PenningTrap trap(std::vector<Particle>{p1});
|
||||
trap.write_simulation_to_dir("output/N_steps/RK4/" + std::to_string(steps) +
|
||||
"_steps",
|
||||
time, steps, "rk4", false);
|
||||
}
|
||||
|
||||
for (int i = 0; i < 4; i++) {
|
||||
int steps = 4000 * (i + 1);
|
||||
PenningTrap trap(std::vector<Particle>{p1});
|
||||
trap.write_simulation_to_dir("output/N_steps/euler/" + std::to_string(steps) +
|
||||
"_steps",
|
||||
time, steps, "euler", false);
|
||||
}
|
||||
}
|
||||
|
||||
void simulate_100_particles()
|
||||
{
|
||||
PenningTrap trap;
|
||||
|
||||
// Add particles inside trap
|
||||
for (int i = 0; i < PARTICLES; i++) {
|
||||
arma::vec r = arma::vec(3).randn() * 0.1 *
|
||||
trap.get_d(); // random initial position
|
||||
arma::vec v = arma::vec(3).randn() * 0.1 *
|
||||
trap.get_d(); // random initial velocity
|
||||
trap.add_particle(Particle(CHARGE, MASS, r, v));
|
||||
}
|
||||
PenningTrap trap(100);
|
||||
|
||||
double time = 50.; // microseconds
|
||||
double dt = time / (double)N;
|
||||
|
||||
sim_arr res = trap.simulate(time, N);
|
||||
|
||||
sim_rows cur_row;
|
||||
arma::vec::fixed<3> cur_elem;
|
||||
|
||||
mkdir("output", 0777);
|
||||
mkdir("output/simulate_100_particles", 0777);
|
||||
|
||||
std::ofstream ofile;
|
||||
|
||||
// Write particle paths to file
|
||||
#pragma omp parallel for private(cur_row, cur_elem, ofile)
|
||||
for (int i = 0; i < PARTICLES; i++) {
|
||||
cur_row = res[i];
|
||||
ofile.open("output/simulate_100_particles/p" + std::to_string(i) +
|
||||
".txt");
|
||||
for (int j = 0; j < N; j++) {
|
||||
cur_elem = cur_row[j];
|
||||
ofile << cur_elem(0) << "," << cur_elem(1) << "," << cur_elem(2)
|
||||
<< "\n";
|
||||
}
|
||||
ofile.close();
|
||||
}
|
||||
trap.write_simulation_to_dir("output/simulate_100_particles", time, N);
|
||||
}
|
||||
|
||||
int main()
|
||||
{
|
||||
double start = omp_get_wtime();
|
||||
|
||||
simulate_single_particle();
|
||||
|
||||
simulate_two_particles();
|
||||
|
||||
simulate_single_particle_with_different_steps();
|
||||
|
||||
simulate_100_particles();
|
||||
|
||||
double end = omp_get_wtime();
|
||||
|
||||
Loading…
Reference in New Issue
Block a user