diff --git a/Doxyfile b/Doxyfile index ea0f35d..591e8ba 100644 --- a/Doxyfile +++ b/Doxyfile @@ -68,7 +68,7 @@ PROJECT_LOGO = # entered, it will be relative to the location where doxygen was started. If # left blank the current directory will be used. -OUTPUT_DIRECTORY = docs +OUTPUT_DIRECTORY = # If the CREATE_SUBDIRS tag is set to YES then doxygen will create up to 4096 # sub-directories (in 2 levels) under the output directory of each output format @@ -1303,7 +1303,7 @@ GENERATE_HTML = YES # The default directory is: html. # This tag requires that the tag GENERATE_HTML is set to YES. -HTML_OUTPUT = html +HTML_OUTPUT = docs # The HTML_FILE_EXTENSION tag can be used to specify the file extension for each # generated HTML page (for example: .htm, .php, .asp). @@ -2183,7 +2183,7 @@ GENERATE_MAN = YES # The default directory is: man. # This tag requires that the tag GENERATE_MAN is set to YES. -MAN_OUTPUT = man +MAN_OUTPUT = man_pages # The MAN_EXTENSION tag determines the extension that is added to the generated # man pages. In case the manual section does not start with a number, the number diff --git a/docs/Particle_8cpp.html b/docs/Particle_8cpp.html index 9f88834..98aaf8d 100644 --- a/docs/Particle_8cpp.html +++ b/docs/Particle_8cpp.html @@ -114,7 +114,6 @@ $(document).ready(function(){initNavTree('Particle_8cpp.html',''); initResizable Janita Ovidie Sandtrøen Willumsen (janitaws)
Version
0.1
Bug:
No known bugs
-
Todo:
Implement constructor

Definition in file Particle.cpp.

diff --git a/docs/Particle_8cpp_source.html b/docs/Particle_8cpp_source.html index 159fd35..2eaf72e 100644 --- a/docs/Particle_8cpp_source.html +++ b/docs/Particle_8cpp_source.html @@ -102,16 +102,24 @@ $(document).ready(function(){initNavTree('Particle_8cpp_source.html',''); initRe
Go to the documentation of this file.
1
-
15#include "Particle.hpp"
-
16
-
17Particle::Particle(double q, double m,
-
18 arma::vec::fixed<3> r_vec,
-
19 arma::vec::fixed<3> v_vec)
-
20{
-
21
-
22}
+
14#include "Particle.hpp"
+
15
+
16Particle::Particle(double q, double m,
+
17 arma::vec::fixed<3> r_vec,
+
18 arma::vec::fixed<3> v_vec)
+
19{
+
20 // Giving the particle its properties
+
21 this->q = q;
+
22 this->m = m;
+
23 this->r_vec = r_vec;
+
24 this->v_vec = v_vec;
+
25}
Particle.hpp
A class that holds the properties of a particle.
-
Particle::Particle
Particle(double q, double m, arma::vec::fixed< 3 > r_vec, arma::vec::fixed< 3 > v_vec)
Initialize the particle.
Definition: Particle.cpp:17
+
Particle::Particle
Particle(double q, double m, arma::vec::fixed< 3 > r_vec, arma::vec::fixed< 3 > v_vec)
Initialize the particle.
Definition: Particle.cpp:16
+
Particle::v_vec
arma::vec::fixed< 3 > v_vec
velocity
Definition: Particle.hpp:24
+
Particle::q
double q
Charge.
Definition: Particle.hpp:21
+
Particle::r_vec
arma::vec::fixed< 3 > r_vec
position
Definition: Particle.hpp:23
+
Particle::m
double m
Mass.
Definition: Particle.hpp:22
diff --git a/docs/PenningTrap_8cpp.html b/docs/PenningTrap_8cpp.html index e572253..ceb79d5 100644 --- a/docs/PenningTrap_8cpp.html +++ b/docs/PenningTrap_8cpp.html @@ -104,7 +104,12 @@ $(document).ready(function(){initNavTree('PenningTrap_8cpp.html',''); initResiza

The implementation of the PenningTrap class. More...

-
#include "PenningTrap.hpp"
+
#include "utils.hpp"
+#include "PenningTrap.hpp"
+#include "constants.hpp"
+#include <algorithm>
+#include <stdexcept>
+#include <omp.h>

Go to the source code of this file.

Detailed Description

@@ -114,15 +119,7 @@ $(document).ready(function(){initNavTree('PenningTrap_8cpp.html',''); initResiza Janita Ovidie Sandtrøen Willumsen (janitaws)
Version
0.1
Bug:
No known bugs
-
Todo:

Implement constructor

-

Implement add_particle

-

Implement external_E_field

-

Implement external_B_field

-

Implement force_on_particle

-

Implement total_force_external

-

Implement total_force_particles

-

Implement total_force

-

Implement evolve_RK4

+
Todo:

Implement evolve_RK4

Implement evolve_forward_euler

diff --git a/docs/PenningTrap_8cpp_source.html b/docs/PenningTrap_8cpp_source.html index a4eb76f..0be7c11 100644 --- a/docs/PenningTrap_8cpp_source.html +++ b/docs/PenningTrap_8cpp_source.html @@ -102,69 +102,159 @@ $(document).ready(function(){initNavTree('PenningTrap_8cpp_source.html',''); ini
Go to the documentation of this file.
1
-
24#include "PenningTrap.hpp"
+
16#include "utils.hpp"
+
17#include "PenningTrap.hpp"
+
18#include "constants.hpp"
+
19#include <algorithm>
+
20#include <stdexcept>
+
21#include <omp.h>
+
22
+
23#pragma omp declare reduction( + : arma::vec : omp_out += omp_in ) \
+
24 initializer( omp_priv = omp_orig )
25
26PenningTrap::PenningTrap(double B_0, double V_0, double d)
27{
-
28
-
29}
-
30
-
31void PenningTrap::add_particle(Particle particle)
-
32{
-
33
-
34}
-
35
-
36arma::vec PenningTrap::external_E_field(arma::vec r)
-
37{
-
38
-
39}
-
40
-
41arma::vec PenningTrap::external_B_field(arma::vec r)
-
42{
-
43
-
44}
+
28 this->B_0 = B_0;
+
29 this->V_0 = V_0;
+
30 this->d = d;
+
31}
+
32
+
33void PenningTrap::add_particle(Particle particle)
+
34{
+
35 this->particles.push_back(particle);
+
36}
+
37
+
38arma::vec PenningTrap::external_E_field(arma::vec r)
+
39{
+
40 arma::vec::fixed<3> res;
+
41 double f = this->V_0/(this->d*this->d);
+
42 res(0) = r(0);
+
43 res(1) = r(1);
+
44 res(2) = -2.*r(2);
45
-
46arma::vec PenningTrap::force_on_particle(int i, int j)
-
47{
+
46 return f*res;
+
47}
48
-
49}
-
50
-
51arma::vec PenningTrap::total_force_external(int i)
-
52{
-
53
-
54}
+
49arma::vec PenningTrap::external_B_field(arma::vec r)
+
50{
+
51 arma::vec::fixed<3> res;
+
52 res(0) = 0.;
+
53 res(1) = 0.;
+
54 res(2) = this->B_0;
55
-
56arma::vec PenningTrap::total_force_particles(int i)
-
57{
+
56 return res;
+
57}
58
-
59}
-
60
-
61arma::vec PenningTrap::total_force(int i)
-
62{
-
63
-
64}
-
65
-
66void PenningTrap::evolve_RK4(double dt)
-
67{
-
68
-
69}
+
59arma::vec PenningTrap::force_on_particle(int i, int j)
+
60{
+
61 // Calculate the difference between the particles' position
+
62 arma::vec::fixed<3> res = this->particles.at(i).r_vec
+
63 - this->particles.at(j).r_vec;
+
64
+
65 // Get the distance between the particles
+
66 double norm = arma::norm(res);
+
67
+
68 // Multiply res with p_j's charge divided by the norm cubed
+
69 res *= this->particles.at(j).q/(norm*norm*norm);
70
-
71void PenningTrap::evolve_forward_euler(double dt)
-
72{
+
71 return res;
+
72}
73
-
74}
+
74arma::vec PenningTrap::total_force_external(int i)
+
75{
+
76 Particle p = this->particles.at(i);
+
77
+
78 arma::vec::fixed<3> v_cross_B;
+
79
+
80 arma::vec::fixed<3> B = this->external_B_field(p.r_vec);
+
81
+
82 v_cross_B(0) = p.v_vec(1)*B(2) - p.v_vec(2)*B(1);
+
83 v_cross_B(1) = p.v_vec(2)*B(0) - p.v_vec(0)*B(2);
+
84 v_cross_B(2) = p.v_vec(0)*B(1) - p.v_vec(1)*B(0);
+
85
+
86 arma::vec force = p.q
+
87 *(this->external_E_field(p.r_vec) + v_cross_B);
+
88
+
89 return force;
+
90}
+
91
+
92arma::vec PenningTrap::total_force_particles(int i)
+
93{
+
94 Particle p = this->particles.at(i);
+
95
+
96 arma::vec res(3);
+
97
+
98 for (int j=0; j < this->particles.size(); j++) {
+
99 if (i == j) {
+
100 continue;
+
101 }
+
102
+
103 res += this->force_on_particle(i, j);
+
104 }
+
105
+
106 res *= K_E*(p.q/p.m);
+
107
+
108 return res;
+
109}
+
110
+
111arma::vec PenningTrap::total_force(int i)
+
112{
+
113 return this->total_force_external(i) - this->total_force_particles(i);
+
114}
+
115
+
116void PenningTrap::evolve_RK4(double dt)
+
117{
+
118
+
119}
+
120
+
121void PenningTrap::evolve_forward_euler(double dt)
+
122{
+
123 std::vector<Particle> new_state = this->particles;
+
124
+
125 Particle *p;
+
126
+
127 #pragma omp parallel for private(p)
+
128 for (int i=0; i < this->particles.size(); i++) {
+
129 p = &new_state.at(i);
+
130 p->v_vec += dt*this->total_force(i)/new_state.at(i).m;
+
131 p->r_vec += dt*this->particles.at(i).v_vec;
+
132 }
+
133
+
134 this->particles = new_state;
+
135}
+
136
+
137arma::vec PenningTrap::get_particle(int i)
+
138{
+
139 return this->particles.at(i).r_vec;
+
140}
+
141
+
142double PenningTrap::get_d()
+
143{
+
144 return this->d;
+
145}
PenningTrap.hpp
A class for simulating a Penning trap.
Particle
A class that holds attributes of a particle.
Definition: Particle.hpp:19
-
PenningTrap::total_force_external
arma::vec total_force_external(int i)
Calculate the total external force on a particle.
Definition: PenningTrap.cpp:51
-
PenningTrap::total_force_particles
arma::vec total_force_particles(int i)
Calculate the total force on a particle from other particles.
Definition: PenningTrap.cpp:56
-
PenningTrap::external_B_field
arma::vec external_B_field(arma::vec r)
Calculate B at point r.
Definition: PenningTrap.cpp:41
-
PenningTrap::force_on_particle
arma::vec force_on_particle(int i, int j)
Calculate the force between 2 particles.
Definition: PenningTrap.cpp:46
-
PenningTrap::evolve_forward_euler
void evolve_forward_euler(double dt)
Go forward one timestep using the forward Euler method.
Definition: PenningTrap.cpp:71
-
PenningTrap::add_particle
void add_particle(Particle particle)
Add a particle to the system.
Definition: PenningTrap.cpp:31
+
Particle::v_vec
arma::vec::fixed< 3 > v_vec
velocity
Definition: Particle.hpp:24
+
Particle::q
double q
Charge.
Definition: Particle.hpp:21
+
Particle::r_vec
arma::vec::fixed< 3 > r_vec
position
Definition: Particle.hpp:23
+
Particle::m
double m
Mass.
Definition: Particle.hpp:22
+
PenningTrap::particles
std::vector< Particle > particles
The particles in the Penning trap.
Definition: PenningTrap.hpp:30
+
PenningTrap::total_force_external
arma::vec total_force_external(int i)
Calculate the total external force on a particle.
Definition: PenningTrap.cpp:74
+
PenningTrap::B_0
double B_0
Magnetic field strength.
Definition: PenningTrap.hpp:27
+
PenningTrap::total_force_particles
arma::vec total_force_particles(int i)
Calculate the total force on a particle from other particles.
Definition: PenningTrap.cpp:92
+
PenningTrap::external_B_field
arma::vec external_B_field(arma::vec r)
Calculate B at point r.
Definition: PenningTrap.cpp:49
+
PenningTrap::force_on_particle
arma::vec force_on_particle(int i, int j)
Calculate the force between 2 particles.
Definition: PenningTrap.cpp:59
+
PenningTrap::evolve_forward_euler
void evolve_forward_euler(double dt)
Go forward one timestep using the forward Euler method.
Definition: PenningTrap.cpp:121
+
PenningTrap::d
double d
Characteristic dimension.
Definition: PenningTrap.hpp:29
+
PenningTrap::add_particle
void add_particle(Particle particle)
Add a particle to the system.
Definition: PenningTrap.cpp:33
+
PenningTrap::V_0
double V_0
Applied potential.
Definition: PenningTrap.hpp:28
PenningTrap::PenningTrap
PenningTrap(double B_0=T, double V_0=25.*V/1000., double d=500.)
Set B_0, V_0 and d.
Definition: PenningTrap.cpp:26
-
PenningTrap::total_force
arma::vec total_force(int i)
calculate the total force on a particle.
Definition: PenningTrap.cpp:61
-
PenningTrap::external_E_field
arma::vec external_E_field(arma::vec r)
Calculate E at point r.
Definition: PenningTrap.cpp:36
-
PenningTrap::evolve_RK4
void evolve_RK4(double dt)
Go forward one timestep using the RK4 method.
Definition: PenningTrap.cpp:66
+
PenningTrap::total_force
arma::vec total_force(int i)
calculate the total force on a particle.
Definition: PenningTrap.cpp:111
+
PenningTrap::external_E_field
arma::vec external_E_field(arma::vec r)
Calculate E at point r.
Definition: PenningTrap.cpp:38
+
PenningTrap::evolve_RK4
void evolve_RK4(double dt)
Go forward one timestep using the RK4 method.
Definition: PenningTrap.cpp:116
+
constants.hpp
Library of constants.
+
K_E
#define K_E
Coulomb constant. unit: .
Definition: constants.hpp:15
+
utils.hpp
Function prototypes and macros that are useful.
diff --git a/docs/PenningTrap_8hpp_source.html b/docs/PenningTrap_8hpp_source.html index 6ed7402..b8a0d43 100644 --- a/docs/PenningTrap_8hpp_source.html +++ b/docs/PenningTrap_8hpp_source.html @@ -137,25 +137,29 @@ $(document).ready(function(){initNavTree('PenningTrap_8hpp_source.html',''); ini
73 void evolve_RK4(double dt);
74
77 void evolve_forward_euler(double dt);
-
78};
-
79
-
80#endif
+
78
+
79 arma::vec get_particle(int i);
+
80
+
81 double get_d();
+
82};
+
83
+
84#endif
Particle.hpp
A class that holds the properties of a particle.
Particle
A class that holds attributes of a particle.
Definition: Particle.hpp:19
PenningTrap
A class that simulates a Penning trap.
Definition: PenningTrap.hpp:25
PenningTrap::particles
std::vector< Particle > particles
The particles in the Penning trap.
Definition: PenningTrap.hpp:30
-
PenningTrap::total_force_external
arma::vec total_force_external(int i)
Calculate the total external force on a particle.
Definition: PenningTrap.cpp:51
+
PenningTrap::total_force_external
arma::vec total_force_external(int i)
Calculate the total external force on a particle.
Definition: PenningTrap.cpp:74
PenningTrap::B_0
double B_0
Magnetic field strength.
Definition: PenningTrap.hpp:27
-
PenningTrap::total_force_particles
arma::vec total_force_particles(int i)
Calculate the total force on a particle from other particles.
Definition: PenningTrap.cpp:56
-
PenningTrap::external_B_field
arma::vec external_B_field(arma::vec r)
Calculate B at point r.
Definition: PenningTrap.cpp:41
-
PenningTrap::force_on_particle
arma::vec force_on_particle(int i, int j)
Calculate the force between 2 particles.
Definition: PenningTrap.cpp:46
-
PenningTrap::evolve_forward_euler
void evolve_forward_euler(double dt)
Go forward one timestep using the forward Euler method.
Definition: PenningTrap.cpp:71
+
PenningTrap::total_force_particles
arma::vec total_force_particles(int i)
Calculate the total force on a particle from other particles.
Definition: PenningTrap.cpp:92
+
PenningTrap::external_B_field
arma::vec external_B_field(arma::vec r)
Calculate B at point r.
Definition: PenningTrap.cpp:49
+
PenningTrap::force_on_particle
arma::vec force_on_particle(int i, int j)
Calculate the force between 2 particles.
Definition: PenningTrap.cpp:59
+
PenningTrap::evolve_forward_euler
void evolve_forward_euler(double dt)
Go forward one timestep using the forward Euler method.
Definition: PenningTrap.cpp:121
PenningTrap::d
double d
Characteristic dimension.
Definition: PenningTrap.hpp:29
-
PenningTrap::add_particle
void add_particle(Particle particle)
Add a particle to the system.
Definition: PenningTrap.cpp:31
+
PenningTrap::add_particle
void add_particle(Particle particle)
Add a particle to the system.
Definition: PenningTrap.cpp:33
PenningTrap::V_0
double V_0
Applied potential.
Definition: PenningTrap.hpp:28
-
PenningTrap::total_force
arma::vec total_force(int i)
calculate the total force on a particle.
Definition: PenningTrap.cpp:61
-
PenningTrap::external_E_field
arma::vec external_E_field(arma::vec r)
Calculate E at point r.
Definition: PenningTrap.cpp:36
-
PenningTrap::evolve_RK4
void evolve_RK4(double dt)
Go forward one timestep using the RK4 method.
Definition: PenningTrap.cpp:66
+
PenningTrap::total_force
arma::vec total_force(int i)
calculate the total force on a particle.
Definition: PenningTrap.cpp:111
+
PenningTrap::external_E_field
arma::vec external_E_field(arma::vec r)
Calculate E at point r.
Definition: PenningTrap.cpp:38
+
PenningTrap::evolve_RK4
void evolve_RK4(double dt)
Go forward one timestep using the RK4 method.
Definition: PenningTrap.cpp:116
constants.hpp
Library of constants.
T
#define T
1 Tesla. unit:
Definition: constants.hpp:17
V
#define V
1 Volt. unit:
Definition: constants.hpp:19
diff --git a/docs/classParticle.html b/docs/classParticle.html index 0fcbb7a..997eb37 100644 --- a/docs/classParticle.html +++ b/docs/classParticle.html @@ -185,7 +185,7 @@ Friends

Initialize the particle.

Initialize the particle with a charge, mass, position and velocity.

-

Definition at line 17 of file Particle.cpp.

+

Definition at line 16 of file Particle.cpp.

diff --git a/docs/classPenningTrap-members.html b/docs/classPenningTrap-members.html index dd46fd6..5a26d52 100644 --- a/docs/classPenningTrap-members.html +++ b/docs/classPenningTrap-members.html @@ -112,6 +112,8 @@ $(document).ready(function(){initNavTree('classPenningTrap.html',''); initResiza external_B_field(arma::vec r)PenningTrap external_E_field(arma::vec r)PenningTrap force_on_particle(int i, int j)PenningTrap + get_d() (defined in PenningTrap)PenningTrap + get_particle(int i) (defined in PenningTrap)PenningTrap particlesPenningTrapprivate PenningTrap(double B_0=T, double V_0=25.*V/1000., double d=500.)PenningTrap total_force(int i)PenningTrap diff --git a/docs/classPenningTrap.html b/docs/classPenningTrap.html index fcbbf79..7809dc2 100644 --- a/docs/classPenningTrap.html +++ b/docs/classPenningTrap.html @@ -143,6 +143,10 @@ Public Member Functions void evolve_forward_euler (double dt)  Go forward one timestep using the forward Euler method.
  +arma::vec get_particle (int i) +  +double get_d () +  @@ -222,7 +226,7 @@ Private Attributes

Add a particle to the system.

-

Definition at line 31 of file PenningTrap.cpp.

+

Definition at line 33 of file PenningTrap.cpp.

@@ -244,7 +248,7 @@ Private Attributes

Go forward one timestep using the forward Euler method.

-

Definition at line 71 of file PenningTrap.cpp.

+

Definition at line 121 of file PenningTrap.cpp.

@@ -266,7 +270,7 @@ Private Attributes

Go forward one timestep using the RK4 method.

-

Definition at line 66 of file PenningTrap.cpp.

+

Definition at line 116 of file PenningTrap.cpp.

@@ -288,7 +292,7 @@ Private Attributes

Calculate B at point r.

-

Definition at line 41 of file PenningTrap.cpp.

+

Definition at line 49 of file PenningTrap.cpp.

@@ -310,7 +314,7 @@ Private Attributes

Calculate E at point r.

-

Definition at line 36 of file PenningTrap.cpp.

+

Definition at line 38 of file PenningTrap.cpp.

@@ -343,7 +347,46 @@ Private Attributes

Calculate the force between 2 particles.

Calculate the force exhibited on particle p_i from particle p_j.

-

Definition at line 46 of file PenningTrap.cpp.

+

Definition at line 59 of file PenningTrap.cpp.

+ + + + +

◆ get_d()

+ +
+
+

Private Attributes

+ + + + + + +
double PenningTrap::get_d ()
+
+ +

Definition at line 142 of file PenningTrap.cpp.

+ +
+ + +

◆ get_particle()

+ +
+
+ + + + + + + + +
arma::vec PenningTrap::get_particle (int i)
+
+ +

Definition at line 137 of file PenningTrap.cpp.

@@ -365,7 +408,7 @@ Private Attributes

calculate the total force on a particle.

-

Definition at line 61 of file PenningTrap.cpp.

+

Definition at line 111 of file PenningTrap.cpp.

@@ -388,7 +431,7 @@ Private Attributes

Calculate the total external force on a particle.

Calculate the total amount of force that E and B exhibits on particle p_i.

-

Definition at line 51 of file PenningTrap.cpp.

+

Definition at line 74 of file PenningTrap.cpp.

@@ -410,7 +453,7 @@ Private Attributes

Calculate the total force on a particle from other particles.

-

Definition at line 56 of file PenningTrap.cpp.

+

Definition at line 92 of file PenningTrap.cpp.

diff --git a/docs/dir_68267d1309a1af8e8297ef4c3efbcdba.html b/docs/dir_68267d1309a1af8e8297ef4c3efbcdba.html index fe23dd4..2f63e67 100644 --- a/docs/dir_68267d1309a1af8e8297ef4c3efbcdba.html +++ b/docs/dir_68267d1309a1af8e8297ef4c3efbcdba.html @@ -113,6 +113,8 @@ Files file  PenningTrap.cpp [code]  The implementation of the PenningTrap class.
  +file  test.py [code] +  file  test_suite.cpp [code]  The test suite for the project.
  diff --git a/docs/dir_68267d1309a1af8e8297ef4c3efbcdba.js b/docs/dir_68267d1309a1af8e8297ef4c3efbcdba.js index 17c90f6..f72184d 100644 --- a/docs/dir_68267d1309a1af8e8297ef4c3efbcdba.js +++ b/docs/dir_68267d1309a1af8e8297ef4c3efbcdba.js @@ -3,6 +3,7 @@ var dir_68267d1309a1af8e8297ef4c3efbcdba = [ "main.cpp", "main_8cpp.html", null ], [ "Particle.cpp", "Particle_8cpp.html", null ], [ "PenningTrap.cpp", "PenningTrap_8cpp.html", null ], + [ "test.py", "test_8py_source.html", null ], [ "test_suite.cpp", "test__suite_8cpp.html", null ], [ "utils.cpp", "utils_8cpp.html", "utils_8cpp" ] ]; \ No newline at end of file diff --git a/docs/files.html b/docs/files.html index a998ae5..f3ab17f 100644 --- a/docs/files.html +++ b/docs/files.html @@ -112,8 +112,9 @@ $(document).ready(function(){initNavTree('files.html',''); initResizable(); });  main.cppThe main program for this project  Particle.cppThe implementation of the Particle class  PenningTrap.cppThe implementation of the PenningTrap class - test_suite.cppThe test suite for the project - utils.cppImplementation of the utils + test.py + test_suite.cppThe test suite for the project + utils.cppImplementation of the utils diff --git a/docs/main_8cpp.html b/docs/main_8cpp.html index cf16990..d86dfaa 100644 --- a/docs/main_8cpp.html +++ b/docs/main_8cpp.html @@ -99,6 +99,7 @@ $(document).ready(function(){initNavTree('main_8cpp.html',''); initResizable();
+Macros | Functions
main.cpp File Reference
@@ -106,11 +107,30 @@ $(document).ready(function(){initNavTree('main_8cpp.html',''); initResizable();

The main program for this project. More...

- +
#include <filesystem>
+#include <fstream>
+#include <string>
+#include <omp.h>
+#include <sys/stat.h>
+#include "PenningTrap.hpp"
+

Go to the source code of this file.

+ + + + + + + + + +

+Macros

#define PARTICLES   100
 
#define N   10000
 
#define CHARGE   1.
 
#define MASS   40.
 
+ +

Functions

void euler_100_particles ()
 
int main ()
 
@@ -123,7 +143,91 @@ Janita Ovidie Sandtrøen Willumsen (janitaws)
Bug:
No known bugs

Definition in file main.cpp.

-

Function Documentation

+

Macro Definition Documentation

+ +

◆ CHARGE

+ +
+
+ + + + +
#define CHARGE   1.
+
+ +

Definition at line 23 of file main.cpp.

+ +
+
+ +

◆ MASS

+ +
+
+ + + + +
#define MASS   40.
+
+ +

Definition at line 24 of file main.cpp.

+ +
+
+ +

◆ N

+ +
+
+ + + + +
#define N   10000
+
+ +

Definition at line 22 of file main.cpp.

+ +
+
+ +

◆ PARTICLES

+ +
+
+ + + + +
#define PARTICLES   100
+
+ +

Definition at line 21 of file main.cpp.

+ +
+
+

Function Documentation

+ +

◆ euler_100_particles()

+ +
+
+ + + + + + + +
void euler_100_particles ()
+
+ +

Definition at line 26 of file main.cpp.

+ +
+

◆ main()

@@ -139,7 +243,7 @@ Janita Ovidie Sandtrøen Willumsen (janitaws)
-

Definition at line 13 of file main.cpp.

+

Definition at line 77 of file main.cpp.

diff --git a/docs/main_8cpp_source.html b/docs/main_8cpp_source.html index bf66ff8..3cc7c06 100644 --- a/docs/main_8cpp_source.html +++ b/docs/main_8cpp_source.html @@ -102,10 +102,87 @@ $(document).ready(function(){initNavTree('main_8cpp_source.html',''); initResiza
Go to the documentation of this file.
1
-
13int main()
-
14{
-
15 return 0;
-
16}
+
13#include <filesystem>
+
14#include <fstream>
+
15#include <string>
+
16#include <omp.h>
+
17#include <sys/stat.h>
+
18
+
19#include "PenningTrap.hpp"
+
20
+
21#define PARTICLES 100
+
22#define N 10000
+
23#define CHARGE 1.
+
24#define MASS 40. // unit: amu
+
25
+
26void euler_100_particles()
+
27{
+
28 PenningTrap trap;
+
29
+
30 // Add particles inside trap
+
31 for (int i=0; i < PARTICLES; i++) {
+
32 arma::vec r = arma::vec(3).randn() * 0.1 * trap.get_d(); // random initial position
+
33 arma::vec v = arma::vec(3).randn() * 0.1 * trap.get_d(); // random initial velocity
+
34 trap.add_particle(Particle(CHARGE, MASS, r, v));
+
35 }
+
36
+
37 double time = 50.; // microseconds
+
38 double dt = time / (double) N;
+
39
+
40 auto res = new arma::vec::fixed<3>[PARTICLES][N];
+
41
+
42 int counter = 0;
+
43
+
44 // Get the path of all particles
+
45 for (int j=0; j < N; j++) {
+
46 #pragma omp parallel for
+
47 for (int i=0; i < PARTICLES; i++) {
+
48 res[i][j] = trap.get_particle(i);
+
49 }
+
50 trap.evolve_forward_euler(dt);
+
51 }
+
52
+
53 std::cout << counter << std::endl;
+
54
+
55 arma::vec::fixed<3>* cur_row;
+
56 arma::vec::fixed<3> cur_elem;
+
57
+
58 mkdir("output", 0777);
+
59
+
60 std::ofstream ofile;
+
61
+
62 // Write particle paths to file
+
63 #pragma omp parallel for private(cur_row, cur_elem, ofile)
+
64 for (int i=0; i < PARTICLES; i++) {
+
65 cur_row = res[i];
+
66 ofile.open("output/p" + std::to_string(i) + ".txt");
+
67 for (int j=0; j < N; j++) {
+
68 cur_elem = cur_row[j];
+
69 ofile << cur_elem(0) << ","
+
70 << cur_elem(1) << ","
+
71 << cur_elem(2) << "\n";
+
72 }
+
73 ofile.close();
+
74 }
+
75}
+
76
+
77int main()
+
78{
+
79 double start = omp_get_wtime();
+
80
+
81 euler_100_particles();
+
82
+
83 double end = omp_get_wtime();
+
84
+
85 std::cout << "Time: " << end - start << " seconds" << std::endl;
+
86
+
87 return 0;
+
88}
+
PenningTrap.hpp
A class for simulating a Penning trap.
+
Particle
A class that holds attributes of a particle.
Definition: Particle.hpp:19
+
PenningTrap
A class that simulates a Penning trap.
Definition: PenningTrap.hpp:25
+
PenningTrap::evolve_forward_euler
void evolve_forward_euler(double dt)
Go forward one timestep using the forward Euler method.
Definition: PenningTrap.cpp:121
+
PenningTrap::add_particle
void add_particle(Particle particle)
Add a particle to the system.
Definition: PenningTrap.cpp:33
diff --git a/docs/navtreeindex0.js b/docs/navtreeindex0.js index 6d8de5d..75ec38f 100644 --- a/docs/navtreeindex0.js +++ b/docs/navtreeindex0.js @@ -53,13 +53,14 @@ var NAVTREEINDEX0 = "main_8cpp.html":[4,0,1,0], "main_8cpp_source.html":[4,0,1,0], "pages.html":[], -"test__suite_8cpp.html":[4,0,1,3], -"test__suite_8cpp_source.html":[4,0,1,3], +"test_8py_source.html":[4,0,1,3], +"test__suite_8cpp.html":[4,0,1,4], +"test__suite_8cpp_source.html":[4,0,1,4], "todo.html":[2], -"utils_8cpp.html":[4,0,1,4], -"utils_8cpp.html#a58565270b643b24e3132f38c653e0199":[4,0,1,4,0], -"utils_8cpp.html#acd2a9c7a7d5a7fe9163be8c4cc110746":[4,0,1,4,1], -"utils_8cpp_source.html":[4,0,1,4], +"utils_8cpp.html":[4,0,1,5], +"utils_8cpp.html#a58565270b643b24e3132f38c653e0199":[4,0,1,5,0], +"utils_8cpp.html#acd2a9c7a7d5a7fe9163be8c4cc110746":[4,0,1,5,1], +"utils_8cpp_source.html":[4,0,1,5], "utils_8hpp.html":[4,0,0,3], "utils_8hpp.html#ad54b96a1074f9df4dc892a41d115b72d":[4,0,0,3,1], "utils_8hpp.html#adfb618b2fdff47ef30a4a2b62c04f384":[4,0,0,3,2], diff --git a/docs/test_8py_source.html b/docs/test_8py_source.html new file mode 100644 index 0000000..17942db --- /dev/null +++ b/docs/test_8py_source.html @@ -0,0 +1,177 @@ + + + + + + + +Penning Trap Simulation: src/test.py Source File + + + + + + + + + + + + + + + + +
+
+ + + + + + +
+
Penning Trap Simulation +
+
Simulate particle behavior inside a Penning Trap
+
+
+ + + + + + + +
+
+ +
+
+
+ +
+ +
+
+ + +
+
+
+
+
+
Loading...
+
Searching...
+
No Matches
+
+
+
+
+ +
+
test.py
+
+
+
1import matplotlib.pyplot as plt
+
2import numpy as np
+
3from mpl_toolkits.mplot3d import Axes3D
+
4from matplotlib import animation
+
5
+
6def get_data(files):
+
7 res = []
+
8 for file in files:
+
9 arr = [[], [], []]
+
10 with open(file) as f:
+
11 lines = f.readlines()
+
12
+
13 for line in lines:
+
14 xi,yi,zi = map(lambda x: float(x), line.strip().split(","))
+
15 arr[0].append(xi)
+
16 arr[1].append(yi)
+
17 arr[2].append(zi)
+
18 res.append(arr)
+
19
+
20 return np.array(res)
+
21
+
22def update(num, lines, arr):
+
23 for line, a in zip(lines, arr):
+
24 line.set_data(a[:2, num])
+
25 line.set_3d_properties(a[2, num])
+
26
+
27
+
28
+
29def animate():
+
30 fig = plt.figure()
+
31 ax = fig.add_subplot(projection="3d")
+
32
+
33
+
34 arr = get_data([f"output/p{i}.txt" for i in range(100)])
+
35
+
36 arr = arr[:,:,::10]
+
37
+
38 N = len(arr[0][0])
+
39
+
40 lines = [ax.plot(*a[:,1], "o")[0] for a in arr]
+
41
+
42 ax.set_xlim3d([-500.0, 500.0])
+
43 ax.set_xlabel('X')
+
44
+
45 ax.set_ylim3d([-500.0, 500.0])
+
46 ax.set_ylabel('Y')
+
47
+
48 ax.set_zlim3d([-500.0, 500.0])
+
49 ax.set_zlabel('Z')
+
50
+
51 ani = animation.FuncAnimation(fig, update, N, fargs=(lines, arr),
+
52 interval=1,
+
53 blit=False)
+
54
+
55
+
56 # ani.save("100_particles.gif", writer=animation.FFMpegFileWriter(fps=30))
+
57 plt.show()
+
58
+
59
+
60if __name__ == "__main__":
+
61 animate()
+
62
+
63
+
+
+ + + + diff --git a/docs/todo.html b/docs/todo.html index 8c5247d..b6363dc 100644 --- a/docs/todo.html +++ b/docs/todo.html @@ -102,26 +102,8 @@ $(document).ready(function(){initNavTree('todo.html',''); initResizable(); });
-
File Particle.cpp
-
Implement constructor
File PenningTrap.cpp
-

Implement constructor

-

-

Implement add_particle

-

-

Implement external_E_field

-

-

Implement external_B_field

-

-

Implement force_on_particle

-

-

Implement total_force_external

-

-

Implement total_force_particles

-

-

Implement total_force

-

-

Implement evolve_RK4

+

Implement evolve_RK4

Implement evolve_forward_euler

diff --git a/include/PenningTrap.hpp b/include/PenningTrap.hpp index b0dbd86..1b1752e 100644 --- a/include/PenningTrap.hpp +++ b/include/PenningTrap.hpp @@ -75,6 +75,10 @@ public: /** @brief Go forward one timestep using the forward Euler method * */ void evolve_forward_euler(double dt); + + arma::vec get_particle(int i); + + double get_d(); }; #endif diff --git a/man_pages/man3/Particle.3 b/man_pages/man3/Particle.3 index 466ecc4..db9056a 100644 --- a/man_pages/man3/Particle.3 +++ b/man_pages/man3/Particle.3 @@ -1,4 +1,4 @@ -.TH "Particle" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "Particle" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME @@ -57,7 +57,7 @@ Definition at line \fB19\fP of file \fBParticle\&.hpp\fP\&. .PP Initialize the particle\&. Initialize the particle with a charge, mass, position and velocity\&. .PP -Definition at line \fB17\fP of file \fBParticle\&.cpp\fP\&. +Definition at line \fB16\fP of file \fBParticle\&.cpp\fP\&. .SH "Friends And Related Function Documentation" .PP .SS "friend class \fBPenningTrap\fP\fC [friend]\fP" diff --git a/man_pages/man3/Particle.cpp.3 b/man_pages/man3/Particle.cpp.3 index 1d94e76..2799e38 100644 --- a/man_pages/man3/Particle.cpp.3 +++ b/man_pages/man3/Particle.cpp.3 @@ -1,4 +1,4 @@ -.TH "src/Particle.cpp" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "src/Particle.cpp" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME @@ -34,12 +34,6 @@ Janita Ovidie Sandtrøen Willumsen (janitaws) No known bugs .RE .PP -.PP -\fBTodo\fP -.RS 4 -Implement constructor -.RE -.PP .PP Definition in file \fBParticle\&.cpp\fP\&. diff --git a/man_pages/man3/Particle.hpp.3 b/man_pages/man3/Particle.hpp.3 index 7cec015..2a5a85e 100644 --- a/man_pages/man3/Particle.hpp.3 +++ b/man_pages/man3/Particle.hpp.3 @@ -1,4 +1,4 @@ -.TH "include/Particle.hpp" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "include/Particle.hpp" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME diff --git a/man_pages/man3/PenningTrap.3 b/man_pages/man3/PenningTrap.3 index 503511e..42c26cf 100644 --- a/man_pages/man3/PenningTrap.3 +++ b/man_pages/man3/PenningTrap.3 @@ -1,4 +1,4 @@ -.TH "PenningTrap" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "PenningTrap" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME @@ -52,6 +52,12 @@ PenningTrap \- A class that simulates a Penning trap\&. .RI "void \fBevolve_forward_euler\fP (double dt)" .br .RI "Go forward one timestep using the forward Euler method\&. " +.ti -1c +.RI "arma::vec \fBget_particle\fP (int i)" +.br +.ti -1c +.RI "double \fBget_d\fP ()" +.br .in -1c .SS "Private Attributes" @@ -95,55 +101,63 @@ Definition at line \fB26\fP of file \fBPenningTrap\&.cpp\fP\&. .PP Add a particle to the system\&. .PP -Definition at line \fB31\fP of file \fBPenningTrap\&.cpp\fP\&. +Definition at line \fB33\fP of file \fBPenningTrap\&.cpp\fP\&. .SS "void PenningTrap::evolve_forward_euler (double dt)" .PP Go forward one timestep using the forward Euler method\&. .PP -Definition at line \fB71\fP of file \fBPenningTrap\&.cpp\fP\&. +Definition at line \fB121\fP of file \fBPenningTrap\&.cpp\fP\&. .SS "void PenningTrap::evolve_RK4 (double dt)" .PP Go forward one timestep using the RK4 method\&. .PP -Definition at line \fB66\fP of file \fBPenningTrap\&.cpp\fP\&. +Definition at line \fB116\fP of file \fBPenningTrap\&.cpp\fP\&. .SS "arma::vec PenningTrap::external_B_field (arma::vec r)" .PP Calculate B at point r\&. .PP -Definition at line \fB41\fP of file \fBPenningTrap\&.cpp\fP\&. +Definition at line \fB49\fP of file \fBPenningTrap\&.cpp\fP\&. .SS "arma::vec PenningTrap::external_E_field (arma::vec r)" .PP Calculate E at point r\&. .PP -Definition at line \fB36\fP of file \fBPenningTrap\&.cpp\fP\&. +Definition at line \fB38\fP of file \fBPenningTrap\&.cpp\fP\&. .SS "arma::vec PenningTrap::force_on_particle (int i, int j)" .PP Calculate the force between 2 particles\&. Calculate the force exhibited on particle p_i from particle p_j\&. .PP -Definition at line \fB46\fP of file \fBPenningTrap\&.cpp\fP\&. +Definition at line \fB59\fP of file \fBPenningTrap\&.cpp\fP\&. +.SS "double PenningTrap::get_d ()" + +.PP +Definition at line \fB142\fP of file \fBPenningTrap\&.cpp\fP\&. +.SS "arma::vec PenningTrap::get_particle (int i)" + +.PP +Definition at line \fB137\fP of file \fBPenningTrap\&.cpp\fP\&. .SS "arma::vec PenningTrap::total_force (int i)" .PP calculate the total force on a particle\&. .PP -Definition at line \fB61\fP of file \fBPenningTrap\&.cpp\fP\&. +Definition at line \fB111\fP of file \fBPenningTrap\&.cpp\fP\&. .SS "arma::vec PenningTrap::total_force_external (int i)" .PP Calculate the total external force on a particle\&. Calculate the total amount of force that E and B exhibits on particle p_i\&. .PP -Definition at line \fB51\fP of file \fBPenningTrap\&.cpp\fP\&. +Definition at line \fB74\fP of file \fBPenningTrap\&.cpp\fP\&. .SS "arma::vec PenningTrap::total_force_particles (int i)" .PP Calculate the total force on a particle from other particles\&. .PP -Definition at line \fB56\fP of file \fBPenningTrap\&.cpp\fP\&. +Definition at line \fB92\fP of file \fBPenningTrap\&.cpp\fP\&. .SH "Member Data Documentation" .PP .SS "double PenningTrap::B_0\fC [private]\fP" diff --git a/man_pages/man3/PenningTrap.cpp.3 b/man_pages/man3/PenningTrap.cpp.3 index 68aaa57..fecf53b 100644 --- a/man_pages/man3/PenningTrap.cpp.3 +++ b/man_pages/man3/PenningTrap.cpp.3 @@ -1,4 +1,4 @@ -.TH "src/PenningTrap.cpp" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "src/PenningTrap.cpp" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME @@ -7,8 +7,18 @@ src/PenningTrap.cpp \- The implementation of the \fBPenningTrap\fP class\&. .SH SYNOPSIS .br .PP +\fC#include 'utils\&.hpp'\fP +.br \fC#include 'PenningTrap\&.hpp'\fP .br +\fC#include 'constants\&.hpp'\fP +.br +\fC#include \fP +.br +\fC#include \fP +.br +\fC#include \fP +.br .SH "Detailed Description" .PP @@ -37,22 +47,6 @@ No known bugs .PP \fBTodo\fP .RS 4 -Implement constructor -.PP -Implement add_particle -.PP -Implement external_E_field -.PP -Implement external_B_field -.PP -Implement force_on_particle -.PP -Implement total_force_external -.PP -Implement total_force_particles -.PP -Implement total_force -.PP Implement evolve_RK4 .PP Implement evolve_forward_euler diff --git a/man_pages/man3/PenningTrap.hpp.3 b/man_pages/man3/PenningTrap.hpp.3 index 2857a95..4d434f2 100644 --- a/man_pages/man3/PenningTrap.hpp.3 +++ b/man_pages/man3/PenningTrap.hpp.3 @@ -1,4 +1,4 @@ -.TH "include/PenningTrap.hpp" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "include/PenningTrap.hpp" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME diff --git a/man_pages/man3/bug.3 b/man_pages/man3/bug.3 index 670c28b..0362f46 100644 --- a/man_pages/man3/bug.3 +++ b/man_pages/man3/bug.3 @@ -1,4 +1,4 @@ -.TH "bug" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "bug" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME diff --git a/man_pages/man3/constants.hpp.3 b/man_pages/man3/constants.hpp.3 index 56030e6..2002a1b 100644 --- a/man_pages/man3/constants.hpp.3 +++ b/man_pages/man3/constants.hpp.3 @@ -1,4 +1,4 @@ -.TH "include/constants.hpp" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "include/constants.hpp" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME diff --git a/man_pages/man3/main.cpp.3 b/man_pages/man3/main.cpp.3 index 294b2fa..4cab7a6 100644 --- a/man_pages/man3/main.cpp.3 +++ b/man_pages/man3/main.cpp.3 @@ -1,4 +1,4 @@ -.TH "src/main.cpp" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "src/main.cpp" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME @@ -7,10 +7,42 @@ src/main.cpp \- The main program for this project\&. .SH SYNOPSIS .br .PP +\fC#include \fP +.br +\fC#include \fP +.br +\fC#include \fP +.br +\fC#include \fP +.br +\fC#include \fP +.br +\fC#include 'PenningTrap\&.hpp'\fP +.br + +.SS "Macros" + +.in +1c +.ti -1c +.RI "#define \fBPARTICLES\fP 100" +.br +.ti -1c +.RI "#define \fBN\fP 10000" +.br +.ti -1c +.RI "#define \fBCHARGE\fP 1\&." +.br +.ti -1c +.RI "#define \fBMASS\fP 40\&." +.br +.in -1c .SS "Functions" .in +1c .ti -1c +.RI "void \fBeuler_100_particles\fP ()" +.br +.ti -1c .RI "int \fBmain\fP ()" .br .in -1c @@ -41,12 +73,34 @@ No known bugs .PP Definition in file \fBmain\&.cpp\fP\&. +.SH "Macro Definition Documentation" +.PP +.SS "#define CHARGE 1\&." + +.PP +Definition at line \fB23\fP of file \fBmain\&.cpp\fP\&. +.SS "#define MASS 40\&." + +.PP +Definition at line \fB24\fP of file \fBmain\&.cpp\fP\&. +.SS "#define N 10000" + +.PP +Definition at line \fB22\fP of file \fBmain\&.cpp\fP\&. +.SS "#define PARTICLES 100" + +.PP +Definition at line \fB21\fP of file \fBmain\&.cpp\fP\&. .SH "Function Documentation" .PP +.SS "void euler_100_particles ()" + +.PP +Definition at line \fB26\fP of file \fBmain\&.cpp\fP\&. .SS "int main ()" .PP -Definition at line \fB13\fP of file \fBmain\&.cpp\fP\&. +Definition at line \fB77\fP of file \fBmain\&.cpp\fP\&. .SH "Author" .PP Generated automatically by Doxygen for Penning Trap Simulation from the source code\&. diff --git a/man_pages/man3/test_suite.cpp.3 b/man_pages/man3/test_suite.cpp.3 index c7bdeea..35cfa4b 100644 --- a/man_pages/man3/test_suite.cpp.3 +++ b/man_pages/man3/test_suite.cpp.3 @@ -1,4 +1,4 @@ -.TH "src/test_suite.cpp" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "src/test_suite.cpp" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME diff --git a/man_pages/man3/todo.3 b/man_pages/man3/todo.3 index e48ea02..c342666 100644 --- a/man_pages/man3/todo.3 +++ b/man_pages/man3/todo.3 @@ -1,37 +1,11 @@ -.TH "todo" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "todo" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME todo \- Todo List .PP -.IP "\fBFile \fBParticle\&.cpp\fP \fP" 1c -Implement constructor .IP "\fBFile \fBPenningTrap\&.cpp\fP \fP" 1c -Implement constructor -.PP -.PP -Implement add_particle -.PP -.PP -Implement external_E_field -.PP -.PP -Implement external_B_field -.PP -.PP -Implement force_on_particle -.PP -.PP -Implement total_force_external -.PP -.PP -Implement total_force_particles -.PP -.PP -Implement total_force -.PP -.PP Implement evolve_RK4 .PP .PP diff --git a/man_pages/man3/utils.cpp.3 b/man_pages/man3/utils.cpp.3 index 1660e32..910bd0a 100644 --- a/man_pages/man3/utils.cpp.3 +++ b/man_pages/man3/utils.cpp.3 @@ -1,4 +1,4 @@ -.TH "src/utils.cpp" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "src/utils.cpp" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME diff --git a/man_pages/man3/utils.hpp.3 b/man_pages/man3/utils.hpp.3 index 3aba776..76556bd 100644 --- a/man_pages/man3/utils.hpp.3 +++ b/man_pages/man3/utils.hpp.3 @@ -1,4 +1,4 @@ -.TH "include/utils.hpp" 3 "Thu Sep 28 2023" "Penning Trap Simulation" \" -*- nroff -*- +.TH "include/utils.hpp" 3 "Mon Oct 2 2023" "Penning Trap Simulation" \" -*- nroff -*- .ad l .nh .SH NAME diff --git a/src/Makefile b/src/Makefile index fd101c7..f0cd401 100644 --- a/src/Makefile +++ b/src/Makefile @@ -8,7 +8,7 @@ CLASSOBJS=$(CLASSSRCS:.cpp=.o) INCLUDE=../include -CFLAGS=-larmadillo -llapack -std=c++11 +CFLAGS=-larmadillo -llapack -std=c++11 -O2 OPENMP=-fopenmp # Add a debug flag when compiling (For the DEBUG macro in utils.hpp) @@ -24,7 +24,7 @@ endif all: test_suite main # Rules for executables -main: main.o $(LIBOBJS) +main: main.o $(LIBOBJS) $(CLASSOBJS) $(CC) $^ -o $@ $(CFLAGS) $(DBGFLAG) -I$(INCLUDE) $(OPENMP) test_suite: test_suite.o $(LIBOBJS) @@ -32,7 +32,7 @@ test_suite: test_suite.o $(LIBOBJS) # Rule for object files %.o: %.cpp - $(CC) -c $^ -o $@ $(DBGFLAG) -I$(INCLUDE) $(OPENMP) + $(CC) -c $^ -o $@ $(CFLAGS) $(DBGFLAG) -I$(INCLUDE) $(OPENMP) clean: rm *.o diff --git a/src/Particle.cpp b/src/Particle.cpp index 1bcdb15..0ff18e8 100644 --- a/src/Particle.cpp +++ b/src/Particle.cpp @@ -9,7 +9,6 @@ * * @bug No known bugs * - * @todo Implement constructor * */ #include "Particle.hpp" @@ -18,5 +17,9 @@ Particle::Particle(double q, double m, arma::vec::fixed<3> r_vec, arma::vec::fixed<3> v_vec) { - + // Giving the particle its properties + this->q = q; + this->m = m; + this->r_vec = r_vec; + this->v_vec = v_vec; } diff --git a/src/PenningTrap.cpp b/src/PenningTrap.cpp index 7a4ff2b..59be75c 100644 --- a/src/PenningTrap.cpp +++ b/src/PenningTrap.cpp @@ -9,58 +9,108 @@ * * @bug No known bugs * - * @todo Implement constructor - * @todo Implement add_particle - * @todo Implement external_E_field - * @todo Implement external_B_field - * @todo Implement force_on_particle - * @todo Implement total_force_external - * @todo Implement total_force_particles - * @todo Implement total_force * @todo Implement evolve_RK4 * @todo Implement evolve_forward_euler * */ +#include "utils.hpp" #include "PenningTrap.hpp" +#include "constants.hpp" +#include +#include +#include + +#pragma omp declare reduction( + : arma::vec : omp_out += omp_in ) \ + initializer( omp_priv = omp_orig ) PenningTrap::PenningTrap(double B_0, double V_0, double d) { - + this->B_0 = B_0; + this->V_0 = V_0; + this->d = d; } void PenningTrap::add_particle(Particle particle) { - + this->particles.push_back(particle); } arma::vec PenningTrap::external_E_field(arma::vec r) { + arma::vec::fixed<3> res; + double f = this->V_0/(this->d*this->d); + res(0) = r(0); + res(1) = r(1); + res(2) = -2.*r(2); + return f*res; } arma::vec PenningTrap::external_B_field(arma::vec r) { + arma::vec::fixed<3> res; + res(0) = 0.; + res(1) = 0.; + res(2) = this->B_0; + return res; } arma::vec PenningTrap::force_on_particle(int i, int j) { + // Calculate the difference between the particles' position + arma::vec::fixed<3> res = this->particles.at(i).r_vec + - this->particles.at(j).r_vec; + // Get the distance between the particles + double norm = arma::norm(res); + + // Multiply res with p_j's charge divided by the norm cubed + res *= this->particles.at(j).q/(norm*norm*norm); + + return res; } arma::vec PenningTrap::total_force_external(int i) { + Particle p = this->particles.at(i); + arma::vec::fixed<3> v_cross_B; + + arma::vec::fixed<3> B = this->external_B_field(p.r_vec); + + v_cross_B(0) = p.v_vec(1)*B(2) - p.v_vec(2)*B(1); + v_cross_B(1) = p.v_vec(2)*B(0) - p.v_vec(0)*B(2); + v_cross_B(2) = p.v_vec(0)*B(1) - p.v_vec(1)*B(0); + + arma::vec force = p.q + *(this->external_E_field(p.r_vec) + v_cross_B); + + return force; } arma::vec PenningTrap::total_force_particles(int i) { + Particle p = this->particles.at(i); + arma::vec res(3); + + for (int j=0; j < this->particles.size(); j++) { + if (i == j) { + continue; + } + + res += this->force_on_particle(i, j); + } + + res *= K_E*(p.q/p.m); + + return res; } arma::vec PenningTrap::total_force(int i) { - + return this->total_force_external(i) - this->total_force_particles(i); } void PenningTrap::evolve_RK4(double dt) @@ -70,5 +120,26 @@ void PenningTrap::evolve_RK4(double dt) void PenningTrap::evolve_forward_euler(double dt) { + std::vector new_state = this->particles; + Particle *p; + + #pragma omp parallel for private(p) + for (int i=0; i < this->particles.size(); i++) { + p = &new_state.at(i); + p->v_vec += dt*this->total_force(i)/new_state.at(i).m; + p->r_vec += dt*this->particles.at(i).v_vec; + } + + this->particles = new_state; +} + +arma::vec PenningTrap::get_particle(int i) +{ + return this->particles.at(i).r_vec; +} + +double PenningTrap::get_d() +{ + return this->d; } diff --git a/src/main.cpp b/src/main.cpp index 4c46902..b479649 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -10,7 +10,79 @@ * @bug No known bugs * */ -int main() +#include +#include +#include +#include +#include + +#include "PenningTrap.hpp" + +#define PARTICLES 100 +#define N 10000 +#define CHARGE 1. +#define MASS 40. // unit: amu + +void euler_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)); + } + + double time = 50.; // microseconds + double dt = time / (double) N; + + auto res = new arma::vec::fixed<3>[PARTICLES][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 < PARTICLES; i++) { + res[i][j] = trap.get_particle(i); + } + trap.evolve_forward_euler(dt); + } + + std::cout << counter << std::endl; + + arma::vec::fixed<3>* cur_row; + arma::vec::fixed<3> cur_elem; + + mkdir("output", 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/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(); + } +} + +int main() +{ + double start = omp_get_wtime(); + + euler_100_particles(); + + double end = omp_get_wtime(); + + std::cout << "Time: " << end - start << " seconds" << std::endl; + return 0; } diff --git a/src/test.py b/src/test.py new file mode 100644 index 0000000..a215cad --- /dev/null +++ b/src/test.py @@ -0,0 +1,63 @@ +import matplotlib.pyplot as plt +import numpy as np +from mpl_toolkits.mplot3d import Axes3D +from matplotlib import animation + +def get_data(files): + res = [] + for file in files: + arr = [[], [], []] + with open(file) as f: + lines = f.readlines() + + for line in lines: + xi,yi,zi = map(lambda x: float(x), line.strip().split(",")) + arr[0].append(xi) + arr[1].append(yi) + arr[2].append(zi) + res.append(arr) + + return np.array(res) + +def update(num, lines, arr): + for line, a in zip(lines, arr): + line.set_data(a[:2, num]) + line.set_3d_properties(a[2, num]) + + + +def animate(): + fig = plt.figure() + ax = fig.add_subplot(projection="3d") + + + arr = get_data([f"output/p{i}.txt" for i in range(100)]) + + arr = arr[:,:,::10] + + N = len(arr[0][0]) + + lines = [ax.plot(*a[:,1], "o")[0] for a in arr] + + ax.set_xlim3d([-500.0, 500.0]) + ax.set_xlabel('X') + + ax.set_ylim3d([-500.0, 500.0]) + ax.set_ylabel('Y') + + ax.set_zlim3d([-500.0, 500.0]) + ax.set_zlabel('Z') + + ani = animation.FuncAnimation(fig, update, N, fargs=(lines, arr), + interval=1, + blit=False) + + + # ani.save("100_particles.gif", writer=animation.FFMpegFileWriter(fps=30)) + plt.show() + + +if __name__ == "__main__": + animate() + +