From 6275765e5e8d1ed522a75fd5ccbffeb30154fb18 Mon Sep 17 00:00:00 2001 From: nilsuatlan Date: Thu, 1 Aug 2024 14:08:37 +0300 Subject: [PATCH] Addition of integration with parallel model solution --- benchmark/benchmark_Cell_ECM.hpp | 2 +- src/main.cpp | 4 +- src/modules/Module_p.cpp | 1 + src/modules/Module_p.hpp | 8 +- tests/integration/parallel_model_sln.cpp | 107 +++++++++++++++++++++++ 5 files changed, 115 insertions(+), 7 deletions(-) create mode 100644 tests/integration/parallel_model_sln.cpp diff --git a/benchmark/benchmark_Cell_ECM.hpp b/benchmark/benchmark_Cell_ECM.hpp index 0869b58..3de7ac0 100644 --- a/benchmark/benchmark_Cell_ECM.hpp +++ b/benchmark/benchmark_Cell_ECM.hpp @@ -123,7 +123,7 @@ inline void run_Cell_ECM_parallel_3_default_pulse() }; auto module = Module_p("Par3", T_ENV, true, false, std::size(cs), 1, 1); - module.setSUs(cs, false); + module.setSUs(cs, false); // set storage unit module.setBlockDegAndTherm(true); diff --git a/src/main.cpp b/src/main.cpp index 23545d9..7ce182e 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -206,7 +206,7 @@ int main() //!< Benchmarks: // slide::benchmarks::run_Cell_Bucket(); - // slide::benchmarks::run_Cell_ECM(); + slide::benchmarks::run_Cell_ECM(); // slide::benchmarks::run_Cell_SPM_1(1); // slide::benchmarks::run_Cell_SPM_2(1); // slide::benchmarks::run_LP_case_SmallPack(); @@ -223,7 +223,7 @@ int main() // slide::benchmarks::run_Cell_ECM_2_RC_single_default_pulse(); // slide::benchmarks::run_Cell_ECM_2_RC_single_default_CCCV(); - // slide::benchmarks::run_Cell_ECM_parallel_3_default_pulse(); + slide::benchmarks::run_Cell_ECM_parallel_3_default_pulse(); // slide::benchmarks::run_Cell_ECM_parallel_3_default_CCCV(); diff --git a/src/modules/Module_p.cpp b/src/modules/Module_p.cpp index 8f7e15f..9fe3239 100644 --- a/src/modules/Module_p.cpp +++ b/src/modules/Module_p.cpp @@ -204,6 +204,7 @@ Status Module_p::setVoltage(double Vnew, bool checkI, bool print) return StatusNow; // #TODO add some voltage/current etc. Also return max iter condition!!!! } +// NİLSU Status Module_p::setCurrent(double Inew, bool checkV, bool print) { /* diff --git a/src/modules/Module_p.hpp b/src/modules/Module_p.hpp index da2d7c7..e041cc3 100644 --- a/src/modules/Module_p.hpp +++ b/src/modules/Module_p.hpp @@ -25,10 +25,10 @@ class Module_p : public Module : Module(ID_, Ti, print, pari, Ncells_, coolControl, cooltype) {} //!< the voltage limits are the most constraining limits of all cells ie the highest Vmin of SUs is the Vmin of the module - double Vmin() const override { return transform_max(SUs, free::get_Vmin); } - double VMIN() const override { return transform_max(SUs, free::get_VMIN); } - double Vmax() const override { return transform_min(SUs, free::get_Vmax); } - double VMAX() const override { return transform_min(SUs, free::get_VMAX); } + double Vmin() const override { return transform_max(SUs, free::get_Vmin); } // soft limit + double VMIN() const override { return transform_max(SUs, free::get_VMIN); } // hard limit + double Vmax() const override { return transform_min(SUs, free::get_Vmax); } // soft limit + double VMAX() const override { return transform_min(SUs, free::get_VMAX); } // hard limit: to stop the simulation completely double I() const override { return transform_sum(SUs, free::get_I); } //!< the current is the sum of the current of each cell. Returns 0 if empty. double Cap() const override { return transform_sum(SUs, free::get_Cap); } //!< module capacity (sum of cells) diff --git a/tests/integration/parallel_model_sln.cpp b/tests/integration/parallel_model_sln.cpp new file mode 100644 index 0000000..73effb3 --- /dev/null +++ b/tests/integration/parallel_model_sln.cpp @@ -0,0 +1,107 @@ +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +using namespace std; // bu yerine fonksiyonun içine ekle +using namespace Eigen; + +double ocv_eval(const vector& ocv_coefs, double z) { + double ocv = 0.0; + for (size_t i = 0; i < ocv_coefs.size(); ++i) { + ocv += ocv_coefs[i] * pow(z, ocv_coefs.size() - 1 - i); + } + return ocv; +} // already existing function --> otomatik olarak değer fonksiyona gelecek + +void new_compute_A11_A12_A21_A22(int n_par, VectorXd R, VectorXd C, VectorXd Q, VectorXd tau, VectorXd r, MatrixXd &A11, MatrixXd &A12, MatrixXd &A21, MatrixXd &A22, MatrixXd &m) { + // Initialize A11 and A12 + A11 = MatrixXd::Zero(2 * n_par, 2 * n_par); + A12 = MatrixXd::Zero(2 * n_par, 2 * n_par); + + A11.block(0, 0, 2, 2) << 0, 0, 0, -tau(0); + A12.block(0, 0, 2, 2) << 1/Q(0), 0, 0, 1/C(0); + + for (int i = 1; i < n_par; ++i) { + A11.block(2 * i, 2 * i, 2, 2) << 0, 0, 0, -tau(i); + A12.block(2 * i, 2 * i, 2, 2) << 1/Q(i), 0, 0, 1/C(i); + } + + // Initialize A22 + A22 = MatrixXd::Zero(n_par, n_par); + VectorXd main_diag = -(r + R); + main_diag(0) = r(0); + + VectorXd lower_diag = r.head(n_par - 1); + lower_diag(0) = r(0); + + A22.diagonal() = main_diag; + A22.diagonal(-1) = lower_diag; + + for (int i = 1; i < n_par; ++i) { + for (int j = i + 1; j < n_par; ++j) { + A22(i, j) = -R(i); + } + } + + A22.row(0).setOnes(); + cout << "A22:\n" << A22 << endl; + + // Initialize A21 + A21 = MatrixXd::Zero((n_par - 1) * 2, n_par); + for (int i = 0; i < n_par - 1; ++i) { + A21.block(2 * i, 0, 2, n_par) = MatrixXd::Identity(2, n_par); + } + + // Calculate m (inverse of A22) + m = MatrixXd::Zero(n_par, n_par); + double Rsum_inv = 1 / R.sum(); + + for (int ktil = 0; ktil < n_par; ++ktil) { + for (int i = 0; i < n_par - 1; ++i) { + if (ktil - i - 1 == 0) { + m(i + 1, i) = pow(1 / R(i + 1), 2) * Rsum_inv - 1 / R(i + 1); + } else { + m(ktil, i) = 1 / (R(ktil) * R(i + 1)) * Rsum_inv; + } + } + } + + m(n_par - 1, n_par - 1) = 1 / (R(n_par - 1) * Rsum_inv); + for (int i = 0; i < n_par - 1; ++i) { + m(i, n_par - 1) = 1 / (R(i) * Rsum_inv); + } + + MatrixXd inv_A22 = A22.inverse(); + double error = (inv_A22 - m).norm(); + + cout << "Error: " << error << endl; +} + +int main() { + int n_par = 3; + VectorXd R(n_par), C(n_par), Q(n_par), tau(n_par), r(n_par); + + R << 0, 0, 0; + C << 634.0, 634.0, 634.0; + Q << 9000.0, 9000.0, 9000.0; + tau << 0.04, 0.04, 0.04; + r << 0.029, 0.029, 0.029; + + MatrixXd A11, A12, A21, A22, m; + new_compute_A11_A12_A21_A22(n_par, R, C, Q, tau, r, A11, A12, A21, A22, m); + + cout << "A11:\n" << A11 << endl; + cout << "A12:\n" << A12 << endl; + cout << "A21:\n" << A21 << endl; + cout << "A22:\n" << A22 << endl; + cout << "m:\n" << m << endl; + + return 0; +}