GSI_Toolbox 1.0.0
A toolbox for Gas-Surface Interaction simulations
Loading...
Searching...
No Matches
Raytracer.h
1
10template <typename T>
11class Raytracer {
12
13 private:
14
18 Surface<T> surface;
19
23 Gas<T> gas;
24
28 void (*local_kernel)(T*, T*, T*, T*, T*);
29
33 Matrix<T> incident_velocity;
34
38 std::string sim_name;
39
43 unsigned long num_particles;
44
45 public:
46
53
67 Raytracer(Surface<T> surface, Gas<T> gas, void (*local_kernel)(T*, T*, T*, T*, T*), Matrix<T> incident_velocity, unsigned long num_particles, std::string sim_name);
68
73
82 std::vector<trajectory<T>> simulate();
83
92 std::vector<trajectory<T>> import_data(std::string filename);
93
102 void save(std::vector<trajectory<T>> trajectories, std::string filename);
103
113 void combine_files(std::string filename, int num_procs);
114
120 void set_num_particles(unsigned long num_particles) { this->num_particles = num_particles; }
121
127 std::string& get_sim_name() { return sim_name; }
128
134 unsigned long get_num_particles() { return num_particles; }
135
136};
137
138template <typename T>
139Raytracer<T>::Raytracer(Surface<T> surface, Gas<T> gas, void (*local_kernel)(T*, T*, T*, T*, T*), Matrix<T> incident_velocity, unsigned long num_particles, std::string sim_name) {
140
141 this->surface = surface;
142 this->gas = gas;
143 this->local_kernel = local_kernel;
144 this->incident_velocity = incident_velocity;
145 this->num_particles = num_particles;
146 this->sim_name = sim_name;
147
148}
149
150template <typename T>
151std::vector<trajectory<T>> Raytracer<T>::simulate() {
152
153 Matrix<T> unit_incident = this->incident_velocity / this->incident_velocity.norm() * (-1.0);
154
155 T theta_i = std::acos(- this->incident_velocity(2, 0) / this->incident_velocity.norm());
156
157 long num_col = 0;
158 bool collided = true;
159 unsigned int collision_index;
160
161 std::vector<trajectory<T>> trajectories;
162
163 std::random_device rd;
164 std::seed_seq fullSeed{rd(), rd(), rd(), rd(), rd(), rd(), rd(), rd(), rd(), rd()};
165 std::mt19937 rng(fullSeed);
166 std::uniform_real_distribution<T> uniformDist(0.0f, 1.0f);
167 std::normal_distribution<T> normDist(0.0f, 1.0f);
168
169 T offset_size = 15.0;
170
171 for(long i = 0; i < this->num_particles; i ++) {
172
173 std::vector<Matrix<T>> p_positions;
174 std::vector<Matrix<T>> p_velocities;
175
176 T m_dot = 1.0;
177 num_col = 0;
178
179 collided = true;
180
181 T rotation_angle = uniformDist(rng) * 2.0 * M_PI;
182
183 surface.get_geometry().rotate(2, rotation_angle);
184
185 Matrix<T> origin = unit_incident / std::cos(theta_i) * 100.0;
186 Matrix<T> offset(3, 1, {2.0 * (uniformDist(rng) - 0.5) * offset_size, 2.0 * (uniformDist(rng) - 0.5) * offset_size, 0.0});
187 origin += offset;
188 Ray<T> ray(origin, incident_velocity, m_dot, 1, -1);
189 while (collided) {
190
191 std::vector<Ray<T>> new_rays = ray.propagate(this->surface.get_geometry(), this->gas, this->local_kernel, &collision_index);
192 if(new_rays.size() == 0) {
193 collided = false;
194 }
195 else {
196 ray = new_rays.back();
197 p_velocities.push_back(ray.get_velocity());
198 p_positions.push_back(ray.get_origin());
199 }
200 num_col ++;
201 }
202
203 trajectory<T> trajectory = {p_velocities, p_positions};
204 trajectories.push_back(trajectory);
205 }
206 return trajectories;
207}
208
209template<typename T>
210void Raytracer<T>::save(std::vector<trajectory<T>> trajectories, std::string filename) {
211
212 std::ofstream myfile;
213 myfile.open(filename);
214
215 unsigned long num_particles = trajectories.size();
216
217 for(auto i = 0; i < num_particles; i ++) {
218
219 trajectory<T> p_trajectory = trajectories[i];
220 unsigned long num_bounces = p_trajectory.positions.size();
221
222 myfile << num_bounces << "\n";
223
224 for(auto j = 0; j < 3; j ++) {
225 for(auto k = 0; k < num_bounces; k ++) {
226 myfile << p_trajectory.positions[k](j, 0) << " ";
227 }
228 myfile << "\n";
229 }
230
231 for(auto j = 0; j < 3; j ++) {
232 for(auto k = 0; k < num_bounces; k ++) {
233 myfile << p_trajectory.velocities[k](j, 0) << " ";
234 }
235 myfile << "\n";
236 }
237
238 }
239 myfile.close();
240
241}
242
243template <typename T>
244std::vector<trajectory<T>> Raytracer<T>::import_data(std::string filename) {
245
246 std::ifstream linecounter(filename);
247 std::string line;
248 std::vector<trajectory<T>> trajectories;
249
250 long num_lines = 0;
251
252 while( std::getline(linecounter, line)) {
253 ++ num_lines;
254 }
255
256 long num_particles = num_lines / 7;
257 long num_bounces;
258
259 std::ifstream myfile(filename);
260
261 for(long i = 0; i < num_particles; i ++) {
262
263 myfile >> num_bounces;
264 Matrix<T> position(3, 1);
265 Matrix<T> velocity(3, 1);
266 std::vector<Matrix<T>> p_positions(num_bounces, Matrix<T>(3, 1));
267 std::vector<Matrix<T>> p_velocities(num_bounces, Matrix<T>(3, 1));
268
269 for(long k = 0; k < 3; k ++) {
270 for(long j = 0; j < num_bounces; j ++) {
271 myfile >> p_positions[j](k, 0);
272 }
273 }
274 for(long k = 0; k < 3; k ++) {
275 for(long j = 0; j < num_bounces; j ++) {
276 myfile >> p_velocities[j](k, 0);
277 }
278 }
279 trajectory<T> trajectory = {p_velocities, p_positions};
280 trajectories.push_back(trajectory);
281 }
282
283 return trajectories;
284}
285
286template <typename T>
287void Raytracer<T>::combine_files(std::string filename, int num_procs) {
288
289 std::vector<trajectory<T>> trajectories_total;
290
291 for( long i = 0; i < num_procs; i ++) {
292 std::vector<trajectory<T>> trajectories = import_data(filename + "_" + std::to_string(i) + ".dat");
293 trajectories_total.insert(trajectories_total.end(), trajectories.begin(), trajectories.end());
294 std::remove((filename + "_" + std::to_string(i) + ".dat").c_str());
295 }
296 save(trajectories_total, filename + ".dat");
297}
Class representing a Gas with properties such as density, temperature, molar mass,...
Definition Gas.h:10
This class implements a matrix object used for linear algebra and vectorized operations.
Definition Matrix.h:25
T norm()
Calculates the Frobenius norm of the matrix.
Definition Matrix.h:674
Class representing a Ray for simulating particle movement in a geometric domain.
Definition Ray.h:10
Matrix< T > & get_velocity()
Gets the velocity of the ray.
Definition Ray.h:91
std::vector< Ray< T > > propagate(Geometry< T > &geometry, const Gas< T > &gas, void(*kernel)(T *, T *, T *, T *, T *), unsigned int *collision_index)
Propagates the ray through the geometry and handles collision and reflection.
Definition Ray.h:356
Matrix< T > & get_origin()
Gets the origin of the ray.
Definition Ray.h:105
A class that simulates ray tracing through a surface and gas environment.
Definition Raytracer.h:11
void combine_files(std::string filename, int num_procs)
Combines output files from multiple processes into a single file.
Definition Raytracer.h:287
void save(std::vector< trajectory< T > > trajectories, std::string filename)
Saves the simulated trajectories to a file.
Definition Raytracer.h:210
Raytracer()
Default constructor.
Definition Raytracer.h:52
unsigned long get_num_particles()
Retrieves the number of particles in the simulation.
Definition Raytracer.h:134
std::vector< trajectory< T > > simulate()
Simulates ray tracing and generates ray trajectories.
Definition Raytracer.h:151
void set_num_particles(unsigned long num_particles)
Sets the number of particles to simulate.
Definition Raytracer.h:120
std::string & get_sim_name()
Retrieves the name of the simulation.
Definition Raytracer.h:127
std::vector< trajectory< T > > import_data(std::string filename)
Imports trajectory data from a file.
Definition Raytracer.h:244
~Raytracer()
Destructor for the Raytracer class.
Definition Raytracer.h:72
A class representing a surface with geometry and statistical properties.
Definition Surface.h:11
Geometry< T > & get_geometry()
Retrieves the geometry of the surface.
Definition Surface.h:93
Represents a trajectory of a certain type.
Definition PG_Kernel.h:10