-
Notifications
You must be signed in to change notification settings - Fork 0
/
main.cpp
121 lines (104 loc) · 3.35 KB
/
main.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
#include <omp.h>
#include <iostream>
#include <fstream>
#include "parameter.h"
#include "particle.h"
#include "tree.h"
#include "function.h"
using namespace std;
int main()
{
// Initialize
double Size = 1000.0, // Size of the simulation [kpc]
Origin[3] = {0 - Size / 2, 0 - Size / 2, 0 - Size / 2};
// Read particles positions and velocities from file
ifstream file("halo.txt");
cout << "Reading file..."
<< "\n";
// Create a folder to store the output
ifstream folder("out");
if (!folder)
{
system("mkdir out");
cout << "Created ./out folder"
<< "\n";
}
else
{
system("rm out/*");
cout << "Cleaned ./out folder"
<< "\n";
}
// copy the input file to the output folder
system("cp halo.txt out/halo_0.txt");
int total_particles;
file >> total_particles;
cout << "Number of particles: " << total_particles << "\n";
cout << "Number of threads: " << omp_get_max_threads() << "\n";
// Initialize particles
Particle *particles = new Particle[total_particles];
for (int i = 0; i < total_particles; i++)
{
file >> particles[i].pos[0] >> particles[i].pos[1] >> particles[i].pos[2] >> particles[i].vel[0] >> particles[i].vel[1] >> particles[i].vel[2];
particles[i].mass = Mass;
}
// Timer
double tstart, tend, tmid;
tstart = omp_get_wtime();
tmid = omp_get_wtime();
// Start the iteration
cout << "Set up a tree with size " << Size
<< " and origin at " << Origin[0] << ", " << Origin[1] << ", " << Origin[2] << "\n";
cout << "Start the iteration...(N = " << N << ")"
<< "\n";
for (int i = 0; i < N; i++)
{
// Create the tree
Node root = Node(Origin, Size);
for (int j = 0; j < total_particles; j++)
{
root.add(particles[j]);
}
#pragma omp parallel for
// Update the particles (velocity)
for (int j = 0; j < total_particles; j++)
{
update(dt, root, particles[j]);
}
#pragma omp parallel for
// Update the particles (position)
for (int j = 0; j < total_particles; j++)
{
for (int k = 0; k < 3; k++)
{
particles[j].pos[k] = particles[j].new_pos[k];
}
}
// Write the particles positions to file
if (i % save_iter == 0 && i != 0)
{
cout << "Iteration " << i << " done!"
<< "\n";
cout << "Time elapsed: " << omp_get_wtime() - tmid << " seconds"
<< "\n";
tmid = omp_get_wtime();
ofstream file("out/halo_" + to_string(i) + ".txt");
// file << total_particles << "\n";
for (int j = 0; j < total_particles; j++)
{
file << particles[j].pos[0] << " " << particles[j].pos[1] << " " << particles[j].pos[2] << " "
<< particles[j].vel[0] << " " << particles[j].vel[1] << " " << particles[j].vel[2] << "\n";
}
file.close();
}
}
// Calculate time elapsed
tend = omp_get_wtime();
cout << "Time elapsed: " << tend - tstart << " seconds"
<< "\n";
cout << "Time per iteration: " << (tend - tstart) / N << " seconds"
<< "\n";
cout << "Done!"
<< "\n";
return 0;
}