diff --git a/libnufi_lib.a b/libnufi_lib.a index d996a15..0b5fb6c 100644 Binary files a/libnufi_lib.a and b/libnufi_lib.a differ diff --git a/nufi/nufi_solver.h b/nufi/nufi_solver.h index 756cf2b..3d8669a 100644 --- a/nufi/nufi_solver.h +++ b/nufi/nufi_solver.h @@ -1,6 +1,7 @@ #ifndef NUFI_SOLVER_H #define NUFI_SOLVER_H +#include #include #include #include @@ -30,6 +31,9 @@ private: double Lx = Parameters::LX; + double x_min = Parameters::X_DOMAIN_LEFT; + double x_max = Parameters::X_DOMAIN_RIGHT; + std::vector rho; unsigned int order; diff --git a/nufi/save_results.h b/nufi/save_results.h index 8e6258b..e1a3247 100644 --- a/nufi/save_results.h +++ b/nufi/save_results.h @@ -23,4 +23,7 @@ void save_Efield(unsigned int n, unsigned int Nx_out, const std::string &filename); + +void save_space_vector(const std::vector& vals, const std::string& filename, size_t it); + #endif diff --git a/rho_E.mp4 b/rho_E.mp4 index dd19f84..0f4af77 100644 Binary files a/rho_E.mp4 and b/rho_E.mp4 differ diff --git a/src/nufi_solver.cc b/src/nufi_solver.cc index a060bda..de4b9a0 100644 --- a/src/nufi_solver.cc +++ b/src/nufi_solver.cc @@ -28,7 +28,7 @@ double NuFISolver::eval_ftilda(unsigned int n, if ( n == 0 ) return f0(x,u); const size_t stride_x = 1; - const size_t stride_t = stride_x*(Parameters::SPLINE_NX + Parameters::SPLINE_ORDER - 1); + const size_t stride_t = stride_x*(Nx + Parameters::SPLINE_ORDER - 1); double Ex; const double *c; @@ -81,7 +81,7 @@ void NuFISolver::run() if ( rho == nullptr ) throw std::bad_alloc {}; - Gradient grad(Parameters::X_DOMAIN_LEFT, Parameters::X_DOMAIN_RIGHT, Parameters::SPLINE_NX); + Gradient grad(x_min, x_max, Nx); for (unsigned int it = 0; it < Nt; ++it) { @@ -99,28 +99,17 @@ void NuFISolver::run() rho.get()[i] = ith_rho; } - poisson.set_rhs_function(std::make_unique>(rho.get(), Parameters::SPLINE_NX)); + poisson.set_rhs_function(std::make_unique>(rho.get(), Nx)); poisson.solve_step(); - // std::vector sampled_potential = poisson.sample_electric_potential(Parameters::X_DOMAIN_LEFT, Parameters::X_DOMAIN_RIGHT, Parameters::SPLINE_NX); + std::vector sampled_potential = poisson.sample_electric_potential(x_min, x_max, Nx); - // sampled_potential.erase(sampled_potential.begin()); - // sampled_potential.erase(sampled_potential.end()-1); - // - // std::ofstream file("results/potential_" + std::to_string(it) + ".dat"); - // file << sampled_potential.size() << "\n"; - // file << Parameters::X_DOMAIN_LEFT << " " << Parameters::X_DOMAIN_RIGHT << "\n"; - // file << std::fixed << std::setprecision(8); - // for (double val : sampled_potential) file << val << "\n"; - // - // - // std::vector E_vals = grad.compute(sampled_potential); - std::vector E_vals = poisson.sample_electric_field(Parameters::X_DOMAIN_LEFT, Parameters::X_DOMAIN_RIGHT,Parameters::SPLINE_NX); - std::ofstream file("results/potential_" + std::to_string(it) + ".dat"); - file << E_vals.size() << "\n"; - file << Parameters::X_DOMAIN_LEFT << " " << Parameters::X_DOMAIN_RIGHT << "\n"; - file << std::fixed << std::setprecision(8); - for (double val : E_vals) file << val << "\n"; + save_space_vector(sampled_potential, "potential", it); + + std::vector E_vals = grad.compute(sampled_potential); + // std::vector E_vals = poisson.sample_electric_field(x_min, x_max, Nx); + + save_space_vector(E_vals, "electric", it); double* current_coeffs = coeffs.get() + it*stride_t; interpolate(current_coeffs, E_vals.data()); diff --git a/src/save_results.cc b/src/save_results.cc index 8654d66..dfb8f7d 100644 --- a/src/save_results.cc +++ b/src/save_results.cc @@ -2,7 +2,9 @@ #include "nufi/parameters.h" #include +#include #include +#include #include "nufi/nufi_solver.h" @@ -102,3 +104,15 @@ void save_Efield(unsigned int n, } file.close(); } + +void save_space_vector(const std::vector& vals, const std::string& filename, size_t it) +{ + std::ofstream file("results/" + filename + "_" + std::to_string(it) + ".dat"); + + if (!file) throw std::runtime_error("failed to start file in results/"); + + file << vals.size() << "\n"; + file << Parameters::X_DOMAIN_LEFT << " " << Parameters::X_DOMAIN_RIGHT << "\n"; + file << std::fixed << std::setprecision(8); + for (double val : vals) file << val << "\n"; +}