From f1a4cf4a24a768a30721465c8c6f939b3b68b4f9 Mon Sep 17 00:00:00 2001 From: Cory Date: Tue, 19 Dec 2023 16:08:58 +0100 Subject: [PATCH] Make some changes --- include/WaveSimulation.hpp | 20 +++++---- lib/WaveSimulation.cpp | 83 +++++++++++++++++++++++++++----------- src/main.cpp | 69 ++++++++++++++++--------------- 3 files changed, 107 insertions(+), 65 deletions(-) diff --git a/include/WaveSimulation.hpp b/include/WaveSimulation.hpp index 0ef9acf..ae9ba7e 100644 --- a/include/WaveSimulation.hpp +++ b/include/WaveSimulation.hpp @@ -26,20 +26,26 @@ protected: double dt; double T; + void build_A(); + void build_B(); + void initialize_U(double x_c, double y_c, double sigma_x, double sigma_y, + double p_x, double p_y); + void build_V(double thickness, double pos_x, double aperture_separation, + double aperture, uint32_t slits); public: int32_t N; arma::cx_mat V; arma::cx_mat U; arma::sp_cx_mat A; - WaveSimulation(double h, double dt, double T); - virtual void solve(std::ofstream& ofile); - void build_A(); - void build_B(); - void initialize_U(double x_c, double y_c, double sigma_x, double sigma_y, - double p_x, double p_y); + WaveSimulation(double h, double dt, double T, double x_c, double y_c, + double sigma_x, double sigma_y, double p_x, double p_y, + double thickness, double pos_x, double ap_sep, double ap, + uint32_t slits); + WaveSimulation(double h, double dt, double T, double x_c, double y_c, + double sigma_x, double sigma_y, double p_x, double p_y); + virtual void solve(std::ofstream &ofile); void write_U(std::ofstream &ofile); void step(); - void build_V(double thickness, double pos_x, double aperture_sparation, double aperture, uint32_t slits); }; #endif diff --git a/lib/WaveSimulation.cpp b/lib/WaveSimulation.cpp index 056cb05..7be2fd5 100644 --- a/lib/WaveSimulation.cpp +++ b/lib/WaveSimulation.cpp @@ -18,17 +18,36 @@ #include #include -WaveSimulation::WaveSimulation(double h, double dt, double T) +WaveSimulation::WaveSimulation(double h, double dt, double T, double x_c, + double y_c, double sigma_x, double sigma_y, + double p_x, double p_y, double thickness, + double pos_x, double ap_sep, double ap, + uint32_t slits) { this->dt = dt; this->h = h; this->T = T; this->M = 1. / h; this->N = M - 2; - this->V.set_size(this->N, this->N); - this->V.fill(0.); - this->U.set_size(this->N, this->N); - this->U.fill(0.); + this->build_V(thickness, pos_x, ap_sep, ap, slits); + this->initialize_U(x_c, y_c, sigma_x, sigma_y, p_x, p_y); + this->build_A(); + this->build_B(); +} + +WaveSimulation::WaveSimulation(double h, double dt, double T, double x_c, double y_c, + double sigma_x, double sigma_y, double p_x, double p_y) +{ + this->dt = dt; + this->h = h; + this->T = T; + this->M = 1. / h; + this->N = M - 2; + this->build_V(0.,0.,0.,0., 0); + this->initialize_U(x_c, y_c, sigma_x, sigma_y, p_x, p_y); + this->build_A(); + this->build_B(); + } void WaveSimulation::solve(std::ofstream &ofile) @@ -43,7 +62,6 @@ void WaveSimulation::solve(std::ofstream &ofile) void WaveSimulation::step() { - DEBUG("Inside step"); arma::cx_vec tmp = this->B * this->U.as_col(); arma::spsolve(this->U, this->A, tmp); } @@ -111,6 +129,7 @@ void WaveSimulation::build_B() void WaveSimulation::initialize_U(double x_c, double y_c, double sigma_x, double sigma_y, double p_x, double p_y) { + this->U.set_size(this->N, this->N); double x, y, diff_x, diff_y; std::complex sum = 0.; for (size_t j = 0; j < this->U.n_cols; j++) { @@ -145,25 +164,42 @@ void WaveSimulation::build_V(double thickness, double pos_x, double aperture_separation, double aperture, uint32_t slits) { - uint32_t mid_y = this->N / 2 - (this->N % 2 == 0); + this->V.set_size(this->N, this->N); + this->V.fill(0.); + if (slits == 0) { + return; + } arma::cx_vec res; if (slits % 2 == 0) { - res = arma::cx_vec(aperture_separation/this->h,arma::fill::value(1e10)); - for (size_t i=0; i < slits; i+=2) { - res = arma::join_cols(res, arma::cx_vec(aperture/this->h,arma::fill::zeros)); - res = arma::join_cols(arma::cx_vec(aperture/this->h,arma::fill::zeros), res); - res = arma::join_cols(res, arma::cx_vec(aperture_separation/this->h,arma::fill::value(1e10))); - res = arma::join_cols(arma::cx_vec(aperture_separation/this->h,arma::fill::value(1e10)), res); + res = arma::cx_vec(aperture_separation / this->h, + arma::fill::value(1e10)); + for (size_t i = 0; i < slits; i += 2) { + res = arma::join_cols( + res, arma::cx_vec(aperture / this->h, arma::fill::zeros)); + res = arma::join_cols( + arma::cx_vec(aperture / this->h, arma::fill::zeros), res); + res = + arma::join_cols(res, arma::cx_vec(aperture_separation / this->h, + arma::fill::value(1e10))); + res = arma::join_cols(arma::cx_vec(aperture_separation / this->h, + arma::fill::value(1e10)), + res); } } else { - res = arma::cx_vec(aperture/this->h,arma::fill::value(0)); - for (size_t i=0; i < slits-1; i+=2) { - res = arma::join_cols(res, arma::cx_vec(aperture_separation/this->h,arma::fill::value(1e10))); - res = arma::join_cols(arma::cx_vec(aperture_separation/this->h,arma::fill::value(1e10)), res); - res = arma::join_cols(res, arma::cx_vec(aperture/this->h,arma::fill::zeros)); - res = arma::join_cols(arma::cx_vec(aperture/this->h,arma::fill::zeros), res); + res = arma::cx_vec(aperture / this->h, arma::fill::value(0)); + for (size_t i = 0; i < slits - 1; i += 2) { + res = + arma::join_cols(res, arma::cx_vec(aperture_separation / this->h, + arma::fill::value(1e10))); + res = arma::join_cols(arma::cx_vec(aperture_separation / this->h, + arma::fill::value(1e10)), + res); + res = arma::join_cols( + res, arma::cx_vec(aperture / this->h, arma::fill::zeros)); + res = arma::join_cols( + arma::cx_vec(aperture / this->h, arma::fill::zeros), res); } } if (res.size() > this->N) { @@ -171,10 +207,11 @@ void WaveSimulation::build_V(double thickness, double pos_x, } uint32_t fill = (this->N - res.size()) / 2; res = arma::join_cols(arma::cx_vec(fill, arma::fill::value(1e10)), res); - res = arma::join_cols(res, arma::cx_vec(fill + ((this->N - res.size()) % 2), arma::fill::value(1e10))); + res = arma::join_cols(res, arma::cx_vec(fill + ((this->N - res.size()) % 2), + arma::fill::value(1e10))); - uint32_t start = pos_x/this->h - thickness/this->h/2; - for (size_t i=0; i < thickness/this->h; i++) { - this->V.col(start+i) = res; + uint32_t start = pos_x / this->h - thickness / this->h / 2; + for (size_t i = 0; i < thickness / this->h; i++) { + this->V.col(start + i) = res; } } diff --git a/src/main.cpp b/src/main.cpp index c6891b4..8c2cb6a 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -11,19 +11,14 @@ * */ #include "WaveSimulation.hpp" #include "utils.hpp" +#include void probability_deviation() { - WaveSimulation sim_narrow(.005, 2.5e-5, .008); - sim_narrow.initialize_U(.25, .5, .05, .05, 200., 0.); - sim_narrow.build_A(); - sim_narrow.build_B(); + WaveSimulation sim_narrow(.005, 2.5e-5, .008, .25, .5, .05, .05, 200., 0.); - WaveSimulation sim_wide(.005, 2.5e-5, .008); - sim_wide.build_V(.02, .5, .05, .05, 2); - sim_wide.initialize_U(.25, .5, .05, .10, 200., 0.); - sim_wide.build_A(); - sim_wide.build_B(); + WaveSimulation sim_wide(.005, 2.5e-5, .008, .25, .5, .05, .10, 200., 0., + 0.02, .5, .05, .05, 2); std::ofstream ofile; utils::mkpath("data"); @@ -41,7 +36,8 @@ void probability_deviation() sim_narrow.step(); sim_wide.step(); - ofile << i << ',' << sum_narrow << ',' << sum_wide << '\n'; + ofile << i << ',' << utils::scientific_format(sum_narrow) << ',' + << utils::scientific_format(sum_wide) << '\n'; } ofile.close(); @@ -49,11 +45,8 @@ void probability_deviation() void color_map() { - WaveSimulation sim(.005, 2.5e-5, .008); - sim.build_V(.02, .5, .05, .05, 2); - sim.initialize_U(.25, .5, .05, .10, 200., 0.); - sim.build_A(); - sim.build_B(); + WaveSimulation sim(.005, 2.5e-5, .002, .25, .5, .05, .20, 200., 0., 0.02, + .5, .05, .05, 2); std::ofstream ofile; ofile.open("data/color_map.txt"); @@ -72,50 +65,56 @@ void color_map() void detector_screen() { - WaveSimulation sim(.005, 2.5e-5, .008); - sim.build_V(.02, .5, .05, .05, 1); - sim.initialize_U(.25, .5, .05, .10, 200., 0.); - sim.build_A(); - sim.build_B(); + WaveSimulation *sim = new WaveSimulation( + .005, 2.5e-5, .002, .25, .5, .05, .20, 200., 0., 0.02, .5, .05, .05, 1); std::ofstream ofile; utils::mkpath("data/screen"); ofile.open("data/screen/single_slit.txt"); - ofile << sim.N << '\n'; + ofile << sim->N << '\n'; for (size_t i = 0; i < 80; i++) { - sim.step(); + sim->step(); } - sim.write_U(ofile); + sim->write_U(ofile); ofile.close(); - sim.build_V(.02, .5, .05, .05, 2); - sim.initialize_U(.25, .5, .05, .10, 200., 0.); + delete sim; + sim = new WaveSimulation(.005, 2.5e-5, .002, .25, .5, .05, .20, 200., 0., + 0.02, .5, .05, .05, 2); ofile.open("data/screen/double_slit.txt"); - ofile << sim.N << '\n'; + ofile << sim->N << '\n'; for (size_t i = 0; i < 80; i++) { - sim.step(); + sim->step(); } - sim.write_U(ofile); + sim->write_U(ofile); ofile.close(); - sim.build_V(.02, .5, .05, .05, 3); - sim.initialize_U(.25, .5, .05, .10, 200., 0.); + delete sim; + sim = new WaveSimulation(.005, 2.5e-5, .002, .25, .5, .05, .20, 200., 0., + 0.02, .5, .05, .05, 3); ofile.open("data/screen/triple_slit.txt"); - ofile << sim.N << '\n'; + ofile << sim->N << '\n'; for (size_t i = 0; i < 80; i++) { - sim.step(); + sim->step(); } - sim.write_U(ofile); + sim->write_U(ofile); ofile.close(); } int main() { - probability_deviation(); - color_map(); + //probability_deviation(); + //color_map(); detector_screen(); + + //std::ofstream ofile; + //ofile.open("test.txt"); + //WaveSimulation sim(.005, 2.5e-5, .008, .25, .5, .05, .10, 200., 0., + //0.02, .5, .05, .05, 2); + //sim.solve(ofile); + return 0; }