From 0a79604e5513afbe6d6ff8f44bad99fa6025a80d Mon Sep 17 00:00:00 2001 From: superwhiskers Date: Mon, 2 Jun 2025 10:39:32 -0400 Subject: [PATCH 1/3] address the remaining warnings resulting from the new linting flags --- src/Solver/Optimization/DynamicConstraint.cpp | 88 ++++++++++++++----- src/Solver/Optimization/DynamicObjective.cpp | 82 ++++++++++++----- 2 files changed, 125 insertions(+), 45 deletions(-) diff --git a/src/Solver/Optimization/DynamicConstraint.cpp b/src/Solver/Optimization/DynamicConstraint.cpp index 784d4e3f..363d42fa 100644 --- a/src/Solver/Optimization/DynamicConstraint.cpp +++ b/src/Solver/Optimization/DynamicConstraint.cpp @@ -1,9 +1,9 @@ -#include "DynamicConstraint.hpp" - #include #include +#include "DynamicConstraint.hpp" + namespace AnalysisManager { namespace IpoptInterface @@ -32,13 +32,13 @@ namespace AnalysisManager // Number of parameters is size of the system plus 1 fictitious parameter // to store the objective value. - n = model_->sizeParams() + 1; + n = static_cast(model_->sizeParams() + 1); // There is one constraint m = 1; // Jacobian is a dense row matrix of length n+1. - nnz_jac_g = model_->sizeParams() + 1; + nnz_jac_g = static_cast(model_->sizeParams() + 1); // Using numerical Hessian. nnz_h_lag = 0; @@ -50,7 +50,12 @@ namespace AnalysisManager } template - bool DynamicConstraint::get_bounds_info(Index n, Number* x_l, Number* x_u, Index m, Number* g_l, Number* g_u) + bool DynamicConstraint::get_bounds_info([[maybe_unused]] Index n, + Number* x_l, + Number* x_u, + [[maybe_unused]] Index m, + Number* g_l, + Number* g_u) { // Check if sizes are set correctly assert(n == (Index) (model_->sizeParams() + 1)); @@ -75,7 +80,15 @@ namespace AnalysisManager } template - bool DynamicConstraint::get_starting_point(Index n, bool init_x, Number* x, bool init_z, Number* z_L, Number* z_U, Index m, bool init_lambda, Number* lambda) + bool DynamicConstraint::get_starting_point([[maybe_unused]] Index n, + [[maybe_unused]] bool init_x, + Number* x, + [[maybe_unused]] bool init_z, + [[maybe_unused]] Number* z_L, + [[maybe_unused]] Number* z_U, + [[maybe_unused]] Index m, + [[maybe_unused]] bool init_lambda, + [[maybe_unused]] Number* lambda) { // Only initial values for x provided. assert(init_x == true); @@ -93,7 +106,10 @@ namespace AnalysisManager } template - bool DynamicConstraint::eval_f(Index n, const Number* x, bool new_x, Number& obj_value) + bool DynamicConstraint::eval_f([[maybe_unused]] Index n, + const Number* x, + [[maybe_unused]] bool new_x, + Number& obj_value) { // Set objective to fictitious optimization parameter x[n-1] obj_value = x[model_->sizeParams()]; @@ -102,7 +118,10 @@ namespace AnalysisManager } template - bool DynamicConstraint::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f) + bool DynamicConstraint::eval_grad_f([[maybe_unused]] Index n, + [[maybe_unused]] const Number* x, + [[maybe_unused]] bool new_x, + Number* grad_f) { // Objective function equals to the fictitious parameter x[n-1]. // Gradient, then assumes the simple form: @@ -114,7 +133,11 @@ namespace AnalysisManager } template - bool DynamicConstraint::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g) + bool DynamicConstraint::eval_g([[maybe_unused]] Index n, + const Number* x, + [[maybe_unused]] bool new_x, + [[maybe_unused]] Index m, + Number* g) { // Update optimization parameters for (IdxT i = 0; i < model_->sizeParams(); ++i) @@ -143,7 +166,14 @@ namespace AnalysisManager } template - bool DynamicConstraint::eval_jac_g(Index n, const Number* x, bool new_x, Index m, Index nele_jac, Index* iRow, Index* jCol, Number* values) + bool DynamicConstraint::eval_jac_g([[maybe_unused]] Index n, + const Number* x, + [[maybe_unused]] bool new_x, + [[maybe_unused]] Index m, + [[maybe_unused]] Index nele_jac, + Index* iRow, + Index* jCol, + Number* values) { // Set Jacobian sparsity pattern ... if (!values) @@ -151,10 +181,10 @@ namespace AnalysisManager for (IdxT i = 0; i < model_->sizeParams(); ++i) { iRow[i] = 0; - jCol[i] = i; + jCol[i] = static_cast(i); } iRow[model_->sizeParams()] = 0; - jCol[model_->sizeParams()] = model_->sizeParams(); + jCol[model_->sizeParams()] = static_cast(model_->sizeParams()); } // ... or compute Jacobian derivatives else @@ -202,23 +232,33 @@ namespace AnalysisManager } template - bool DynamicConstraint::eval_h(Index n, const Number* x, bool new_x, Number obj_factor, Index m, const Number* lambda, bool new_lambda, Index nele_hess, Index* iRow, Index* jCol, Number* values) + bool DynamicConstraint::eval_h([[maybe_unused]] Index n, + [[maybe_unused]] const Number* x, + [[maybe_unused]] bool new_x, + [[maybe_unused]] Number obj_factor, + [[maybe_unused]] Index m, + [[maybe_unused]] const Number* lambda, + [[maybe_unused]] bool new_lambda, + [[maybe_unused]] Index nele_hess, + [[maybe_unused]] Index* iRow, + [[maybe_unused]] Index* jCol, + [[maybe_unused]] Number* values) { return true; } template - void DynamicConstraint::finalize_solution(SolverReturn status, - Index n, - const Number* x, - const Number* z_L, - const Number* z_U, - Index m, - const Number* g, - const Number* lambda, - Number obj_value, - const IpoptData* ip_data, - IpoptCalculatedQuantities* ip_cq) + void DynamicConstraint::finalize_solution([[maybe_unused]] SolverReturn status, + [[maybe_unused]] Index n, + [[maybe_unused]] const Number* x, + [[maybe_unused]] const Number* z_L, + [[maybe_unused]] const Number* z_U, + [[maybe_unused]] Index m, + [[maybe_unused]] const Number* g, + [[maybe_unused]] const Number* lambda, + [[maybe_unused]] Number obj_value, + [[maybe_unused]] const IpoptData* ip_data, + [[maybe_unused]] IpoptCalculatedQuantities* ip_cq) { } diff --git a/src/Solver/Optimization/DynamicObjective.cpp b/src/Solver/Optimization/DynamicObjective.cpp index 81ad25cd..f667a7e8 100644 --- a/src/Solver/Optimization/DynamicObjective.cpp +++ b/src/Solver/Optimization/DynamicObjective.cpp @@ -1,9 +1,9 @@ -#include "DynamicObjective.hpp" - #include #include +#include "DynamicObjective.hpp" + namespace AnalysisManager { namespace IpoptInterface @@ -31,7 +31,7 @@ namespace AnalysisManager assert(model_->sizeQuadrature() == 1); // Number of optimization variables. - n = model_->sizeParams(); + n = static_cast(model_->sizeParams()); // There are no constraints m = 0; @@ -49,7 +49,12 @@ namespace AnalysisManager } template - bool DynamicObjective::get_bounds_info(Index n, Number* x_l, Number* x_u, Index m, Number* g_l, Number* g_u) + bool DynamicObjective::get_bounds_info([[maybe_unused]] Index n, + Number* x_l, + Number* x_u, + [[maybe_unused]] Index m, + [[maybe_unused]] Number* g_l, + [[maybe_unused]] Number* g_u) { // Check if sizes are set correctly assert(n == (Index) model_->sizeParams()); @@ -66,7 +71,15 @@ namespace AnalysisManager } template - bool DynamicObjective::get_starting_point(Index n, bool init_x, Number* x, bool init_z, Number* z_L, Number* z_U, Index m, bool init_lambda, Number* lambda) + bool DynamicObjective::get_starting_point([[maybe_unused]] Index n, + [[maybe_unused]] bool init_x, + Number* x, + [[maybe_unused]] bool init_z, + [[maybe_unused]] Number* z_L, + [[maybe_unused]] Number* z_U, + [[maybe_unused]] Index m, + [[maybe_unused]] bool init_lambda, + [[maybe_unused]] Number* lambda) { // Only initial values for x provided. assert(init_x == true); @@ -81,7 +94,10 @@ namespace AnalysisManager } template - bool DynamicObjective::eval_f(Index n, const Number* x, bool new_x, Number& obj_value) + bool DynamicObjective::eval_f([[maybe_unused]] Index n, + const Number* x, + [[maybe_unused]] bool new_x, + Number& obj_value) { // Update optimization parameters for (IdxT i = 0; i < model_->sizeParams(); ++i) @@ -100,7 +116,10 @@ namespace AnalysisManager } template - bool DynamicObjective::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f) + bool DynamicObjective::eval_grad_f([[maybe_unused]] Index n, + const Number* x, + [[maybe_unused]] bool new_x, + Number* grad_f) { assert(model_->sizeParams() == static_cast(n)); // Update optimization parameters @@ -130,35 +149,56 @@ namespace AnalysisManager } template - bool DynamicObjective::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g) + bool DynamicObjective::eval_g([[maybe_unused]] Index n, + [[maybe_unused]] const Number* x, + [[maybe_unused]] bool new_x, + [[maybe_unused]] Index m, + [[maybe_unused]] Number* g) { return true; } template - bool DynamicObjective::eval_jac_g(Index n, const Number* x, bool new_x, Index m, Index nele_jac, Index* iRow, Index* jCol, Number* values) + bool DynamicObjective::eval_jac_g([[maybe_unused]] Index n, + [[maybe_unused]] const Number* x, + [[maybe_unused]] bool new_x, + [[maybe_unused]] Index m, + [[maybe_unused]] Index nele_jac, + [[maybe_unused]] Index* iRow, + [[maybe_unused]] Index* jCol, + [[maybe_unused]] Number* values) { return true; } template - bool DynamicObjective::eval_h(Index n, const Number* x, bool new_x, Number obj_factor, Index m, const Number* lambda, bool new_lambda, Index nele_hess, Index* iRow, Index* jCol, Number* values) + bool DynamicObjective::eval_h([[maybe_unused]] Index n, + [[maybe_unused]] const Number* x, + [[maybe_unused]] bool new_x, + [[maybe_unused]] Number obj_factor, + [[maybe_unused]] Index m, + [[maybe_unused]] const Number* lambda, + [[maybe_unused]] bool new_lambda, + [[maybe_unused]] Index nele_hess, + [[maybe_unused]] Index* iRow, + [[maybe_unused]] Index* jCol, + [[maybe_unused]] Number* values) { return true; } template - void DynamicObjective::finalize_solution(SolverReturn status, - Index n, - const Number* x, - const Number* z_L, - const Number* z_U, - Index m, - const Number* g, - const Number* lambda, - Number obj_value, - const IpoptData* ip_data, - IpoptCalculatedQuantities* ip_cq) + void DynamicObjective::finalize_solution([[maybe_unused]] SolverReturn status, + [[maybe_unused]] Index n, + [[maybe_unused]] const Number* x, + [[maybe_unused]] const Number* z_L, + [[maybe_unused]] const Number* z_U, + [[maybe_unused]] Index m, + [[maybe_unused]] const Number* g, + [[maybe_unused]] const Number* lambda, + [[maybe_unused]] Number obj_value, + [[maybe_unused]] const IpoptData* ip_data, + [[maybe_unused]] IpoptCalculatedQuantities* ip_cq) { } From 471286a0825f88c63a1402314bb07056f62259ff Mon Sep 17 00:00:00 2001 From: superwhiskers Date: Mon, 2 Jun 2025 10:40:18 -0400 Subject: [PATCH 2/3] remove unnecessary comment in the root `CMakeLists.txt` now that we've added linting flags --- CMakeLists.txt | 3 --- 1 file changed, 3 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index a65f53af..a9982f3e 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -20,9 +20,6 @@ set(PACKAGE_VERSION_PATCH "0") set(PACKAGE_VERSION "${PACKAGE_VERSION_MAJOR}.${PACKAGE_VERSION_MINOR}.${PACKAGE_VERSION_PATCH}") -# TODO: Probably beter to set a debug interface target -# set(CMAKE_CXX_FLAGS_DEBUG "-Wall -O0 -g -DDEBUG") - set(CMAKE_CXX_STANDARD 17) set(CMAKE_CXX_FLAGS "-Wall -Wextra -Wconversion -Wpedantic") From 12fa7be8483f8738af968fb3778ed69d2fa195b3 Mon Sep 17 00:00:00 2001 From: superwhiskers Date: Mon, 2 Jun 2025 16:49:27 -0400 Subject: [PATCH 3/3] finish making all implicit casts explicit --- .../PowerElectronics/Microgrid/Microgrid.cpp | 2 +- .../PowerElectronics/RLCircuit/RLCircuit.cpp | 4 +- .../ScaleMicrogrid/ScaleMicrogrid.cpp | 2 +- src/LinearAlgebra/SparseMatrix/COO_Matrix.hpp | 57 +- .../SynchronousMachine/GENROUwS/Genrou.cpp | 825 +++++++++--------- .../PowerElectronics/Capacitor/Capacitor.cpp | 6 +- src/Model/PowerElectronics/CircuitGraph.hpp | 2 +- .../DistributedGenerator.cpp | 8 +- .../InductionMotor/InductionMotor.cpp | 6 +- .../PowerElectronics/Inductor/Inductor.cpp | 6 +- .../LinearTransformer/LinearTransformer.cpp | 6 +- .../MicrogridBusDQ/MicrogridBusDQ.cpp | 6 +- .../MicrogridLine/MicrogridLine.cpp | 6 +- .../MicrogridLoad/MicrogridLoad.cpp | 6 +- .../PowerElectronics/Resistor/Resistor.cpp | 6 +- .../SynchronousMachine/SynchronousMachine.cpp | 6 +- .../TransmissionLine/TransmissionLine.cpp | 6 +- .../VoltageSource/VoltageSource.cpp | 6 +- src/Model/PowerFlow/Bus/BusPQ.cpp | 16 +- src/Model/PowerFlow/Bus/BusPV.cpp | 16 +- src/Model/PowerFlow/Generator2/Generator2.cpp | 2 +- src/Model/PowerFlow/Generator4/Generator4.cpp | 11 +- .../Generator4Governor/Generator4Governor.cpp | 111 ++- .../Generator4Governor/Generator4Governor.hpp | 20 +- .../Generator4Param/Generator4Param.cpp | 27 +- src/Model/PowerFlow/ModelEvaluatorImpl.hpp | 22 +- src/Solver/Dynamic/Ida.cpp | 33 +- src/Solver/Optimization/DynamicConstraint.cpp | 10 +- src/Solver/Optimization/DynamicObjective.cpp | 10 +- src/Solver/SteadyState/Kinsol.cpp | 9 +- 30 files changed, 631 insertions(+), 622 deletions(-) diff --git a/examples/PowerElectronics/Microgrid/Microgrid.cpp b/examples/PowerElectronics/Microgrid/Microgrid.cpp index 5ca5f35c..3475d789 100644 --- a/examples/PowerElectronics/Microgrid/Microgrid.cpp +++ b/examples/PowerElectronics/Microgrid/Microgrid.cpp @@ -106,7 +106,7 @@ int main(int /* argc */, char const** /* argv */) dg1->setExternalConnectionNodes(1, dqbus1); dg1->setExternalConnectionNodes(2, dqbus1 + 1); //"grounding" of the difference - dg1->setExternalConnectionNodes(3, -1); + dg1->setExternalConnectionNodes(3, static_cast(-1)); // internal connections for (size_t i = 0; i < 12; i++) { diff --git a/examples/PowerElectronics/RLCircuit/RLCircuit.cpp b/examples/PowerElectronics/RLCircuit/RLCircuit.cpp index 1bea6597..4c80d0a8 100644 --- a/examples/PowerElectronics/RLCircuit/RLCircuit.cpp +++ b/examples/PowerElectronics/RLCircuit/RLCircuit.cpp @@ -37,7 +37,7 @@ int main(int /* argc */, char const** /* argv */) // input induct->setExternalConnectionNodes(0, 1); // output - induct->setExternalConnectionNodes(1, -1); + induct->setExternalConnectionNodes(1, static_cast(-1)); // internal induct->setExternalConnectionNodes(2, 2); // add component @@ -59,7 +59,7 @@ int main(int /* argc */, char const** /* argv */) GridKit::VoltageSource* vsource = new GridKit::VoltageSource(idoff, vinit); // Form index to node uid realations // input - vsource->setExternalConnectionNodes(0, -1); + vsource->setExternalConnectionNodes(0, static_cast(-1)); // output vsource->setExternalConnectionNodes(1, 0); // internal diff --git a/examples/PowerElectronics/ScaleMicrogrid/ScaleMicrogrid.cpp b/examples/PowerElectronics/ScaleMicrogrid/ScaleMicrogrid.cpp index 9864235a..ebc573d4 100644 --- a/examples/PowerElectronics/ScaleMicrogrid/ScaleMicrogrid.cpp +++ b/examples/PowerElectronics/ScaleMicrogrid/ScaleMicrogrid.cpp @@ -192,7 +192,7 @@ int test(index_type Nsize, real_type error_tol, bool debug_output) dg_ref->setExternalConnectionNodes(1, vdqbus_index[0]); dg_ref->setExternalConnectionNodes(2, vdqbus_index[0] + 1); //"grounding" of the difference - dg_ref->setExternalConnectionNodes(3, -1); + dg_ref->setExternalConnectionNodes(3, static_cast(-1)); // internal connections for (index_type i = 0; i < 12; i++) { diff --git a/src/LinearAlgebra/SparseMatrix/COO_Matrix.hpp b/src/LinearAlgebra/SparseMatrix/COO_Matrix.hpp index 33d28324..7e05ab25 100644 --- a/src/LinearAlgebra/SparseMatrix/COO_Matrix.hpp +++ b/src/LinearAlgebra/SparseMatrix/COO_Matrix.hpp @@ -189,14 +189,15 @@ inline std::vector COO_Matrix::getCSRRowData() { if (!this->isSorted()) this->sortSparse(); - std::vector row_size_vec(this->rows_size_ + 1, 0); + std::vector row_size_vec(static_cast(this->rows_size_ + 1), 0); IdxT counter = 0; for (IdxT i = 0; i < static_cast(row_size_vec.size() - 1); i++) { - row_size_vec[i + 1] = row_size_vec[i]; - while (counter < static_cast(this->row_indices_.size()) && i == this->row_indices_[counter]) + row_size_vec[static_cast(i + 1)] = row_size_vec[static_cast(i)]; + while (counter < static_cast(this->row_indices_.size()) + && i == this->row_indices_[static_cast(counter)]) { - row_size_vec[i + 1]++; + row_size_vec[static_cast(i + 1)]++; counter++; } } @@ -231,12 +232,12 @@ inline void COO_Matrix::setValues(std::vector r, std::vecto for (IdxT i = 0; i < static_cast(this->row_indices_.size()); i++) { // pushback values_ when they are not in current matrix - while (a_iter < static_cast(r.size()) && (r[a_iter] < this->row_indices_[i] || (r[a_iter] == this->row_indices_[i] && c[a_iter] < this->column_indices_[i]))) + while (a_iter < static_cast(r.size()) && (r[static_cast(a_iter)] < this->row_indices_[static_cast(i)] || (r[static_cast(a_iter)] == this->row_indices_[static_cast(i)] && c[static_cast(a_iter)] < this->column_indices_[static_cast(i)]))) { - this->row_indices_.push_back(r[a_iter]); - this->column_indices_.push_back(c[a_iter]); - this->values_.push_back(v[a_iter]); - this->checkIncreaseSize(r[a_iter], c[a_iter]); + this->row_indices_.push_back(r[static_cast(a_iter)]); + this->column_indices_.push_back(c[static_cast(a_iter)]); + this->values_.push_back(v[static_cast(a_iter)]); + this->checkIncreaseSize(r[static_cast(a_iter)], c[static_cast(a_iter)]); a_iter++; } if (a_iter >= static_cast(r.size())) @@ -244,20 +245,20 @@ inline void COO_Matrix::setValues(std::vector r, std::vecto break; } - if (r[a_iter] == this->row_indices_[i] && c[a_iter] == this->column_indices_[i]) + if (r[static_cast(a_iter)] == this->row_indices_[static_cast(i)] && c[static_cast(a_iter)] == this->column_indices_[static_cast(i)]) { - this->values_[i] = v[a_iter]; + this->values_[static_cast(i)] = v[static_cast(a_iter)]; a_iter++; } } // push back rest that was not found sorted for (IdxT i = a_iter; i < static_cast(r.size()); i++) { - this->row_indices_.push_back(r[i]); - this->column_indices_.push_back(c[i]); - this->values_.push_back(v[i]); + this->row_indices_.push_back(r[static_cast(i)]); + this->column_indices_.push_back(c[static_cast(i)]); + this->values_.push_back(v[static_cast(i)]); - this->checkIncreaseSize(r[i], c[i]); + this->checkIncreaseSize(r[static_cast(i)], c[static_cast(i)]); } this->sorted_ = false; @@ -304,13 +305,13 @@ inline void COO_Matrix::axpy(ScalarT alpha, COO_Matrix(this->row_indices_.size()); i++) { // pushback values when they are not in current matrix - while (a_iter < static_cast(r.size()) && (r[a_iter] < this->row_indices_[i] || (r[a_iter] == this->row_indices_[i] && c[a_iter] < this->column_indices_[i]))) + while (a_iter < static_cast(r.size()) && (r[static_cast(a_iter)] < this->row_indices_[static_cast(i)] || (r[static_cast(a_iter)] == this->row_indices_[static_cast(i)] && c[static_cast(a_iter)] < this->column_indices_[static_cast(i)]))) { - this->row_indices_.push_back(r[a_iter]); - this->column_indices_.push_back(c[a_iter]); - this->values_.push_back(alpha * val[a_iter]); + this->row_indices_.push_back(r[static_cast(a_iter)]); + this->column_indices_.push_back(c[static_cast(a_iter)]); + this->values_.push_back(alpha * val[static_cast(a_iter)]); - this->checkIncreaseSize(r[a_iter], c[a_iter]); + this->checkIncreaseSize(r[static_cast(a_iter)], c[static_cast(a_iter)]); a_iter++; } if (a_iter >= static_cast(r.size())) @@ -318,20 +319,20 @@ inline void COO_Matrix::axpy(ScalarT alpha, COO_Matrixrow_indices_[i] && c[a_iter] == this->column_indices_[i]) + if (r[static_cast(a_iter)] == this->row_indices_[static_cast(i)] && c[static_cast(a_iter)] == this->column_indices_[static_cast(i)]) { - this->values_[i] += alpha * val[a_iter]; + this->values_[static_cast(i)] += alpha * val[static_cast(a_iter)]; a_iter++; } } // push back rest that was not found sorted_ for (IdxT i = a_iter; i < static_cast(r.size()); i++) { - this->row_indices_.push_back(r[i]); - this->column_indices_.push_back(c[i]); - this->values_.push_back(alpha * val[i]); + this->row_indices_.push_back(r[static_cast(i)]); + this->column_indices_.push_back(c[static_cast(i)]); + this->values_.push_back(alpha * val[static_cast(i)]); - this->checkIncreaseSize(r[i], c[i]); + this->checkIncreaseSize(r[static_cast(i)], c[static_cast(i)]); } this->sorted_ = false; @@ -780,7 +781,9 @@ inline void COO_Matrix::sortSparseCOO(std::vector& rows, st std::sort(std::begin(ordervec), std::end(ordervec), [&](int i1, int i2) - { return (rows[i1] < rows[i2]) || (rows[i1] == rows[i2] && columns[i1] < columns[i2]); }); + { return (rows[static_cast(i1)] < rows[static_cast(i2)]) + || (rows[static_cast(i1)] == rows[static_cast(i2)] + && columns[static_cast(i1)] < columns[static_cast(i2)]); }); // reorder based of index-sorting. Only swap cost no extra memory. // @todo see if extra memory creation is fine diff --git a/src/Model/PhasorDynamics/SynchronousMachine/GENROUwS/Genrou.cpp b/src/Model/PhasorDynamics/SynchronousMachine/GENROUwS/Genrou.cpp index 1edd2dee..d55f38ca 100644 --- a/src/Model/PhasorDynamics/SynchronousMachine/GENROUwS/Genrou.cpp +++ b/src/Model/PhasorDynamics/SynchronousMachine/GENROUwS/Genrou.cpp @@ -1,413 +1,412 @@ -/** - * @file Genrou.cpp - * @author Adam Birchfield (abirchfield@tamu.edu) - * @author Slaven Peles (peless@ornl.gov) - * @brief Definition of a GENROU generator model. - * - * - */ - -#define _USE_MATH_DEFINES -#include "Genrou.hpp" - -#include -#include - -#include - -#define _USE_MATH_DEFINES - -namespace GridKit -{ - namespace PhasorDynamics - { - /*! - * @brief Constructor for a pi-model branch - * - * Arguments passed to ModelEvaluatorImpl: - * - Number of equations = 0 - * - Number of independent variables = 0 - * - Number of quadratures = 0 - * - Number of optimization parameters = 0 - */ - template - Genrou::Genrou(bus_type* bus, int unit_id) - : bus_(bus), - busID_(0), - unit_id_(unit_id), - p0_(0.), - q0_(0.), - H_(3.), - D_(0.), - Ra_(0.), - Tdop_(7.), - Tdopp_(.04), - Tqopp_(.05), - Tqop_(.75), - Xd_(2.1), - Xdp_(0.2), - Xdpp_(0.18), - Xq_(.5), - Xqp_(.5), - Xqpp_(.18), - Xl_(.15), - S10_(0.), - S12_(0.) - { - size_ = 21; - setDerivedParams(); - - // Temporary, to eliminate compiler warnings - (void) busID_; - (void) unit_id_; - } - - /*! - * @brief Constructor for a pi-model branch - * - * Arguments passed to ModelEvaluatorImpl: - * - Number of equations = 0 - * - Number of independent variables = 0 - * - Number of quadratures = 0 - * - Number of optimization parameters = 0 - */ - template - Genrou::Genrou(bus_type* bus, - int unit_id, - ScalarT p0, - ScalarT q0, - real_type H, - real_type D, - real_type Ra, - real_type Tdop, - real_type Tdopp, - real_type Tqopp, - real_type Tqop, - real_type Xd, - real_type Xdp, - real_type Xdpp, - real_type Xq, - real_type Xqp, - real_type Xqpp, - real_type Xl, - real_type S10, - real_type S12) - : bus_(bus), - busID_(0), - unit_id_(unit_id), - p0_(p0), - q0_(q0), - H_(H), - D_(D), - Ra_(Ra), - Tdop_(Tdop), - Tdopp_(Tdopp), - Tqopp_(Tqopp), - Tqop_(Tqop), - Xd_(Xd), - Xdp_(Xdp), - Xdpp_(Xdpp), - Xq_(Xq), - Xqp_(Xqp), - Xqpp_(Xqpp), - Xl_(Xl), - S10_(S10), - S12_(S12) - { - size_ = 21; - setDerivedParams(); - } - - // /** - // * @brief Destroy the Genrou - // * - // * @tparam ScalarT - // * @tparam IdxT - // */ - // template - // Genrou::~Genrou() - // { - // // std::cout << "Destroy Genrou..." << std::endl; - // } - - /*! - * @brief allocate method computes sparsity pattern of the Jacobian. - */ - template - int Genrou::allocate() - { - f_.resize(size_); - y_.resize(size_); - yp_.resize(size_); - tag_.resize(size_); - fB_.resize(size_); - yB_.resize(size_); - ypB_.resize(size_); - return 0; - } - - /** - * Initialization of the branch model - * - */ - template - int Genrou::initialize() - { - /* Initialization tricks -- assuming NO saturation */ - ScalarT vr = Vr(); - ScalarT vi = Vi(); - ScalarT p = p0_; - ScalarT q = q0_; - ScalarT vm2 = vr * vr + vi * vi; - ScalarT Er = vr + (Ra_ * p * vr + Ra_ * q * vi - Xq_ * p * vi + Xq_ * q * vr) / vm2; - ScalarT Ei = vi + (Ra_ * p * vi - Ra_ * q * vr + Xq_ * p * vr + Xq_ * q * vi) / vm2; - ScalarT delta = atan2(Ei, Er); - ScalarT omega = 0; - ScalarT ir = (p * vr + q * vi) / vm2; - ScalarT ii = (p * vi - q * vr) / vm2; - ScalarT id = ir * sin(delta) - ii * cos(delta); - ScalarT iq = ir * cos(delta) + ii * sin(delta); - ScalarT vd = vr * sin(delta) - vi * cos(delta) + id * Ra_ - iq * Xqpp_; - ScalarT vq = vr * cos(delta) + vi * sin(delta) + id * Xqpp_ - iq * Ra_; - ScalarT psiqpp = -vd / (1 + omega); - ScalarT psidpp = vq / (1 + omega); - ScalarT Te = (psidpp - id * Xdpp_) * iq - (psiqpp - iq * Xdpp_) * id; - ScalarT psiqp = -(-(Xqp_ - Xl_) * iq + psiqpp * (Xqp_ - Xl_) / (Xqpp_ - Xl_)) - / (1 + (Xqp_ - Xqpp_) / (Xqpp_ - Xl_)); - ScalarT Edp = psiqp - (Xqp_ - Xl_) * iq; - ScalarT psidp = -((Xdp_ - Xl_) * id - psidpp * (Xdp_ - Xl_) / (Xdpp_ - Xl_)) - / (1 + (Xdp_ - Xdpp_) / (Xdpp_ - Xl_)); - ScalarT Eqp = psidp + (Xdp_ - Xl_) * id; - - /* Now we have the state variables, solve for alg. variables */ - ScalarT ksat; - ScalarT psipp; - - y_[0] = delta; //= 0.55399038; - y_[1] = omega; // = 0; - y_[2] = Eqp; // = 0.995472581; - y_[3] = psidp; // = 0.971299567; - y_[4] = psiqp; // = 0.306880069; - y_[5] = Edp; // = 0; - - y_[6] = psiqpp = -psiqp * Xq4_ - Edp * Xq5_; - y_[7] = psidpp = psidp * Xd4_ + Eqp * Xd5_; - y_[8] = psipp = sqrt(psiqpp * psiqpp + psidpp * psidpp); - y_[9] = ksat = SB_ * pow(psipp - SA_, 2); - y_[10] = vd = -psiqpp * (1 + omega); - y_[11] = vq = psidpp * (1 + omega); - y_[12] = Te = (psidpp - id * Xdpp_) * iq - (psiqpp - iq * Xdpp_) * id; - y_[13] = id; - y_[14] = iq; - y_[15] = ir; - y_[16] = ii; - y_[17] = pmech_set_ = Te; - y_[18] = efd_set_ = Eqp + Xd1_ * (id + Xd3_ * (Eqp - psidp - Xd2_ * id)) + psidpp * ksat; - y_[19] = G_ * (vd * sin(delta) + vq * cos(delta)) - - B_ * (vd * -cos(delta) + vq * sin(delta)); /* inort, real */ - y_[20] = B_ * (vd * sin(delta) + vq * cos(delta)) - + G_ * (vd * -cos(delta) + vq * sin(delta)); /* inort, imag */ - - for (IdxT i = 0; i < size_; ++i) - yp_[i] = 0.0; - - return 0; - } - - /** - * \brief Identify differential variables. - */ - template - int Genrou::tagDifferentiable() - { - for (IdxT i = 0; i < size_; ++i) - { - tag_[i] = i < 6; - } - return 0; - } - - /** - * \brief Residual contribution of the branch is pushed to the - * two terminal buses. - * - */ - template - int Genrou::evaluateResidual() - { - /* Read variables */ - ScalarT delta = y_[0]; - ScalarT omega = y_[1]; - ScalarT Eqp = y_[2]; - ScalarT psidp = y_[3]; - ScalarT psiqp = y_[4]; - ScalarT Edp = y_[5]; - ScalarT psiqpp = y_[6]; - ScalarT psidpp = y_[7]; - ScalarT psipp = y_[8]; - ScalarT ksat = y_[9]; - ScalarT vd = y_[10]; - ScalarT vq = y_[11]; - ScalarT telec = y_[12]; - ScalarT id = y_[13]; - ScalarT iq = y_[14]; - ScalarT ir = y_[15]; - ScalarT ii = y_[16]; - ScalarT pmech = y_[17]; - ScalarT efd = y_[18]; - ScalarT inr = y_[19]; - ScalarT ini = y_[20]; - ScalarT vr = Vr(); - ScalarT vi = Vi(); - - /* Read derivatives */ - ScalarT delta_dot = yp_[0]; - ScalarT omega_dot = yp_[1]; - ScalarT Eqp_dot = yp_[2]; - ScalarT psidp_dot = yp_[3]; - ScalarT psiqp_dot = yp_[4]; - ScalarT Edp_dot = yp_[5]; - - /* 6 Genrou differential equations */ - f_[0] = delta_dot - omega * (2 * M_PI * 60); - f_[1] = omega_dot - (1 / (2 * H_)) * ((pmech - D_ * omega) / (1 + omega) - telec); - f_[2] = Eqp_dot - (1 / Tdop_) * (efd - (Eqp + Xd1_ * (id + Xd3_ * (Eqp - psidp - Xd2_ * id)) + psidpp * ksat)); - f_[3] = psidp_dot - (1 / Tdopp_) * (Eqp - psidp - Xd2_ * id); - f_[4] = psiqp_dot - (1 / Tqopp_) * (Edp - psiqp + Xq2_ * iq); - f_[5] = Edp_dot - (1 / Tqop_) * (-Edp + Xqd_ * psiqpp * ksat + Xq1_ * (iq - Xq3_ * (Edp + iq * Xq2_ - psiqp))); - - /* 11 Genrou algebraic equations */ - f_[6] = psiqpp - (-psiqp * Xq4_ - Edp * Xq5_); - f_[7] = psidpp - (psidp * Xd4_ + Eqp * Xd5_); - f_[8] = psipp - sqrt(pow(psidpp, 2.0) + pow(psiqpp, 2.0)); - f_[9] = ksat - SB_ * pow(psipp - SA_, 2.0); - f_[10] = vd + psiqpp * (1 + omega); - f_[11] = vq - psidpp * (1 + omega); - f_[12] = telec - ((psidpp - id * Xdpp_) * iq - (psiqpp - iq * Xdpp_) * id); - f_[13] = id - (ir * sin(delta) - ii * cos(delta)); - f_[14] = iq - (ir * cos(delta) + ii * sin(delta)); - f_[15] = ir + G_ * vr - B_ * vi - inr; - f_[16] = ii + B_ * vr + G_ * vi - ini; - - /* 2 Genrou control inputs are set to constant for this example */ - f_[17] = pmech - pmech_set_; - f_[18] = efd - efd_set_; - - /* 2 Genrou current source definitions */ - f_[19] = inr - (G_ * (sin(delta) * vd + cos(delta) * vq) - B_ * (-cos(delta) * vd + sin(delta) * vq)); - f_[20] = ini - (B_ * (sin(delta) * vd + cos(delta) * vq) + G_ * (-cos(delta) * vd + sin(delta) * vq)); - - /* Current balance */ - Ir() += inr - Vr() * G_ + Vi() * B_; - Ii() += ini - Vr() * B_ - Vi() * G_; - - return 0; - } - - /** - * @brief Jacobian evaluation not implemented yet - * - * @tparam ScalarT - scalar data type - * @tparam IdxT - matrix index data type - * @return int - error code, 0 = success - */ - template - int Genrou::evaluateJacobian() - { - std::cout << "Evaluate Jacobian for Genrou..." << std::endl; - std::cout << "Jacobian evaluation not implemented!" << std::endl; - return 0; - } - - /** - * @brief Integrand (objective) evaluation not implemented yet - * - * @tparam ScalarT - scalar data type - * @tparam IdxT - matrix index data type - * @return int - error code, 0 = success - */ - template - int Genrou::evaluateIntegrand() - { - // std::cout << "Evaluate Integrand for Genrou..." << std::endl; - return 0; - } - - /** - * @brief Adjoint initialization not implemented yet - * - * @tparam ScalarT - scalar data type - * @tparam IdxT - matrix index data type - * @return int - error code, 0 = success - */ - template - int Genrou::initializeAdjoint() - { - // std::cout << "Initialize adjoint for Genrou..." << std::endl; - return 0; - } - - /** - * @brief Adjoint residual evaluation not implemented yet - * - * @tparam ScalarT - scalar data type - * @tparam IdxT - matrix index data type - * @return int - error code, 0 = success - */ - template - int Genrou::evaluateAdjointResidual() - { - // std::cout << "Evaluate adjoint residual for Genrou..." << std::endl; - return 0; - } - - /** - * @brief Adjoint integrand (objective) evaluation not implemented yet - * - * @tparam ScalarT - scalar data type - * @tparam IdxT - matrix index data type - * @return int - error code, 0 = success - */ - template - int Genrou::evaluateAdjointIntegrand() - { - // std::cout << "Evaluate adjoint Integrand for Genrou..." << std::endl; - return 0; - } - - template - void Genrou::setDerivedParams() - { - SA_ = 0; - SB_ = 0; - if (S12_ != 0) - { - real_type s112 = sqrt(S10_ / S12_); - - SA_ = (1.2 * s112 + 1) / (s112 + 1); - SB_ = (1.2 * s112 - 1) / (s112 - 1); - if (SB_ < SA_) - SA_ = SB_; - SB_ = S12_ / pow(SA_ - 1.2, 2); - } - Xd1_ = Xd_ - Xdp_; - Xd2_ = Xdp_ - Xl_; - Xd3_ = (Xdp_ - Xdpp_) / (Xd2_ * Xd2_); - Xd4_ = (Xdp_ - Xdpp_) / Xd2_; - Xd5_ = (Xdpp_ - Xl_) / Xd2_; - Xq1_ = Xq_ - Xqp_; - Xq2_ = Xqp_ - Xl_; - Xq3_ = (Xqp_ - Xqpp_) / (Xq2_ * Xq2_); - Xq4_ = (Xqp_ - Xqpp_) / Xq2_; - Xq5_ = (Xqpp_ - Xl_) / Xq2_; - Xqd_ = (Xq_ - Xl_) / (Xd_ - Xl_); - G_ = Ra_ / (Ra_ * Ra_ + Xqpp_ * Xqpp_); - B_ = -Xqpp_ / (Ra_ * Ra_ + Xqpp_ * Xqpp_); - } - - // Available template instantiations - template class Genrou; - template class Genrou; - - } // namespace PhasorDynamics -} // namespace GridKit +/** + * @file Genrou.cpp + * @author Adam Birchfield (abirchfield@tamu.edu) + * @author Slaven Peles (peless@ornl.gov) + * @brief Definition of a GENROU generator model. + * + * + */ + +#define _USE_MATH_DEFINES +#include +#include + +#include "Genrou.hpp" +#include + +#define _USE_MATH_DEFINES + +namespace GridKit +{ + namespace PhasorDynamics + { + /*! + * @brief Constructor for a pi-model branch + * + * Arguments passed to ModelEvaluatorImpl: + * - Number of equations = 0 + * - Number of independent variables = 0 + * - Number of quadratures = 0 + * - Number of optimization parameters = 0 + */ + template + Genrou::Genrou(bus_type* bus, int unit_id) + : bus_(bus), + busID_(0), + unit_id_(unit_id), + p0_(0.), + q0_(0.), + H_(3.), + D_(0.), + Ra_(0.), + Tdop_(7.), + Tdopp_(.04), + Tqopp_(.05), + Tqop_(.75), + Xd_(2.1), + Xdp_(0.2), + Xdpp_(0.18), + Xq_(.5), + Xqp_(.5), + Xqpp_(.18), + Xl_(.15), + S10_(0.), + S12_(0.) + { + size_ = 21; + setDerivedParams(); + + // Temporary, to eliminate compiler warnings + (void) busID_; + (void) unit_id_; + } + + /*! + * @brief Constructor for a pi-model branch + * + * Arguments passed to ModelEvaluatorImpl: + * - Number of equations = 0 + * - Number of independent variables = 0 + * - Number of quadratures = 0 + * - Number of optimization parameters = 0 + */ + template + Genrou::Genrou(bus_type* bus, + int unit_id, + ScalarT p0, + ScalarT q0, + real_type H, + real_type D, + real_type Ra, + real_type Tdop, + real_type Tdopp, + real_type Tqopp, + real_type Tqop, + real_type Xd, + real_type Xdp, + real_type Xdpp, + real_type Xq, + real_type Xqp, + real_type Xqpp, + real_type Xl, + real_type S10, + real_type S12) + : bus_(bus), + busID_(0), + unit_id_(unit_id), + p0_(p0), + q0_(q0), + H_(H), + D_(D), + Ra_(Ra), + Tdop_(Tdop), + Tdopp_(Tdopp), + Tqopp_(Tqopp), + Tqop_(Tqop), + Xd_(Xd), + Xdp_(Xdp), + Xdpp_(Xdpp), + Xq_(Xq), + Xqp_(Xqp), + Xqpp_(Xqpp), + Xl_(Xl), + S10_(S10), + S12_(S12) + { + size_ = 21; + setDerivedParams(); + } + + // /** + // * @brief Destroy the Genrou + // * + // * @tparam ScalarT + // * @tparam IdxT + // */ + // template + // Genrou::~Genrou() + // { + // // std::cout << "Destroy Genrou..." << std::endl; + // } + + /*! + * @brief allocate method computes sparsity pattern of the Jacobian. + */ + template + int Genrou::allocate() + { + f_.resize(static_cast(size_)); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + tag_.resize(static_cast(size_)); + fB_.resize(static_cast(size_)); + yB_.resize(static_cast(size_)); + ypB_.resize(static_cast(size_)); + return 0; + } + + /** + * Initialization of the branch model + * + */ + template + int Genrou::initialize() + { + /* Initialization tricks -- assuming NO saturation */ + ScalarT vr = Vr(); + ScalarT vi = Vi(); + ScalarT p = p0_; + ScalarT q = q0_; + ScalarT vm2 = vr * vr + vi * vi; + ScalarT Er = vr + (Ra_ * p * vr + Ra_ * q * vi - Xq_ * p * vi + Xq_ * q * vr) / vm2; + ScalarT Ei = vi + (Ra_ * p * vi - Ra_ * q * vr + Xq_ * p * vr + Xq_ * q * vi) / vm2; + ScalarT delta = atan2(Ei, Er); + ScalarT omega = 0; + ScalarT ir = (p * vr + q * vi) / vm2; + ScalarT ii = (p * vi - q * vr) / vm2; + ScalarT id = ir * sin(delta) - ii * cos(delta); + ScalarT iq = ir * cos(delta) + ii * sin(delta); + ScalarT vd = vr * sin(delta) - vi * cos(delta) + id * Ra_ - iq * Xqpp_; + ScalarT vq = vr * cos(delta) + vi * sin(delta) + id * Xqpp_ - iq * Ra_; + ScalarT psiqpp = -vd / (1 + omega); + ScalarT psidpp = vq / (1 + omega); + ScalarT Te = (psidpp - id * Xdpp_) * iq - (psiqpp - iq * Xdpp_) * id; + ScalarT psiqp = -(-(Xqp_ - Xl_) * iq + psiqpp * (Xqp_ - Xl_) / (Xqpp_ - Xl_)) + / (1 + (Xqp_ - Xqpp_) / (Xqpp_ - Xl_)); + ScalarT Edp = psiqp - (Xqp_ - Xl_) * iq; + ScalarT psidp = -((Xdp_ - Xl_) * id - psidpp * (Xdp_ - Xl_) / (Xdpp_ - Xl_)) + / (1 + (Xdp_ - Xdpp_) / (Xdpp_ - Xl_)); + ScalarT Eqp = psidp + (Xdp_ - Xl_) * id; + + /* Now we have the state variables, solve for alg. variables */ + ScalarT ksat; + ScalarT psipp; + + y_[0] = delta; //= 0.55399038; + y_[1] = omega; // = 0; + y_[2] = Eqp; // = 0.995472581; + y_[3] = psidp; // = 0.971299567; + y_[4] = psiqp; // = 0.306880069; + y_[5] = Edp; // = 0; + + y_[6] = psiqpp = -psiqp * Xq4_ - Edp * Xq5_; + y_[7] = psidpp = psidp * Xd4_ + Eqp * Xd5_; + y_[8] = psipp = sqrt(psiqpp * psiqpp + psidpp * psidpp); + y_[9] = ksat = SB_ * pow(psipp - SA_, 2); + y_[10] = vd = -psiqpp * (1 + omega); + y_[11] = vq = psidpp * (1 + omega); + y_[12] = Te = (psidpp - id * Xdpp_) * iq - (psiqpp - iq * Xdpp_) * id; + y_[13] = id; + y_[14] = iq; + y_[15] = ir; + y_[16] = ii; + y_[17] = pmech_set_ = Te; + y_[18] = efd_set_ = Eqp + Xd1_ * (id + Xd3_ * (Eqp - psidp - Xd2_ * id)) + psidpp * ksat; + y_[19] = G_ * (vd * sin(delta) + vq * cos(delta)) + - B_ * (vd * -cos(delta) + vq * sin(delta)); /* inort, real */ + y_[20] = B_ * (vd * sin(delta) + vq * cos(delta)) + + G_ * (vd * -cos(delta) + vq * sin(delta)); /* inort, imag */ + + for (IdxT i = 0; i < size_; ++i) + yp_[static_cast(i)] = 0.0; + + return 0; + } + + /** + * \brief Identify differential variables. + */ + template + int Genrou::tagDifferentiable() + { + for (IdxT i = 0; i < size_; ++i) + { + tag_[static_cast(i)] = i < 6; + } + return 0; + } + + /** + * \brief Residual contribution of the branch is pushed to the + * two terminal buses. + * + */ + template + int Genrou::evaluateResidual() + { + /* Read variables */ + ScalarT delta = y_[0]; + ScalarT omega = y_[1]; + ScalarT Eqp = y_[2]; + ScalarT psidp = y_[3]; + ScalarT psiqp = y_[4]; + ScalarT Edp = y_[5]; + ScalarT psiqpp = y_[6]; + ScalarT psidpp = y_[7]; + ScalarT psipp = y_[8]; + ScalarT ksat = y_[9]; + ScalarT vd = y_[10]; + ScalarT vq = y_[11]; + ScalarT telec = y_[12]; + ScalarT id = y_[13]; + ScalarT iq = y_[14]; + ScalarT ir = y_[15]; + ScalarT ii = y_[16]; + ScalarT pmech = y_[17]; + ScalarT efd = y_[18]; + ScalarT inr = y_[19]; + ScalarT ini = y_[20]; + ScalarT vr = Vr(); + ScalarT vi = Vi(); + + /* Read derivatives */ + ScalarT delta_dot = yp_[0]; + ScalarT omega_dot = yp_[1]; + ScalarT Eqp_dot = yp_[2]; + ScalarT psidp_dot = yp_[3]; + ScalarT psiqp_dot = yp_[4]; + ScalarT Edp_dot = yp_[5]; + + /* 6 Genrou differential equations */ + f_[0] = delta_dot - omega * (2 * M_PI * 60); + f_[1] = omega_dot - (1 / (2 * H_)) * ((pmech - D_ * omega) / (1 + omega) - telec); + f_[2] = Eqp_dot - (1 / Tdop_) * (efd - (Eqp + Xd1_ * (id + Xd3_ * (Eqp - psidp - Xd2_ * id)) + psidpp * ksat)); + f_[3] = psidp_dot - (1 / Tdopp_) * (Eqp - psidp - Xd2_ * id); + f_[4] = psiqp_dot - (1 / Tqopp_) * (Edp - psiqp + Xq2_ * iq); + f_[5] = Edp_dot - (1 / Tqop_) * (-Edp + Xqd_ * psiqpp * ksat + Xq1_ * (iq - Xq3_ * (Edp + iq * Xq2_ - psiqp))); + + /* 11 Genrou algebraic equations */ + f_[6] = psiqpp - (-psiqp * Xq4_ - Edp * Xq5_); + f_[7] = psidpp - (psidp * Xd4_ + Eqp * Xd5_); + f_[8] = psipp - sqrt(pow(psidpp, 2.0) + pow(psiqpp, 2.0)); + f_[9] = ksat - SB_ * pow(psipp - SA_, 2.0); + f_[10] = vd + psiqpp * (1 + omega); + f_[11] = vq - psidpp * (1 + omega); + f_[12] = telec - ((psidpp - id * Xdpp_) * iq - (psiqpp - iq * Xdpp_) * id); + f_[13] = id - (ir * sin(delta) - ii * cos(delta)); + f_[14] = iq - (ir * cos(delta) + ii * sin(delta)); + f_[15] = ir + G_ * vr - B_ * vi - inr; + f_[16] = ii + B_ * vr + G_ * vi - ini; + + /* 2 Genrou control inputs are set to constant for this example */ + f_[17] = pmech - pmech_set_; + f_[18] = efd - efd_set_; + + /* 2 Genrou current source definitions */ + f_[19] = inr - (G_ * (sin(delta) * vd + cos(delta) * vq) - B_ * (-cos(delta) * vd + sin(delta) * vq)); + f_[20] = ini - (B_ * (sin(delta) * vd + cos(delta) * vq) + G_ * (-cos(delta) * vd + sin(delta) * vq)); + + /* Current balance */ + Ir() += inr - Vr() * G_ + Vi() * B_; + Ii() += ini - Vr() * B_ - Vi() * G_; + + return 0; + } + + /** + * @brief Jacobian evaluation not implemented yet + * + * @tparam ScalarT - scalar data type + * @tparam IdxT - matrix index data type + * @return int - error code, 0 = success + */ + template + int Genrou::evaluateJacobian() + { + std::cout << "Evaluate Jacobian for Genrou..." << std::endl; + std::cout << "Jacobian evaluation not implemented!" << std::endl; + return 0; + } + + /** + * @brief Integrand (objective) evaluation not implemented yet + * + * @tparam ScalarT - scalar data type + * @tparam IdxT - matrix index data type + * @return int - error code, 0 = success + */ + template + int Genrou::evaluateIntegrand() + { + // std::cout << "Evaluate Integrand for Genrou..." << std::endl; + return 0; + } + + /** + * @brief Adjoint initialization not implemented yet + * + * @tparam ScalarT - scalar data type + * @tparam IdxT - matrix index data type + * @return int - error code, 0 = success + */ + template + int Genrou::initializeAdjoint() + { + // std::cout << "Initialize adjoint for Genrou..." << std::endl; + return 0; + } + + /** + * @brief Adjoint residual evaluation not implemented yet + * + * @tparam ScalarT - scalar data type + * @tparam IdxT - matrix index data type + * @return int - error code, 0 = success + */ + template + int Genrou::evaluateAdjointResidual() + { + // std::cout << "Evaluate adjoint residual for Genrou..." << std::endl; + return 0; + } + + /** + * @brief Adjoint integrand (objective) evaluation not implemented yet + * + * @tparam ScalarT - scalar data type + * @tparam IdxT - matrix index data type + * @return int - error code, 0 = success + */ + template + int Genrou::evaluateAdjointIntegrand() + { + // std::cout << "Evaluate adjoint Integrand for Genrou..." << std::endl; + return 0; + } + + template + void Genrou::setDerivedParams() + { + SA_ = 0; + SB_ = 0; + if (S12_ != 0) + { + real_type s112 = sqrt(S10_ / S12_); + + SA_ = (1.2 * s112 + 1) / (s112 + 1); + SB_ = (1.2 * s112 - 1) / (s112 - 1); + if (SB_ < SA_) + SA_ = SB_; + SB_ = S12_ / pow(SA_ - 1.2, 2); + } + Xd1_ = Xd_ - Xdp_; + Xd2_ = Xdp_ - Xl_; + Xd3_ = (Xdp_ - Xdpp_) / (Xd2_ * Xd2_); + Xd4_ = (Xdp_ - Xdpp_) / Xd2_; + Xd5_ = (Xdpp_ - Xl_) / Xd2_; + Xq1_ = Xq_ - Xqp_; + Xq2_ = Xqp_ - Xl_; + Xq3_ = (Xqp_ - Xqpp_) / (Xq2_ * Xq2_); + Xq4_ = (Xqp_ - Xqpp_) / Xq2_; + Xq5_ = (Xqpp_ - Xl_) / Xq2_; + Xqd_ = (Xq_ - Xl_) / (Xd_ - Xl_); + G_ = Ra_ / (Ra_ * Ra_ + Xqpp_ * Xqpp_); + B_ = -Xqpp_ / (Ra_ * Ra_ + Xqpp_ * Xqpp_); + } + + // Available template instantiations + template class Genrou; + template class Genrou; + + } // namespace PhasorDynamics +} // namespace GridKit diff --git a/src/Model/PowerElectronics/Capacitor/Capacitor.cpp b/src/Model/PowerElectronics/Capacitor/Capacitor.cpp index 497c3b19..a7b0d5bb 100644 --- a/src/Model/PowerElectronics/Capacitor/Capacitor.cpp +++ b/src/Model/PowerElectronics/Capacitor/Capacitor.cpp @@ -39,9 +39,9 @@ namespace GridKit template int Capacitor::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/CircuitGraph.hpp b/src/Model/PowerElectronics/CircuitGraph.hpp index 26d5d6a2..5a57104f 100644 --- a/src/Model/PowerElectronics/CircuitGraph.hpp +++ b/src/Model/PowerElectronics/CircuitGraph.hpp @@ -105,7 +105,7 @@ size_t CircuitGraph::amountHyperEdges() */ template -void CircuitGraph::printBiPartiteGraph(bool verbose) +void CircuitGraph::printBiPartiteGraph([[maybe_unused]] bool verbose) { std::cout << "Amount of HyperNodes: " << this->amountHyperNodes() << std::endl; diff --git a/src/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.cpp b/src/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.cpp index 5be1802d..bbc16b52 100644 --- a/src/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.cpp +++ b/src/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.cpp @@ -57,9 +57,9 @@ namespace GridKit template int DistributedGenerator::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } @@ -186,7 +186,7 @@ namespace GridKit std::vector valsder(13, -1.0); for (int i = 0; i < 13; i++) { - rcordder[i] = i + 3; + rcordder[static_cast(i)] = static_cast(i + 3); } COO_Matrix Jacder = COO_Matrix(rcordder, rcordder, valsder, 16, 16); diff --git a/src/Model/PowerElectronics/InductionMotor/InductionMotor.cpp b/src/Model/PowerElectronics/InductionMotor/InductionMotor.cpp index 34e1f743..98cecf2b 100644 --- a/src/Model/PowerElectronics/InductionMotor/InductionMotor.cpp +++ b/src/Model/PowerElectronics/InductionMotor/InductionMotor.cpp @@ -51,9 +51,9 @@ namespace GridKit template int InductionMotor::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/Inductor/Inductor.cpp b/src/Model/PowerElectronics/Inductor/Inductor.cpp index fa9d0732..83a704f6 100644 --- a/src/Model/PowerElectronics/Inductor/Inductor.cpp +++ b/src/Model/PowerElectronics/Inductor/Inductor.cpp @@ -38,9 +38,9 @@ namespace GridKit int Inductor::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/LinearTransformer/LinearTransformer.cpp b/src/Model/PowerElectronics/LinearTransformer/LinearTransformer.cpp index 189dc1d6..600df629 100644 --- a/src/Model/PowerElectronics/LinearTransformer/LinearTransformer.cpp +++ b/src/Model/PowerElectronics/LinearTransformer/LinearTransformer.cpp @@ -53,9 +53,9 @@ namespace GridKit template int LinearTransformer::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.cpp b/src/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.cpp index 3de0b9e1..9706fd76 100644 --- a/src/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.cpp +++ b/src/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.cpp @@ -42,9 +42,9 @@ namespace GridKit template int MicrogridBusDQ::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/MicrogridLine/MicrogridLine.cpp b/src/Model/PowerElectronics/MicrogridLine/MicrogridLine.cpp index 030668ea..b5446f77 100644 --- a/src/Model/PowerElectronics/MicrogridLine/MicrogridLine.cpp +++ b/src/Model/PowerElectronics/MicrogridLine/MicrogridLine.cpp @@ -46,9 +46,9 @@ namespace GridKit template int MicrogridLine::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.cpp b/src/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.cpp index f36dd889..4a034191 100644 --- a/src/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.cpp +++ b/src/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.cpp @@ -44,9 +44,9 @@ namespace GridKit template int MicrogridLoad::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/Resistor/Resistor.cpp b/src/Model/PowerElectronics/Resistor/Resistor.cpp index 6549ec12..832297d0 100644 --- a/src/Model/PowerElectronics/Resistor/Resistor.cpp +++ b/src/Model/PowerElectronics/Resistor/Resistor.cpp @@ -37,9 +37,9 @@ namespace GridKit template int Resistor::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.cpp b/src/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.cpp index 9dc535cd..16b217ec 100644 --- a/src/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.cpp +++ b/src/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.cpp @@ -69,9 +69,9 @@ namespace GridKit template int SynchronousMachine::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/TransmissionLine/TransmissionLine.cpp b/src/Model/PowerElectronics/TransmissionLine/TransmissionLine.cpp index 85ec9151..399c1a8a 100644 --- a/src/Model/PowerElectronics/TransmissionLine/TransmissionLine.cpp +++ b/src/Model/PowerElectronics/TransmissionLine/TransmissionLine.cpp @@ -50,9 +50,9 @@ namespace GridKit template int TransmissionLine::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerElectronics/VoltageSource/VoltageSource.cpp b/src/Model/PowerElectronics/VoltageSource/VoltageSource.cpp index 8a2adaeb..8f45852c 100644 --- a/src/Model/PowerElectronics/VoltageSource/VoltageSource.cpp +++ b/src/Model/PowerElectronics/VoltageSource/VoltageSource.cpp @@ -37,9 +37,9 @@ namespace GridKit template int VoltageSource::allocate() { - y_.resize(size_); - yp_.resize(size_); - f_.resize(size_); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + f_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerFlow/Bus/BusPQ.cpp b/src/Model/PowerFlow/Bus/BusPQ.cpp index 207c304d..7b53a9c8 100644 --- a/src/Model/PowerFlow/Bus/BusPQ.cpp +++ b/src/Model/PowerFlow/Bus/BusPQ.cpp @@ -70,14 +70,14 @@ namespace GridKit int BusPQ::allocate() { // std::cout << "Allocate PQ bus ..." << std::endl; - f_.resize(size_); - y_.resize(size_); - yp_.resize(size_); - tag_.resize(size_); - - fB_.resize(size_); - yB_.resize(size_); - ypB_.resize(size_); + f_.resize(static_cast(size_)); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + tag_.resize(static_cast(size_)); + + fB_.resize(static_cast(size_)); + yB_.resize(static_cast(size_)); + ypB_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerFlow/Bus/BusPV.cpp b/src/Model/PowerFlow/Bus/BusPV.cpp index 1552fb06..c3312c70 100644 --- a/src/Model/PowerFlow/Bus/BusPV.cpp +++ b/src/Model/PowerFlow/Bus/BusPV.cpp @@ -68,14 +68,14 @@ namespace GridKit int BusPV::allocate() { // std::cout << "Allocate PV bus ..." << std::endl; - f_.resize(size_); - y_.resize(size_); - yp_.resize(size_); - tag_.resize(size_); - - fB_.resize(size_); - yB_.resize(size_); - ypB_.resize(size_); + f_.resize(static_cast(size_)); + y_.resize(static_cast(size_)); + yp_.resize(static_cast(size_)); + tag_.resize(static_cast(size_)); + + fB_.resize(static_cast(size_)); + yB_.resize(static_cast(size_)); + ypB_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerFlow/Generator2/Generator2.cpp b/src/Model/PowerFlow/Generator2/Generator2.cpp index 2d6b6538..454df1dc 100644 --- a/src/Model/PowerFlow/Generator2/Generator2.cpp +++ b/src/Model/PowerFlow/Generator2/Generator2.cpp @@ -48,7 +48,7 @@ namespace GridKit template int Generator2::allocate() { - tag_.resize(size_); + tag_.resize(static_cast(size_)); return 0; } diff --git a/src/Model/PowerFlow/Generator4/Generator4.cpp b/src/Model/PowerFlow/Generator4/Generator4.cpp index f5d5c2b1..e336a30e 100644 --- a/src/Model/PowerFlow/Generator4/Generator4.cpp +++ b/src/Model/PowerFlow/Generator4/Generator4.cpp @@ -1,10 +1,9 @@ #define _USE_MATH_DEFINES -#include "Generator4.hpp" - #include #include +#include "Generator4.hpp" #include namespace GridKit @@ -57,7 +56,7 @@ namespace GridKit int Generator4::allocate() { // std::cout << "Allocate Generator4..." << std::endl; - tag_.resize(size_); + tag_.resize(static_cast(size_)); return 0; } @@ -143,7 +142,7 @@ namespace GridKit for (IdxT i = 4; i < size_; ++i) { - tag_[i] = false; + tag_[static_cast(i)] = false; } return 0; @@ -217,8 +216,8 @@ namespace GridKit // std::cout << "Initialize adjoint for Generator4..." << std::endl; for (IdxT i = 0; i < size_; ++i) { - yB_[i] = 0.0; - ypB_[i] = 0.0; + yB_[static_cast(i)] = 0.0; + ypB_[static_cast(i)] = 0.0; } ypB_[1] = frequencyPenaltyDer(y_[1]); diff --git a/src/Model/PowerFlow/Generator4Governor/Generator4Governor.cpp b/src/Model/PowerFlow/Generator4Governor/Generator4Governor.cpp index 3736eb77..cc65f32e 100644 --- a/src/Model/PowerFlow/Generator4Governor/Generator4Governor.cpp +++ b/src/Model/PowerFlow/Generator4Governor/Generator4Governor.cpp @@ -1,10 +1,9 @@ #define _USE_MATH_DEFINES -#include "Generator4Governor.hpp" - #include #include +#include "Generator4Governor.hpp" #include "Model/PowerFlow/Bus/BaseBus.hpp" namespace GridKit @@ -66,7 +65,7 @@ namespace GridKit int Generator4Governor::allocate() { // std::cout << "Allocate Gen2..." << std::endl; - tag_.resize(size_); + tag_.resize(static_cast(size_)); return 0; } @@ -111,29 +110,29 @@ namespace GridKit const ScalarT Eqp = V() * cos(delta - theta()) + Rs_ * Iq + Xdp_ * Id; // Initialize generator - y_[offsetGen_ + 0] = delta; - y_[offsetGen_ + 1] = omega_s_ + 0.2; // <~ this is hack to perturb omega - y_[offsetGen_ + 2] = Edp; - y_[offsetGen_ + 3] = Eqp; - y_[offsetGen_ + 4] = Id; - y_[offsetGen_ + 5] = Iq; - yp_[offsetGen_ + 0] = 0.0; - yp_[offsetGen_ + 1] = 0.0; - yp_[offsetGen_ + 2] = 0.0; - yp_[offsetGen_ + 3] = 0.0; - yp_[offsetGen_ + 4] = 0.0; - yp_[offsetGen_ + 5] = 0.0; + y_[static_cast(offsetGen_ + 0)] = delta; + y_[static_cast(offsetGen_ + 1)] = omega_s_ + 0.2; // <~ this is hack to perturb omega + y_[static_cast(offsetGen_ + 2)] = Edp; + y_[static_cast(offsetGen_ + 3)] = Eqp; + y_[static_cast(offsetGen_ + 4)] = Id; + y_[static_cast(offsetGen_ + 5)] = Iq; + yp_[static_cast(offsetGen_ + 0)] = 0.0; + yp_[static_cast(offsetGen_ + 1)] = 0.0; + yp_[static_cast(offsetGen_ + 2)] = 0.0; + yp_[static_cast(offsetGen_ + 3)] = 0.0; + yp_[static_cast(offsetGen_ + 4)] = 0.0; + yp_[static_cast(offsetGen_ + 5)] = 0.0; Pm0_ = Edp * Id + Eqp * Iq + (Xqp_ - Xdp_) * Id * Iq; // <~ set to steady state value Ef0_ = Eqp + (Xd_ - Xdp_) * Id; // <~ set to steady state value // Initialize governor - y_[offsetGov_ + 0] = Pm0_; - y_[offsetGov_ + 1] = 0.0; - y_[offsetGov_ + 2] = 0.0; - yp_[offsetGov_ + 0] = 0.0; - yp_[offsetGov_ + 1] = 0.0; - yp_[offsetGov_ + 2] = 0.0; + y_[static_cast(offsetGov_ + 0)] = Pm0_; + y_[static_cast(offsetGov_ + 1)] = 0.0; + y_[static_cast(offsetGov_ + 2)] = 0.0; + yp_[static_cast(offsetGov_ + 0)] = 0.0; + yp_[static_cast(offsetGov_ + 1)] = 0.0; + yp_[static_cast(offsetGov_ + 2)] = 0.0; param_[1] = K_; param_up_[1] = 20.0; @@ -153,16 +152,16 @@ namespace GridKit int Generator4Governor::tagDifferentiable() { // std::cout << "size of tag vector is " << tag_.size() << "\n"; - tag_[offsetGen_ + 0] = true; - tag_[offsetGen_ + 1] = true; - tag_[offsetGen_ + 2] = true; - tag_[offsetGen_ + 3] = true; - tag_[offsetGen_ + 4] = false; - tag_[offsetGen_ + 5] = false; + tag_[static_cast(offsetGen_ + 0)] = true; + tag_[static_cast(offsetGen_ + 1)] = true; + tag_[static_cast(offsetGen_ + 2)] = true; + tag_[static_cast(offsetGen_ + 3)] = true; + tag_[static_cast(offsetGen_ + 4)] = false; + tag_[static_cast(offsetGen_ + 5)] = false; - tag_[offsetGov_ + 0] = true; - tag_[offsetGov_ + 1] = true; - tag_[offsetGov_ + 2] = false; + tag_[static_cast(offsetGov_ + 0)] = true; + tag_[static_cast(offsetGov_ + 1)] = true; + tag_[static_cast(offsetGov_ + 2)] = false; return 0; } @@ -207,21 +206,21 @@ namespace GridKit int Generator4Governor::evaluateResidual() { // Generator equations - f_[offsetGen_ + 0] = dotDelta() - omega_b_ * (omega() - omega_s_); - f_[offsetGen_ + 1] = (2.0 * H_) / omega_s_ * dotOmega() - Lm(y_[offsetGov_ + 0]) + Eqp() * Iq() + Edp() * Id() + (-Xdp_ + Xqp_) * Id() * Iq() + D_ * (omega() - omega_s_); - f_[offsetGen_ + 2] = Tq0p_ * dotEdp() + Edp() - (Xq_ - Xqp_) * Iq(); - f_[offsetGen_ + 3] = Td0p_ * dotEqp() + Eqp() + (Xd_ - Xdp_) * Id() - Ef0_; - f_[offsetGen_ + 4] = Rs_ * Id() - Xqp_ * Iq() + V() * sin(delta() - theta()) - Edp(); - f_[offsetGen_ + 5] = Xdp_ * Id() + Rs_ * Iq() + V() * cos(delta() - theta()) - Eqp(); + f_[static_cast(offsetGen_ + 0)] = dotDelta() - omega_b_ * (omega() - omega_s_); + f_[static_cast(offsetGen_ + 1)] = (2.0 * H_) / omega_s_ * dotOmega() - Lm(y_[static_cast(offsetGov_ + 0)]) + Eqp() * Iq() + Edp() * Id() + (-Xdp_ + Xqp_) * Id() * Iq() + D_ * (omega() - omega_s_); + f_[static_cast(offsetGen_ + 2)] = Tq0p_ * dotEdp() + Edp() - (Xq_ - Xqp_) * Iq(); + f_[static_cast(offsetGen_ + 3)] = Td0p_ * dotEqp() + Eqp() + (Xd_ - Xdp_) * Id() - Ef0_; + f_[static_cast(offsetGen_ + 4)] = Rs_ * Id() - Xqp_ * Iq() + V() * sin(delta() - theta()) - Edp(); + f_[static_cast(offsetGen_ + 5)] = Xdp_ * Id() + Rs_ * Iq() + V() * cos(delta() - theta()) - Eqp(); // Bus equations P() += Pg(); Q() += Qg(); // Governor equations - f_[offsetGov_ + 0] = yp_[offsetGov_ + 0] - Ln(y_[offsetGov_ + 2]); - f_[offsetGov_ + 1] = T1() * yp_[offsetGov_ + 1] + y_[offsetGov_ + 1] - (1.0 - T2() / T1()) * (omega() - omega_s_); - f_[offsetGov_ + 2] = T3() * y_[offsetGov_ + 2] - Pm0_ + Lm(y_[offsetGov_ + 0]) + K() * y_[offsetGov_ + 1] + K() * T2() / T1() * (omega() - omega_s_); + f_[static_cast(offsetGov_ + 0)] = yp_[static_cast(offsetGov_ + 0)] - Ln(y_[static_cast(offsetGov_ + 2)]); + f_[static_cast(offsetGov_ + 1)] = T1() * yp_[static_cast(offsetGov_ + 1)] + y_[static_cast(offsetGov_ + 1)] - (1.0 - T2() / T1()) * (omega() - omega_s_); + f_[static_cast(offsetGov_ + 2)] = T3() * y_[static_cast(offsetGov_ + 2)] - Pm0_ + Lm(y_[static_cast(offsetGov_ + 0)]) + K() * y_[static_cast(offsetGov_ + 1)] + K() * T2() / T1() * (omega() - omega_s_); return 0; } @@ -253,10 +252,10 @@ namespace GridKit // std::cout << "Initialize adjoint for Generator4Governor..." << std::endl; for (IdxT i = 0; i < size_; ++i) { - yB_[i] = 0.0; - ypB_[i] = 0.0; + yB_[static_cast(i)] = 0.0; + ypB_[static_cast(i)] = 0.0; } - ypB_[offsetGen_ + 1] = frequencyPenaltyDer(omega()); + ypB_[static_cast(offsetGen_ + 1)] = frequencyPenaltyDer(omega()); return 0; } @@ -296,26 +295,26 @@ namespace GridKit ScalarT cosPhi = cos(delta() - theta()); // Generator adjoint - fB_[offsetGen_ + 0] = ypB_[offsetGen_ + 0] - yB_[offsetGen_ + 4] * V() * cosPhi + yB_[offsetGen_ + 5] * V() * sinPhi; - fB_[offsetGen_ + 1] = 2.0 * H_ / omega_s_ * ypB_[offsetGen_ + 1] + yB_[offsetGen_ + 0] * omega_b_ - yB_[offsetGen_ + 1] * D_ + frequencyPenaltyDer(omega()) - + yB_[offsetGov_ + 1] * (1.0 - T2() / T1()) - yB_[offsetGov_ + 2] * K() * T2() / T1(); - fB_[offsetGen_ + 2] = Tq0p_ * ypB_[offsetGen_ + 2] - yB_[offsetGen_ + 1] * Id() - yB_[offsetGen_ + 2] + yB_[offsetGen_ + 4] + fB_[static_cast(offsetGen_ + 0)] = ypB_[static_cast(offsetGen_ + 0)] - yB_[static_cast(offsetGen_ + 4)] * V() * cosPhi + yB_[static_cast(offsetGen_ + 5)] * V() * sinPhi; + fB_[static_cast(offsetGen_ + 1)] = 2.0 * H_ / omega_s_ * ypB_[static_cast(offsetGen_ + 1)] + yB_[static_cast(offsetGen_ + 0)] * omega_b_ - yB_[static_cast(offsetGen_ + 1)] * D_ + frequencyPenaltyDer(omega()) + + yB_[static_cast(offsetGov_ + 1)] * (1.0 - T2() / T1()) - yB_[static_cast(offsetGov_ + 2)] * K() * T2() / T1(); + fB_[static_cast(offsetGen_ + 2)] = Tq0p_ * ypB_[static_cast(offsetGen_ + 2)] - yB_[static_cast(offsetGen_ + 1)] * Id() - yB_[static_cast(offsetGen_ + 2)] + yB_[static_cast(offsetGen_ + 4)] + lambdaP() * Id() - lambdaQ() * Iq(); - fB_[offsetGen_ + 3] = Td0p_ * ypB_[offsetGen_ + 3] - yB_[offsetGen_ + 1] * Iq() - yB_[offsetGen_ + 3] + yB_[offsetGen_ + 5] + fB_[static_cast(offsetGen_ + 3)] = Td0p_ * ypB_[static_cast(offsetGen_ + 3)] - yB_[static_cast(offsetGen_ + 1)] * Iq() - yB_[static_cast(offsetGen_ + 3)] + yB_[static_cast(offsetGen_ + 5)] + lambdaP() * Iq() + lambdaQ() * Id(); - fB_[offsetGen_ + 4] = -yB_[offsetGen_ + 1] * (Edp() + (Xqp_ - Xdp_) * Iq()) - yB_[offsetGen_ + 3] * (Xd_ - Xdp_) - yB_[offsetGen_ + 4] * Rs_ - yB_[offsetGen_ + 5] * Xdp_ + fB_[static_cast(offsetGen_ + 4)] = -yB_[static_cast(offsetGen_ + 1)] * (Edp() + (Xqp_ - Xdp_) * Iq()) - yB_[static_cast(offsetGen_ + 3)] * (Xd_ - Xdp_) - yB_[static_cast(offsetGen_ + 4)] * Rs_ - yB_[static_cast(offsetGen_ + 5)] * Xdp_ + lambdaP() * (Edp() + (Xqp_ - Xdp_) * Iq() - 2.0 * Rs_ * Id()) + lambdaQ() * (Eqp() - 2.0 * Xdp_ * Id()); - fB_[offsetGen_ + 5] = -yB_[offsetGen_ + 1] * (Eqp() + (Xqp_ - Xdp_) * Id()) + yB_[offsetGen_ + 2] * (Xq_ - Xqp_) + yB_[offsetGen_ + 4] * Xqp_ - yB_[offsetGen_ + 5] * Rs_ + fB_[static_cast(offsetGen_ + 5)] = -yB_[static_cast(offsetGen_ + 1)] * (Eqp() + (Xqp_ - Xdp_) * Id()) + yB_[static_cast(offsetGen_ + 2)] * (Xq_ - Xqp_) + yB_[static_cast(offsetGen_ + 4)] * Xqp_ - yB_[static_cast(offsetGen_ + 5)] * Rs_ + lambdaP() * (Eqp() + (Xqp_ - Xdp_) * Id() - 2.0 * Rs_ * Iq()) - lambdaQ() * (Edp() + 2.0 * Xqp_ * Iq()); // Bus adjoint - PB() += (-yB_[offsetGen_ + 4] * sinPhi - yB_[offsetGen_ + 5] * cosPhi); - QB() += (yB_[offsetGen_ + 4] * V() * cosPhi - yB_[offsetGen_ + 5] * V() * sinPhi); + PB() += (-yB_[static_cast(offsetGen_ + 4)] * sinPhi - yB_[static_cast(offsetGen_ + 5)] * cosPhi); + QB() += (yB_[static_cast(offsetGen_ + 4)] * V() * cosPhi - yB_[static_cast(offsetGen_ + 5)] * V() * sinPhi); // Governor adjoint - fB_[offsetGov_ + 0] = ypB_[offsetGov_ + 0] - yB_[offsetGov_ + 2] * dLm(y_[offsetGov_ + 0]) + yB_[offsetGen_ + 1] * dLm(y_[offsetGov_ + 0]); - fB_[offsetGov_ + 1] = ypB_[offsetGov_ + 1] * T1() - yB_[offsetGov_ + 1] - yB_[offsetGov_ + 2] * K(); - fB_[offsetGov_ + 2] = yB_[offsetGov_ + 0] * dLn(y_[offsetGov_ + 2]) - yB_[offsetGov_ + 2] * T3(); + fB_[static_cast(offsetGov_ + 0)] = ypB_[static_cast(offsetGov_ + 0)] - yB_[static_cast(offsetGov_ + 2)] * dLm(y_[static_cast(offsetGov_ + 0)]) + yB_[static_cast(offsetGen_ + 1)] * dLm(y_[static_cast(offsetGov_ + 0)]); + fB_[static_cast(offsetGov_ + 1)] = ypB_[static_cast(offsetGov_ + 1)] * T1() - yB_[static_cast(offsetGov_ + 1)] - yB_[static_cast(offsetGov_ + 2)] * K(); + fB_[static_cast(offsetGov_ + 2)] = yB_[static_cast(offsetGov_ + 0)] * dLn(y_[static_cast(offsetGov_ + 2)]) - yB_[static_cast(offsetGov_ + 2)] * T3(); return 0; } @@ -334,10 +333,10 @@ namespace GridKit // std::cout << "Evaluate adjoint Integrand for Gen2..." << std::endl; // K adjoint - gB_[1] = -yB_[offsetGov_ + 2] * (y_[offsetGov_ + 1] + T2() / T1() * (omega() - omega_s_)); + gB_[1] = -yB_[static_cast(offsetGov_ + 2)] * (y_[static_cast(offsetGov_ + 1)] + T2() / T1() * (omega() - omega_s_)); // T2 adjoint - gB_[0] = -yB_[offsetGov_ + 1] * (omega() - omega_s_) / T1() - yB_[offsetGov_ + 2] * K() / T1() * (omega() - omega_s_); + gB_[0] = -yB_[static_cast(offsetGov_ + 1)] * (omega() - omega_s_) / T1() - yB_[static_cast(offsetGov_ + 2)] * K() / T1() * (omega() - omega_s_); return 0; } diff --git a/src/Model/PowerFlow/Generator4Governor/Generator4Governor.hpp b/src/Model/PowerFlow/Generator4Governor/Generator4Governor.hpp index 948f6f67..5464526a 100644 --- a/src/Model/PowerFlow/Generator4Governor/Generator4Governor.hpp +++ b/src/Model/PowerFlow/Generator4Governor/Generator4Governor.hpp @@ -151,52 +151,52 @@ namespace GridKit const ScalarT dotDelta() const { - return yp_[offsetGen_ + 0]; + return yp_[static_cast(offsetGen_ + 0)]; } const ScalarT dotOmega() const { - return yp_[offsetGen_ + 1]; + return yp_[static_cast(offsetGen_ + 1)]; } const ScalarT dotEdp() const { - return yp_[offsetGen_ + 2]; + return yp_[static_cast(offsetGen_ + 2)]; } const ScalarT dotEqp() const { - return yp_[offsetGen_ + 3]; + return yp_[static_cast(offsetGen_ + 3)]; } const ScalarT delta() const { - return y_[offsetGen_ + 0]; + return y_[static_cast(offsetGen_ + 0)]; } const ScalarT omega() const { - return y_[offsetGen_ + 1]; + return y_[static_cast(offsetGen_ + 1)]; } const ScalarT Edp() const { - return y_[offsetGen_ + 2]; + return y_[static_cast(offsetGen_ + 2)]; } const ScalarT Eqp() const { - return y_[offsetGen_ + 3]; + return y_[static_cast(offsetGen_ + 3)]; } const ScalarT Id() const { - return y_[offsetGen_ + 4]; + return y_[static_cast(offsetGen_ + 4)]; } const ScalarT Iq() const { - return y_[offsetGen_ + 5]; + return y_[static_cast(offsetGen_ + 5)]; } const ScalarT K() const diff --git a/src/Model/PowerFlow/Generator4Param/Generator4Param.cpp b/src/Model/PowerFlow/Generator4Param/Generator4Param.cpp index b1dd73f2..c7b81b10 100644 --- a/src/Model/PowerFlow/Generator4Param/Generator4Param.cpp +++ b/src/Model/PowerFlow/Generator4Param/Generator4Param.cpp @@ -1,10 +1,9 @@ #define _USE_MATH_DEFINES -#include "Generator4Param.hpp" - #include #include +#include "Generator4Param.hpp" #include namespace GridKit @@ -53,7 +52,7 @@ namespace GridKit int Generator4Param::allocate() { // std::cout << "Allocate Generator4Param..." << std::endl; - tag_.resize(size_); + tag_.resize(static_cast(size_)); return 0; } @@ -143,7 +142,7 @@ namespace GridKit for (IdxT i = 4; i < size_; ++i) { - tag_[i] = false; + tag_[static_cast(i)] = false; } return 0; @@ -219,8 +218,8 @@ namespace GridKit // std::cout << "Initialize adjoint for Generator4Param..." << std::endl; for (IdxT i = 0; i < size_; ++i) { - yB_[i] = 0.0; - ypB_[i] = 0.0; + yB_[static_cast(i)] = 0.0; + ypB_[static_cast(i)] = 0.0; } // ypB_[1] = frequencyPenaltyDer(y_[1]); ypB_[2] = -trajectoryPenaltyDerEdp(time_) / Tq0p_; @@ -325,16 +324,16 @@ namespace GridKit if (t >= ti && t < tf) { // Interpolate from look-up table - Edp_est = (table_[n + 1][3] - table_[n][3]) / (table_[n + 1][0] - table_[n][0]) * (t - table_[n][0]) + table_[n][3]; - Eqp_est = (table_[n + 1][4] - table_[n][4]) / (table_[n + 1][0] - table_[n][0]) * (t - table_[n][0]) + table_[n][4]; + Edp_est = (table_[static_cast(n + 1)][3] - table_[static_cast(n)][3]) / (table_[static_cast(n + 1)][0] - table_[static_cast(n)][0]) * (t - table_[static_cast(n)][0]) + table_[static_cast(n)][3]; + Eqp_est = (table_[static_cast(n + 1)][4] - table_[static_cast(n)][4]) / (table_[static_cast(n + 1)][0] - table_[static_cast(n)][0]) * (t - table_[static_cast(n)][0]) + table_[static_cast(n)][4]; } else { if (tf <= t && t < tf + dt) { // Extrapolate from look-up table - Edp_est = (table_[n][3] - table_[n - 1][3]) / (table_[n][0] - table_[n - 1][0]) * (t - table_[n - 1][0]) + table_[n - 1][3]; - Eqp_est = (table_[n][4] - table_[n - 1][4]) / (table_[n][0] - table_[n - 1][0]) * (t - table_[n - 1][0]) + table_[n - 1][4]; + Edp_est = (table_[static_cast(n)][3] - table_[static_cast(n - 1)][3]) / (table_[static_cast(n)][0] - table_[static_cast(n - 1)][0]) * (t - table_[static_cast(n - 1)][0]) + table_[static_cast(n - 1)][3]; + Eqp_est = (table_[static_cast(n)][4] - table_[static_cast(n - 1)][4]) / (table_[static_cast(n)][0] - table_[static_cast(n - 1)][0]) * (t - table_[static_cast(n - 1)][0]) + table_[static_cast(n - 1)][4]; } else { @@ -360,13 +359,13 @@ namespace GridKit if (t >= ti && t < tf) { - Edp_est = (table_[n + 1][3] - table_[n][3]) / (table_[n + 1][0] - table_[n][0]) * (t - table_[n][0]) + table_[n][3]; + Edp_est = (table_[static_cast(n + 1)][3] - table_[static_cast(n)][3]) / (table_[static_cast(n + 1)][0] - table_[static_cast(n)][0]) * (t - table_[static_cast(n)][0]) + table_[static_cast(n)][3]; } else { if (tf <= t && t < tf + dt) { - Edp_est = (table_[n][3] - table_[n - 1][3]) / (table_[n][0] - table_[n - 1][0]) * (t - table_[n - 1][0]) + table_[n - 1][3]; + Edp_est = (table_[static_cast(n)][3] - table_[static_cast(n - 1)][3]) / (table_[static_cast(n)][0] - table_[static_cast(n - 1)][0]) * (t - table_[static_cast(n - 1)][0]) + table_[static_cast(n - 1)][3]; } else { @@ -391,13 +390,13 @@ namespace GridKit if (t >= ti && t < tf) { - Eqp_est = (table_[n + 1][4] - table_[n][4]) / (table_[n + 1][0] - table_[n][0]) * (t - table_[n][0]) + table_[n][4]; + Eqp_est = (table_[static_cast(n + 1)][4] - table_[static_cast(n)][4]) / (table_[static_cast(n + 1)][0] - table_[static_cast(n)][0]) * (t - table_[static_cast(n)][0]) + table_[static_cast(n)][4]; } else { if (tf <= t && t < tf + dt) { - Eqp_est = (table_[n][4] - table_[n - 1][4]) / (table_[n][0] - table_[n - 1][0]) * (t - table_[n - 1][0]) + table_[n - 1][4]; + Eqp_est = (table_[static_cast(n)][4] - table_[static_cast(n - 1)][4]) / (table_[static_cast(n)][0] - table_[static_cast(n - 1)][0]) * (t - table_[static_cast(n - 1)][0]) + table_[static_cast(n - 1)][4]; } else { diff --git a/src/Model/PowerFlow/ModelEvaluatorImpl.hpp b/src/Model/PowerFlow/ModelEvaluatorImpl.hpp index 3b7bb482..a8d26fba 100644 --- a/src/Model/PowerFlow/ModelEvaluatorImpl.hpp +++ b/src/Model/PowerFlow/ModelEvaluatorImpl.hpp @@ -30,18 +30,18 @@ namespace GridKit : size_(size), size_quad_(size_quad), size_opt_(size_opt), - y_(size_), - yp_(size_), - f_(size_), - g_(size_quad_), - yB_(size_), - ypB_(size_), - fB_(size_), - gB_(size_opt_), + y_(static_cast(size_)), + yp_(static_cast(size_)), + f_(static_cast(size_)), + g_(static_cast(size_quad_)), + yB_(static_cast(size_)), + ypB_(static_cast(size_)), + fB_(static_cast(size_)), + gB_(static_cast(size_opt_)), jac_(COO_Matrix()), - param_(size_opt_), - param_up_(size_opt_), - param_lo_(size_opt_) + param_(static_cast(size_opt_)), + param_up_(static_cast(size_opt_)), + param_lo_(static_cast(size_opt_)) { } diff --git a/src/Solver/Dynamic/Ida.cpp b/src/Solver/Dynamic/Ida.cpp index 5bde5db4..b8b4fb6e 100644 --- a/src/Solver/Dynamic/Ida.cpp +++ b/src/Solver/Dynamic/Ida.cpp @@ -1,12 +1,11 @@ -#include "Ida.hpp" - #include #include #include #include +#include "Ida.hpp" #include "Model/Evaluator.hpp" namespace AnalysisManager @@ -51,7 +50,7 @@ namespace AnalysisManager int retval = 0; // Allocate solution vectors - yy_ = N_VNew_Serial(model_->size(), context_); + yy_ = N_VNew_Serial(static_cast(model_->size()), context_); checkAllocation((void*) yy_, "N_VNew_Serial"); yp_ = N_VClone(yy_); checkAllocation((void*) yp_, "N_VClone"); @@ -88,7 +87,7 @@ namespace AnalysisManager model_->setMaxSteps(msa); /// \todo Need to set max number of steps based on user input! - retval = IDASetMaxNumSteps(solver_, msa); + retval = IDASetMaxNumSteps(solver_, static_cast(msa)); checkOutput(retval, "IDASetMaxNumSteps"); // Tag differential variables @@ -116,7 +115,11 @@ namespace AnalysisManager int retval = 0; if (model_->hasJacobian()) { - JacobianMat_ = SUNSparseMatrix(model_->size(), model_->size(), model_->size() * model_->size(), CSR_MAT, context_); + JacobianMat_ = SUNSparseMatrix(static_cast(model_->size()), + static_cast(model_->size()), + static_cast(model_->size() * model_->size()), + CSR_MAT, + context_); checkAllocation((void*) JacobianMat_, "SUNSparseMatrix"); linearSolver_ = SUNLinSol_KLU(yy_, JacobianMat_, context_); @@ -130,7 +133,9 @@ namespace AnalysisManager } else { - JacobianMat_ = SUNDenseMatrix(model_->size(), model_->size(), context_); + JacobianMat_ = SUNDenseMatrix(static_cast(model_->size()), + static_cast(model_->size()), + context_); checkAllocation((void*) JacobianMat_, "SUNDenseMatrix"); linearSolver_ = SUNLinSol_Dense(yy_, JacobianMat_, context_); @@ -278,7 +283,7 @@ namespace AnalysisManager int retval = 0; // Create and initialize quadratures - q_ = N_VNew_Serial(model_->sizeQuadrature(), context_); + q_ = N_VNew_Serial(static_cast(model_->sizeQuadrature()), context_); checkAllocation((void*) q_, "N_VNew_Serial"); // Set integrand function and allocate quadrature workspace @@ -364,13 +369,13 @@ namespace AnalysisManager int Ida::configureAdjoint() { // Allocate adjoint vector, derivatives and quadrature - yyB_ = N_VNew_Serial(model_->size(), context_); + yyB_ = N_VNew_Serial(static_cast(model_->size()), context_); checkAllocation((void*) yyB_, "N_VNew_Serial"); ypB_ = N_VClone(yyB_); checkAllocation((void*) ypB_, "N_VClone"); - qB_ = N_VNew_Serial(model_->sizeParams(), context_); + qB_ = N_VNew_Serial(static_cast(model_->sizeParams()), context_); checkAllocation((void*) qB_, "N_VNew_Serial"); return 0; @@ -382,7 +387,7 @@ namespace AnalysisManager int retval = 0; // Create adjoint workspace - retval = IDAAdjInit(solver_, steps, IDA_HERMITE); + retval = IDAAdjInit(solver_, static_cast(steps), IDA_HERMITE); checkOutput(retval, "IDAAdjInit"); return retval; @@ -420,7 +425,9 @@ namespace AnalysisManager checkOutput(retval, "IDASetMaxNumSteps"); // Set up linear solver - JacobianMatB_ = SUNDenseMatrix(model_->size(), model_->size(), context_); + JacobianMatB_ = SUNDenseMatrix(static_cast(model_->size()), + static_cast(model_->size()), + context_); checkAllocation((void*) JacobianMatB_, "SUNDenseMatrix"); linearSolverB_ = SUNLinSol_Dense(yyB_, JacobianMatB_, context_); @@ -450,7 +457,9 @@ namespace AnalysisManager int retval = 0; // Create Jacobian matrix - JacobianMatB_ = SUNDenseMatrix(model_->size(), model_->size(), context_); + JacobianMatB_ = SUNDenseMatrix(static_cast(model_->size()), + static_cast(model_->size()), + context_); checkAllocation((void*) JacobianMatB_, "SUNDenseMatrix"); // Create linear solver diff --git a/src/Solver/Optimization/DynamicConstraint.cpp b/src/Solver/Optimization/DynamicConstraint.cpp index 363d42fa..89e4d6af 100644 --- a/src/Solver/Optimization/DynamicConstraint.cpp +++ b/src/Solver/Optimization/DynamicConstraint.cpp @@ -64,8 +64,8 @@ namespace AnalysisManager // Get boundaries for the optimization parameters for (IdxT i = 0; i < model_->sizeParams(); ++i) { - x_l[i] = model_->param_lo()[i]; - x_u[i] = model_->param_up()[i]; + x_l[i] = model_->param_lo()[static_cast(i)]; + x_u[i] = model_->param_up()[static_cast(i)]; } // No boundaries for fictitious parameter x[n] @@ -97,7 +97,7 @@ namespace AnalysisManager // Initialize optimization parameters x for (IdxT i = 0; i < model_->sizeParams(); ++i) - x[i] = model_->param()[i]; + x[i] = model_->param()[static_cast(i)]; // Initialize fictitious parameter x[n-1] to zero x[model_->sizeParams()] = 0.0; @@ -142,7 +142,7 @@ namespace AnalysisManager // Update optimization parameters for (IdxT i = 0; i < model_->sizeParams(); ++i) { - model_->param()[i] = x[i]; + model_->param()[static_cast(i)] = x[i]; // std::cout << "x[" << i << "] = " << x[i] << "\n"; } @@ -191,7 +191,7 @@ namespace AnalysisManager { // Update optimization parameters for (IdxT i = 0; i < model_->sizeParams(); ++i) - model_->param()[i] = x[i]; + model_->param()[static_cast(i)] = x[i]; // evaluate the gradient of the objective function grad_{x} f(x) // This is creating and deleting adjoint system for each iteration! diff --git a/src/Solver/Optimization/DynamicObjective.cpp b/src/Solver/Optimization/DynamicObjective.cpp index f667a7e8..c27882c0 100644 --- a/src/Solver/Optimization/DynamicObjective.cpp +++ b/src/Solver/Optimization/DynamicObjective.cpp @@ -63,8 +63,8 @@ namespace AnalysisManager // Get boundaries for the optimization parameters for (IdxT i = 0; i < model_->sizeParams(); ++i) { - x_l[i] = model_->param_lo()[i]; - x_u[i] = model_->param_up()[i]; + x_l[i] = model_->param_lo()[static_cast(i)]; + x_u[i] = model_->param_up()[static_cast(i)]; } return true; @@ -88,7 +88,7 @@ namespace AnalysisManager // Initialize optimization parameters x for (IdxT i = 0; i < model_->sizeParams(); ++i) - x[i] = model_->param()[i]; + x[i] = model_->param()[static_cast(i)]; return true; } @@ -101,7 +101,7 @@ namespace AnalysisManager { // Update optimization parameters for (IdxT i = 0; i < model_->sizeParams(); ++i) - model_->param()[i] = x[i]; + model_->param()[static_cast(i)] = x[i]; // Evaluate objective function integrator_->getSavedInitialCondition(); @@ -124,7 +124,7 @@ namespace AnalysisManager assert(model_->sizeParams() == static_cast(n)); // Update optimization parameters for (IdxT i = 0; i < model_->sizeParams(); ++i) - model_->param()[i] = x[i]; + model_->param()[static_cast(i)] = x[i]; // evaluate the gradient of the objective function grad_{x} f(x) // This is creating and deleting adjoint system for each iteration! diff --git a/src/Solver/SteadyState/Kinsol.cpp b/src/Solver/SteadyState/Kinsol.cpp index 67f33982..b8651d20 100644 --- a/src/Solver/SteadyState/Kinsol.cpp +++ b/src/Solver/SteadyState/Kinsol.cpp @@ -8,8 +8,6 @@ * */ -#include "Kinsol.hpp" - #include #include @@ -18,6 +16,7 @@ #include // access to dense SUNLinearSolver #include // access to dense SUNMatrix +#include "Kinsol.hpp" #include "Model/Evaluator.hpp" namespace AnalysisManager @@ -53,7 +52,7 @@ namespace AnalysisManager int retval = 0; // Allocate solution vectors - yy_ = N_VNew_Serial(model_->size(), context_); + yy_ = N_VNew_Serial(static_cast(model_->size()), context_); checkAllocation((void*) yy_, "N_VNew_Serial"); // Allocate scaling vector @@ -93,7 +92,9 @@ namespace AnalysisManager int retval = 0; // Set up linear solver - JacobianMat_ = SUNDenseMatrix(model_->size(), model_->size(), context_); + JacobianMat_ = SUNDenseMatrix(static_cast(model_->size()), + static_cast(model_->size()), + context_); checkAllocation((void*) JacobianMat_, "SUNDenseMatrix"); linearSolver_ = SUNLinSol_Dense(yy_, JacobianMat_, context_);