From 023cf2d13b3f1809c0ca3765943a2bf72fc4338f Mon Sep 17 00:00:00 2001 From: Slaven Peles Date: Thu, 20 Feb 2025 17:26:48 -0500 Subject: [PATCH 1/8] Replace ModelEvaluator class. --- .../AdjointSensitivity/AdjointSensitivity.cpp | 6 +- .../DynamicConstrainedOpt.cpp | 6 +- examples/GenInfiniteBus/GenInfiniteBus.cpp | 6 +- .../ParameterEstimation.cpp | 6 +- .../{EvaluatorDynamics.hpp => Evaluator.hpp} | 6 +- src/Model/PhasorDynamics/BusBase.hpp | 6 +- src/Model/PhasorDynamics/Component.hpp | 6 +- src/ModelEvaluator.hpp | 157 ------------------ src/ModelEvaluatorImpl.hpp | 10 +- src/Solver/Dynamic/DynamicSolver.hpp | 8 +- src/Solver/Dynamic/Ida.cpp | 18 +- src/Solver/Dynamic/Ida.hpp | 6 +- src/Solver/Optimization/DynamicConstraint.cpp | 38 ++--- src/Solver/Optimization/DynamicObjective.cpp | 16 +- .../Optimization/OptimizationSolver.hpp | 4 +- src/Solver/SteadyState/Kinsol.cpp | 8 +- src/Solver/SteadyState/Kinsol.hpp | 6 +- src/Solver/SteadyState/SteadyStateSolver.hpp | 8 +- src/SystemModel.hpp | 72 ++++---- 19 files changed, 118 insertions(+), 275 deletions(-) rename src/Model/{EvaluatorDynamics.hpp => Evaluator.hpp} (96%) delete mode 100644 src/ModelEvaluator.hpp diff --git a/examples/AdjointSensitivity/AdjointSensitivity.cpp b/examples/AdjointSensitivity/AdjointSensitivity.cpp index 11cbe8ba..1516b255 100644 --- a/examples/AdjointSensitivity/AdjointSensitivity.cpp +++ b/examples/AdjointSensitivity/AdjointSensitivity.cpp @@ -142,9 +142,9 @@ int main() const double eps = 2e-3; // Compute gradient of the objective function numerically - std::vector dGdp(model->size_opt()); + std::vector dGdp(model->sizeParams()); - for (unsigned i=0; isize_opt(); ++i) + for (unsigned i=0; isizeParams(); ++i) { model->param()[i] += eps; idas->getSavedInitialCondition(); @@ -185,7 +185,7 @@ int main() int retval = 0; std::cout << "\n\nComparison of numerical and adjoint results:\n\n"; double* neg_dGdp = idas->getAdjointIntegral(); - for (unsigned i=0; isize_opt(); ++i) + for (unsigned i=0; isizeParams(); ++i) { std::cout << "dG/dp" << i << " (numerical) = " << dGdp[i] << "\n"; std::cout << "dG/dp" << i << " (adjoint) = " << -neg_dGdp[i] << "\n\n"; diff --git a/examples/DynamicConstrainedOpt/DynamicConstrainedOpt.cpp b/examples/DynamicConstrainedOpt/DynamicConstrainedOpt.cpp index 4a90a4b1..9f3fc632 100644 --- a/examples/DynamicConstrainedOpt/DynamicConstrainedOpt.cpp +++ b/examples/DynamicConstrainedOpt/DynamicConstrainedOpt.cpp @@ -165,8 +165,8 @@ int main() } // Store dynamic objective optimization results - double* results = new double[model->size_opt()]; - for(unsigned i=0; i size_opt(); ++i) + double* results = new double[model->sizeParams()]; + for(unsigned i=0; i sizeParams(); ++i) { results[i] = model->param()[i]; } @@ -193,7 +193,7 @@ int main() // Compare results of the two optimization methods int retval = 0; - for(unsigned i=0; i size_opt(); ++i) + for(unsigned i=0; i sizeParams(); ++i) { if(!isEqual(results[i], model->param()[i], 10*tol)) --retval; diff --git a/examples/GenInfiniteBus/GenInfiniteBus.cpp b/examples/GenInfiniteBus/GenInfiniteBus.cpp index 63cfd8c1..e60ac888 100644 --- a/examples/GenInfiniteBus/GenInfiniteBus.cpp +++ b/examples/GenInfiniteBus/GenInfiniteBus.cpp @@ -176,8 +176,8 @@ int main() new IpoptInterface::DynamicConstraint(idas); // Store dynamic objective optimization results - double* results = new double[model->size_opt()]; - for(unsigned i=0; i size_opt(); ++i) + double* results = new double[model->sizeParams()]; + for(unsigned i=0; i sizeParams(); ++i) { results[i] = model->param()[i]; } @@ -202,7 +202,7 @@ int main() // Compare results of the two optimization methods int retval = 0; - for(unsigned i=0; i size_opt(); ++i) + for(unsigned i=0; i sizeParams(); ++i) { if(!isEqual(results[i], model->param()[i], 100*tol)) --retval; diff --git a/examples/ParameterEstimation/ParameterEstimation.cpp b/examples/ParameterEstimation/ParameterEstimation.cpp index 54a24313..599e7ff7 100644 --- a/examples/ParameterEstimation/ParameterEstimation.cpp +++ b/examples/ParameterEstimation/ParameterEstimation.cpp @@ -178,8 +178,8 @@ int main() } // Store dynamic objective optimization results - double* results = new double[model->size_opt()]; - for(unsigned i=0; i size_opt(); ++i) + double* results = new double[model->sizeParams()]; + for(unsigned i=0; i sizeParams(); ++i) { results[i] = model->param()[i]; } @@ -207,7 +207,7 @@ int main() // Compare results of the two optimization methods int retval = 0; - for(unsigned i=0; i size_opt(); ++i) + for(unsigned i=0; i sizeParams(); ++i) { if(!isEqual(results[i], model->param()[i], 100*tol)) --retval; diff --git a/src/Model/EvaluatorDynamics.hpp b/src/Model/Evaluator.hpp similarity index 96% rename from src/Model/EvaluatorDynamics.hpp rename to src/Model/Evaluator.hpp index ffa7ad7c..669834dc 100644 --- a/src/Model/EvaluatorDynamics.hpp +++ b/src/Model/Evaluator.hpp @@ -13,13 +13,13 @@ namespace Model * */ template - class EvaluatorDynamics + class Evaluator { public: typedef typename GridKit::ScalarTraits::real_type real_type; - EvaluatorDynamics(){} - virtual ~EvaluatorDynamics(){} + Evaluator(){} + virtual ~Evaluator(){} virtual int allocate() = 0; virtual int initialize() = 0; diff --git a/src/Model/PhasorDynamics/BusBase.hpp b/src/Model/PhasorDynamics/BusBase.hpp index 6a94c647..30f1b687 100644 --- a/src/Model/PhasorDynamics/BusBase.hpp +++ b/src/Model/PhasorDynamics/BusBase.hpp @@ -1,7 +1,7 @@ #pragma once #include -#include +#include namespace GridKit { @@ -12,10 +12,10 @@ namespace PhasorDynamics * */ template - class BusBase : public Model::EvaluatorDynamics + class BusBase : public Model::Evaluator { public: - using real_type = typename Model::EvaluatorDynamics::real_type; + using real_type = typename Model::Evaluator::real_type; enum BusType{DEFAULT=1, SLACK}; diff --git a/src/Model/PhasorDynamics/Component.hpp b/src/Model/PhasorDynamics/Component.hpp index deb83cd8..3a724a02 100644 --- a/src/Model/PhasorDynamics/Component.hpp +++ b/src/Model/PhasorDynamics/Component.hpp @@ -1,7 +1,7 @@ #pragma once #include -#include +#include namespace GridKit { @@ -13,10 +13,10 @@ namespace PhasorDynamics * */ template - class Component : public Model::EvaluatorDynamics + class Component : public Model::Evaluator { public: - using real_type = typename Model::EvaluatorDynamics::real_type; + using real_type = typename Model::Evaluator::real_type; Component() : size_(0), diff --git a/src/ModelEvaluator.hpp b/src/ModelEvaluator.hpp deleted file mode 100644 index dde5130b..00000000 --- a/src/ModelEvaluator.hpp +++ /dev/null @@ -1,157 +0,0 @@ -/* - * - * Copyright (c) 2017, Lawrence Livermore National Security, LLC. - * Produced at the Lawrence Livermore National Laboratory. - * Written by Slaven Peles . - * LLNL-CODE-718378. - * All rights reserved. - * - * This file is part of GridKit™. For details, see github.com/LLNL/GridKit - * Please also read the LICENSE file. - * - * Redistribution and use in source and binary forms, with or without - * modification, are permitted provided that the following conditions are met: - * - Redistributions of source code must retain the above copyright notice, - * this list of conditions and the disclaimer below. - * - Redistributions in binary form must reproduce the above copyright notice, - * this list of conditions and the disclaimer (as noted below) in the - * documentation and/or other materials provided with the distribution. - * - Neither the name of the LLNS/LLNL nor the names of its contributors may - * be used to endorse or promote products derived from this software without - * specific prior written permission. - * - * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND - * CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, - * INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF - * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE - * DISCLAIMED. IN NO EVENT SHALL LAWRENCE LIVERMORE NATIONAL - * SECURITY, LLC, THE U.S. DEPARTMENT OF ENERGY OR CONTRIBUTORS BE - * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, - * OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, - * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, - * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND - * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT - * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISINGIN ANY - * WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF - * THE POSSIBILITY OF SUCH DAMAGE. - * - * Lawrence Livermore National Laboratory is operated by Lawrence Livermore - * National Security, LLC, for the U.S. Department of Energy, National - * Nuclear Security Administration under Contract DE-AC52-07NA27344. - * - * This document was prepared as an account of work sponsored by an agency - * of the United States government. Neither the United States government nor - * Lawrence Livermore National Security, LLC, nor any of their employees - * makes any warranty, expressed or implied, or assumes any legal liability - * or responsibility for the accuracy, completeness, or usefulness of any - * information, apparatus, product, or process disclosed, or represents that - * its use would not infringe privately owned rights. Reference herein to - * any specific commercial product, process, or service by trade name, - * trademark, manufacturer, or otherwise does not necessarily constitute or - * imply its endorsement, recommendation, or favoring by the United States - * government or Lawrence Livermore National Security, LLC. The views and - * opinions of authors expressed herein do not necessarily state or reflect - * those of the United States government or Lawrence Livermore National - * Security, LLC, and shall not be used for advertising or product - * endorsement purposes. - * - */ - -#ifndef _MODEL_EVALUATOR_HPP_ -#define _MODEL_EVALUATOR_HPP_ - -#include -#include -#include - -namespace GridKit -{ - - /*! - * @brief Abstract class describing a model. - * - */ - template - class ModelEvaluator - { - public: - typedef typename GridKit::ScalarTraits::real_type real_type; - - ModelEvaluator(){} - virtual ~ModelEvaluator(){} - - virtual int allocate() = 0; - virtual int initialize() = 0; - virtual int tagDifferentiable() = 0; - virtual int evaluateResidual() = 0; - virtual int evaluateJacobian() = 0; - virtual int evaluateIntegrand() = 0; - - virtual int initializeAdjoint() = 0; - virtual int evaluateAdjointResidual() = 0; - // virtual int evaluateAdjointJacobian() = 0; - virtual int evaluateAdjointIntegrand() = 0; - - virtual IdxT size() = 0; - virtual IdxT nnz() = 0; - - /** - * @brief Is the Jacobian defined. Used in IDA to determine wether DQ is used or not - * - * @return true - * @return false - */ - virtual bool hasJacobian() = 0; - - virtual IdxT size_quad() = 0; - virtual IdxT size_opt() = 0; - virtual void updateTime(real_type t, real_type a) = 0; - virtual void setTolerances(real_type& rel_tol, real_type& abs_tol) const = 0; - virtual void setMaxSteps(IdxT& msa) const = 0; - - virtual std::vector& y() = 0; - virtual const std::vector& y() const = 0; - - virtual std::vector& yp() = 0; - virtual const std::vector& yp() const = 0; - - virtual std::vector& tag() = 0; - virtual const std::vector& tag() const = 0; - - virtual std::vector& yB() = 0; - virtual const std::vector& yB() const = 0; - - virtual std::vector& ypB() = 0; - virtual const std::vector& ypB() const = 0; - - virtual std::vector& param() = 0; - virtual const std::vector& param() const = 0; - - virtual std::vector& param_up() = 0; - virtual const std::vector& param_up() const = 0; - - virtual std::vector& param_lo() = 0; - virtual const std::vector& param_lo() const = 0; - - virtual std::vector& getResidual() = 0; - virtual const std::vector& getResidual() const = 0; - - - virtual COO_Matrix& getJacobian() = 0; - virtual const COO_Matrix& getJacobian() const = 0; - - virtual std::vector& getIntegrand() = 0; - virtual const std::vector& getIntegrand() const = 0; - - virtual std::vector& getAdjointResidual() = 0; - virtual const std::vector& getAdjointResidual() const = 0; - - virtual std::vector& getAdjointIntegrand() = 0; - virtual const std::vector& getAdjointIntegrand() const = 0; - - }; - - -} // namespace GridKit - -#endif // _MODEL_EVALUATOR_HPP_ \ No newline at end of file diff --git a/src/ModelEvaluatorImpl.hpp b/src/ModelEvaluatorImpl.hpp index e3b15863..d08e247d 100644 --- a/src/ModelEvaluatorImpl.hpp +++ b/src/ModelEvaluatorImpl.hpp @@ -61,7 +61,7 @@ #define _MODEL_EVALUATOR_IMPL_HPP_ #include -#include +#include namespace GridKit { @@ -71,10 +71,10 @@ namespace GridKit * */ template - class ModelEvaluatorImpl : public ModelEvaluator + class ModelEvaluatorImpl : public Model::Evaluator { public: - typedef typename ModelEvaluator::real_type real_type; + typedef typename Model::Evaluator::real_type real_type; ModelEvaluatorImpl() : size_(0), @@ -116,12 +116,12 @@ namespace GridKit return false; } - virtual IdxT size_quad() + virtual IdxT sizeQuadrature() { return size_quad_; } - virtual IdxT size_opt() + virtual IdxT sizeParams() { return size_opt_; } diff --git a/src/Solver/Dynamic/DynamicSolver.hpp b/src/Solver/Dynamic/DynamicSolver.hpp index 3f1fb546..ab7387f1 100644 --- a/src/Solver/Dynamic/DynamicSolver.hpp +++ b/src/Solver/Dynamic/DynamicSolver.hpp @@ -61,7 +61,7 @@ #ifndef _DYNAMIC_SOLVER_HPP_ #define _DYNAMIC_SOLVER_HPP_ -#include "ModelEvaluator.hpp" +#include "Model/Evaluator.hpp" namespace AnalysisManager { @@ -69,18 +69,18 @@ namespace AnalysisManager class DynamicSolver { public: - DynamicSolver(GridKit::ModelEvaluator* model) : model_(model) + DynamicSolver(GridKit::Model::Evaluator* model) : model_(model) { } virtual ~DynamicSolver(){} - GridKit::ModelEvaluator* getModel() + GridKit::Model::Evaluator* getModel() { return model_; } protected: - GridKit::ModelEvaluator* model_; + GridKit::Model::Evaluator* model_; }; } diff --git a/src/Solver/Dynamic/Ida.cpp b/src/Solver/Dynamic/Ida.cpp index 0c81ada1..f0e6478b 100644 --- a/src/Solver/Dynamic/Ida.cpp +++ b/src/Solver/Dynamic/Ida.cpp @@ -64,7 +64,7 @@ #include #include -#include "ModelEvaluator.hpp" +#include "Model/Evaluator.hpp" #include "Ida.hpp" @@ -75,7 +75,7 @@ namespace Sundials { template - Ida::Ida(GridKit::ModelEvaluator* model) : DynamicSolver(model) + Ida::Ida(GridKit::Model::Evaluator* model) : DynamicSolver(model) { int retval = 0; @@ -316,7 +316,7 @@ namespace Sundials int retval = 0; // Create and initialize quadratures - q_ = N_VNew_Serial(model_->size_quad(), context_); + q_ = N_VNew_Serial(model_->sizeQuadrature(), context_); checkAllocation((void*) q_, "N_VNew_Serial"); // Set integrand function and allocate quadrature workspace @@ -412,7 +412,7 @@ namespace Sundials ypB_ = N_VClone(yyB_); checkAllocation((void*) ypB_, "N_VClone"); - qB_ = N_VNew_Serial(model_->size_opt(), context_); + qB_ = N_VNew_Serial(model_->sizeParams(), context_); checkAllocation((void*) qB_, "N_VNew_Serial"); return 0; @@ -579,7 +579,7 @@ namespace Sundials template int Ida::Residual(sunrealtype tres, N_Vector yy, N_Vector yp, N_Vector rr, void *user_data) { - GridKit::ModelEvaluator* model = static_cast*>(user_data); + GridKit::Model::Evaluator* model = static_cast*>(user_data); model->updateTime(tres, 0.0); copyVec(yy, model->y()); @@ -596,7 +596,7 @@ namespace Sundials int Ida::Jac(sunrealtype t, sunrealtype cj, N_Vector yy, N_Vector yp, N_Vector resvec, SUNMatrix J, void *user_data, N_Vector tmp1, N_Vector tmp2, N_Vector tmp3) { - GridKit::ModelEvaluator* model = static_cast*>(user_data); + GridKit::Model::Evaluator* model = static_cast*>(user_data); model->updateTime(t, cj); @@ -637,7 +637,7 @@ namespace Sundials template int Ida::Integrand(sunrealtype tt, N_Vector yy, N_Vector yp, N_Vector rhsQ, void *user_data) { - GridKit::ModelEvaluator* model = static_cast*>(user_data); + GridKit::Model::Evaluator* model = static_cast*>(user_data); model->updateTime(tt, 0.0); copyVec(yy, model->y()); @@ -653,7 +653,7 @@ namespace Sundials template int Ida::adjointResidual(sunrealtype tt, N_Vector yy, N_Vector yp, N_Vector yyB, N_Vector ypB, N_Vector rrB, void *user_data) { - GridKit::ModelEvaluator* model = static_cast*>(user_data); + GridKit::Model::Evaluator* model = static_cast*>(user_data); model->updateTime(tt, 0.0); copyVec(yy, model->y()); @@ -672,7 +672,7 @@ namespace Sundials template int Ida::adjointIntegrand(sunrealtype tt, N_Vector yy, N_Vector yp, N_Vector yyB, N_Vector ypB, N_Vector rhsQB, void *user_data) { - GridKit::ModelEvaluator* model = static_cast*>(user_data); + GridKit::Model::Evaluator* model = static_cast*>(user_data); model->updateTime(tt, 0.0); copyVec(yy, model->y()); diff --git a/src/Solver/Dynamic/Ida.hpp b/src/Solver/Dynamic/Ida.hpp index 54011344..272db453 100644 --- a/src/Solver/Dynamic/Ida.hpp +++ b/src/Solver/Dynamic/Ida.hpp @@ -68,7 +68,7 @@ #include /* access to KLU linear solver */ #include /* access to dense linear solver */ -#include "ModelEvaluator.hpp" +#include "Model/Evaluator.hpp" #include "DynamicSolver.hpp" namespace AnalysisManager @@ -83,7 +83,7 @@ namespace AnalysisManager typedef typename GridKit::ScalarTraits::real_type real_type; public: - Ida(GridKit::ModelEvaluator* model); + Ida(GridKit::Model::Evaluator* model); ~Ida(); int configureSimulation(); @@ -212,7 +212,7 @@ namespace AnalysisManager int backwardID_; private: - //static void copyMat(ModelEvaluator::Mat& J, SlsMat Jida); + //static void copyMat(Model::Evaluator::Mat& J, SlsMat Jida); static void copyVec(const N_Vector x, std::vector& y); static void copyVec(const std::vector& x, N_Vector y); static void copyVec(const std::vector& x, N_Vector y); diff --git a/src/Solver/Optimization/DynamicConstraint.cpp b/src/Solver/Optimization/DynamicConstraint.cpp index 05ba82bc..bdec5001 100644 --- a/src/Solver/Optimization/DynamicConstraint.cpp +++ b/src/Solver/Optimization/DynamicConstraint.cpp @@ -91,13 +91,13 @@ bool DynamicConstraint::get_nlp_info(Index& n, Index& m, Index& n // Number of parameters is size of the system plus 1 fictitious parameter // to store the objective value. - n = model_->size_opt() + 1; + n = model_->sizeParams() + 1; // There is one constraint m = 1; // Jacobian is a dense row matrix of length n+1. - nnz_jac_g = model_->size_opt() + 1; + nnz_jac_g = model_->sizeParams() + 1; // Using numerical Hessian. nnz_h_lag = 0; @@ -114,19 +114,19 @@ bool DynamicConstraint::get_bounds_info(Index n, Number* x_l, Num Index m, Number* g_l, Number* g_u) { // Check if sizes are set correctly - assert(n == (Index) (model_->size_opt() + 1)); + assert(n == (Index) (model_->sizeParams() + 1)); assert(m == 1); // Get boundaries for the optimization parameters - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) { x_l[i] = model_->param_lo()[i]; x_u[i] = model_->param_up()[i]; } // No boundaries for fictitious parameter x[n] - x_l[model_->size_opt()] = -1e20; - x_u[model_->size_opt()] = +1e20; + x_l[model_->sizeParams()] = -1e20; + x_u[model_->sizeParams()] = +1e20; // Set constraint g[0] to be equality constraint g[0] = 0 g_l[0] = 0.0; @@ -148,11 +148,11 @@ bool DynamicConstraint::get_starting_point(Index n, bool init_x, assert(init_lambda == false); // Initialize optimization parameters x - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) x[i] = model_->param()[i]; // Initialize fictitious parameter x[n-1] to zero - x[model_->size_opt()] = 0.0; + x[model_->sizeParams()] = 0.0; return true; } @@ -162,7 +162,7 @@ template bool DynamicConstraint::eval_f(Index n, const Number* x, bool new_x, Number& obj_value) { // Set objective to fictitious optimization parameter x[n-1] - obj_value = x[model_->size_opt()]; + obj_value = x[model_->sizeParams()]; return true; } @@ -173,9 +173,9 @@ bool DynamicConstraint::eval_grad_f(Index n, const Number* x, boo { // Objective function equals to the fictitious parameter x[n-1]. // Gradient, then assumes the simple form: - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) grad_f[i] = 0.0; - grad_f[model_->size_opt()] = 1.0; + grad_f[model_->sizeParams()] = 1.0; return true; } @@ -185,7 +185,7 @@ template bool DynamicConstraint::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g) { // Update optimization parameters - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) { model_->param()[i] = x[i]; //std::cout << "x[" << i << "] = " << x[i] << "\n"; @@ -205,7 +205,7 @@ bool DynamicConstraint::eval_g(Index n, const Number* x, bool new } // For now assumes only one forward integrand and multiple optimization parameters. - g[0] = (integrator_->getIntegral())[0] - x[model_->size_opt()]; + g[0] = (integrator_->getIntegral())[0] - x[model_->sizeParams()]; //std::cout << "constraint:" << g[0] << std::endl; return true; } @@ -219,19 +219,19 @@ bool DynamicConstraint::eval_jac_g(Index n, const Number* x, bool // Set Jacobian sparsity pattern ... if(!values) { - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) { iRow[i] = 0; jCol[i] = i; } - iRow[model_->size_opt()] = 0; - jCol[model_->size_opt()] = model_->size_opt(); + iRow[model_->sizeParams()] = 0; + jCol[model_->sizeParams()] = model_->sizeParams(); } // ... or compute Jacobian derivatives else { // Update optimization parameters - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) model_->param()[i] = x[i]; // evaluate the gradient of the objective function grad_{x} f(x) @@ -261,11 +261,11 @@ bool DynamicConstraint::eval_jac_g(Index n, const Number* x, bool } // For now assumes only one forward integrand and multiple optimization parameters. - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) { values[i] = -((integrator_->getAdjointIntegral())[i]); } - values[model_->size_opt()] = -1.0; + values[model_->sizeParams()] = -1.0; integrator_->deleteAdjoint(); } diff --git a/src/Solver/Optimization/DynamicObjective.cpp b/src/Solver/Optimization/DynamicObjective.cpp index adc743ba..6898354c 100644 --- a/src/Solver/Optimization/DynamicObjective.cpp +++ b/src/Solver/Optimization/DynamicObjective.cpp @@ -90,7 +90,7 @@ bool DynamicObjective::get_nlp_info(Index& n, Index& m, Index& nn assert(model_->size_quad() == 1); // Number of optimization variables. - n = model_->size_opt(); + n = model_->sizeParams(); // There are no constraints m = 0; @@ -113,11 +113,11 @@ bool DynamicObjective::get_bounds_info(Index n, Number* x_l, Numb Index m, Number* g_l, Number* g_u) { // Check if sizes are set correctly - assert(n == (Index) model_->size_opt()); + assert(n == (Index) model_->sizeParams()); assert(m == 0); // Get boundaries for the optimization parameters - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) { x_l[i] = model_->param_lo()[i]; x_u[i] = model_->param_up()[i]; @@ -139,7 +139,7 @@ bool DynamicObjective::get_starting_point(Index n, bool init_x, N assert(init_lambda == false); // Initialize optimization parameters x - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) x[i] = model_->param()[i]; return true; @@ -150,7 +150,7 @@ template bool DynamicObjective::eval_f(Index n, const Number* x, bool new_x, Number& obj_value) { // Update optimization parameters - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) model_->param()[i] = x[i]; // Evaluate objective function @@ -169,9 +169,9 @@ bool DynamicObjective::eval_f(Index n, const Number* x, bool new_ template bool DynamicObjective::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f) { - assert(model_->size_opt() == static_cast(n)); + assert(model_->sizeParams() == static_cast(n)); // Update optimization parameters - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) model_->param()[i] = x[i]; // evaluate the gradient of the objective function grad_{x} f(x) @@ -188,7 +188,7 @@ bool DynamicObjective::eval_grad_f(Index n, const Number* x, bool integrator_->runBackwardSimulation(t_init_); // For now assumes only one forward integrand and multiple optimization parameters. - for(IdxT i = 0; i < model_->size_opt(); ++i) + for(IdxT i = 0; i < model_->sizeParams(); ++i) grad_f[i] = -((integrator_->getAdjointIntegral())[i]); integrator_->deleteAdjoint(); diff --git a/src/Solver/Optimization/OptimizationSolver.hpp b/src/Solver/Optimization/OptimizationSolver.hpp index ad7a286d..d89657d4 100644 --- a/src/Solver/Optimization/OptimizationSolver.hpp +++ b/src/Solver/Optimization/OptimizationSolver.hpp @@ -60,7 +60,7 @@ #ifndef _OPTIMIZATION_SOLVER_HPP_ #define _OPTIMIZATION_SOLVER_HPP_ -#include "ModelEvaluator.hpp" +#include "Model/Evaluator.hpp" #include namespace AnalysisManager @@ -83,7 +83,7 @@ namespace AnalysisManager virtual ~OptimizationSolver(){} protected: - GridKit::ModelEvaluator* model_; + GridKit::Model::Evaluator* model_; Sundials::Ida* integrator_; }; diff --git a/src/Solver/SteadyState/Kinsol.cpp b/src/Solver/SteadyState/Kinsol.cpp index 2a764b31..0b14226f 100644 --- a/src/Solver/SteadyState/Kinsol.cpp +++ b/src/Solver/SteadyState/Kinsol.cpp @@ -74,7 +74,7 @@ #include // access to dense SUNMatrix #include // access to dense SUNLinearSolver -#include "ModelEvaluator.hpp" +#include "Model/Evaluator.hpp" #include "Kinsol.hpp" @@ -85,7 +85,7 @@ namespace Sundials { template - Kinsol::Kinsol(GridKit::ModelEvaluator* model) + Kinsol::Kinsol(GridKit::Model::Evaluator* model) : SteadyStateSolver(model) { int retval = 0; @@ -216,8 +216,8 @@ namespace Sundials template int Kinsol::Residual(N_Vector yy, N_Vector rr, void *user_data) { - GridKit::ModelEvaluator* model = - static_cast*>(user_data); + GridKit::Model::Evaluator* model = + static_cast*>(user_data); copyVec(yy, model->y()); diff --git a/src/Solver/SteadyState/Kinsol.hpp b/src/Solver/SteadyState/Kinsol.hpp index eed77cdb..0bd972df 100644 --- a/src/Solver/SteadyState/Kinsol.hpp +++ b/src/Solver/SteadyState/Kinsol.hpp @@ -75,7 +75,7 @@ // #include /* access to KLU linear solver */ #include /* access to dense linear solver */ -#include "ModelEvaluator.hpp" +#include "Model/Evaluator.hpp" #include "SteadyStateSolver.hpp" namespace AnalysisManager @@ -90,7 +90,7 @@ namespace AnalysisManager typedef typename GridKit::ScalarTraits::real_type real_type; public: - Kinsol(GridKit::ModelEvaluator* model); + Kinsol(GridKit::Model::Evaluator* model); ~Kinsol(); int configureSimulation(); @@ -197,7 +197,7 @@ namespace AnalysisManager N_Vector yy0_; ///< Storage for initial values private: - //static void copyMat(ModelEvaluator::Mat& J, SlsMat Jida); + //static void copyMat(Model::Evaluator::Mat& J, SlsMat Jida); static void copyVec(const N_Vector x, std::vector& y); static void copyVec(const std::vector& x, N_Vector y); static void copyVec(const std::vector& x, N_Vector y); diff --git a/src/Solver/SteadyState/SteadyStateSolver.hpp b/src/Solver/SteadyState/SteadyStateSolver.hpp index 78a17255..a17024b4 100644 --- a/src/Solver/SteadyState/SteadyStateSolver.hpp +++ b/src/Solver/SteadyState/SteadyStateSolver.hpp @@ -58,7 +58,7 @@ */ #pragma once -#include "ModelEvaluator.hpp" +#include "Model/Evaluator.hpp" namespace AnalysisManager { @@ -66,18 +66,18 @@ namespace AnalysisManager class SteadyStateSolver { public: - SteadyStateSolver(GridKit::ModelEvaluator* model) : model_(model) + SteadyStateSolver(GridKit::Model::Evaluator* model) : model_(model) { } virtual ~SteadyStateSolver(){} - GridKit::ModelEvaluator* getModel() + GridKit::Model::Evaluator* getModel() { return model_; } protected: - GridKit::ModelEvaluator* model_; + GridKit::Model::Evaluator* model_; }; } diff --git a/src/SystemModel.hpp b/src/SystemModel.hpp index 1b744330..b9139d25 100644 --- a/src/SystemModel.hpp +++ b/src/SystemModel.hpp @@ -74,7 +74,7 @@ namespace GridKit * @brief Prototype for a system model class * * This class maps component data to system data and implements - * ModelEvaluator for the system model. This is still work in + * Model::Evaluator for the system model. This is still work in * progress and code is not optimized. * * @todo Address thread safety for the system model methods. @@ -83,8 +83,8 @@ namespace GridKit template class SystemModel : public ModelEvaluatorImpl { - typedef ModelEvaluator bus_type; - typedef ModelEvaluator component_type; + using bus_type = Model::Evaluator; + using component_type = Model::Evaluator; using real_type = typename ModelEvaluatorImpl::real_type; using ModelEvaluatorImpl::size_; @@ -150,8 +150,8 @@ class SystemModel : public ModelEvaluatorImpl { bus->allocate(); size_ += bus->size(); - size_quad_ += bus->size_quad(); - size_opt_ += bus->size_opt(); + size_quad_ += bus->sizeQuadrature(); + size_opt_ += bus->sizeParams(); } // Allocate all components @@ -159,8 +159,8 @@ class SystemModel : public ModelEvaluatorImpl { component->allocate(); size_ += component->size(); - size_quad_ += component->size_quad(); - size_opt_ += component->size_opt(); + size_quad_ += component->sizeQuadrature(); + size_opt_ += component->sizeParams(); } // Allocate global vectors @@ -230,13 +230,13 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += bus->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { param_[optOffset + j] = bus->param()[j]; param_lo_[optOffset + j] = bus->param_lo()[j]; param_up_[optOffset + j] = bus->param_up()[j]; } - optOffset += bus->size_opt(); + optOffset += bus->sizeParams(); } // Initialize components @@ -254,13 +254,13 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += component->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { param_[optOffset + j] = component->param()[j]; param_lo_[optOffset + j] = component->param_lo()[j]; param_up_[optOffset + j] = component->param_up()[j]; } - optOffset += component->size_opt(); + optOffset += component->sizeParams(); } return 0; @@ -332,11 +332,11 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += bus->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { bus->param()[j] = param_[optOffset + j]; } - optOffset += bus->size_opt(); + optOffset += bus->sizeParams(); bus->evaluateResidual(); } @@ -350,11 +350,11 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += component->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { component->param()[j] = param_[optOffset + j]; } - optOffset += component->size_opt(); + optOffset += component->sizeParams(); component->evaluateResidual(); } @@ -409,11 +409,11 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += bus->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { bus->param()[j] = param_[optOffset + j]; } - optOffset += bus->size_opt(); + optOffset += bus->sizeParams(); bus->evaluateIntegrand(); } @@ -427,11 +427,11 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += component->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { component->param()[j] = param_[optOffset + j]; } - optOffset += component->size_opt(); + optOffset += component->sizeParams(); component->evaluateIntegrand(); } @@ -440,20 +440,20 @@ class SystemModel : public ModelEvaluatorImpl IdxT intOffset = 0; for(const auto& bus: buses_) { - for(IdxT j=0; jsize_quad(); ++j) + for(IdxT j=0; jsizeQuadrature(); ++j) { g_[intOffset + j] = bus->getIntegrand()[j]; } - intOffset += bus->size_quad(); + intOffset += bus->sizeQuadrature(); } for(const auto& component : components_) { - for(IdxT j=0; jsize_quad(); ++j) + for(IdxT j=0; jsizeQuadrature(); ++j) { g_[intOffset + j] = component->getIntegrand()[j]; } - intOffset += component->size_quad(); + intOffset += component->sizeQuadrature(); } return 0; @@ -480,11 +480,11 @@ class SystemModel : public ModelEvaluatorImpl } offset += bus->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { bus->param()[j] = param_[optOffset + j]; } - optOffset += bus->size_opt(); + optOffset += bus->sizeParams(); } // Update component variables and optimization parameters @@ -497,11 +497,11 @@ class SystemModel : public ModelEvaluatorImpl } offset += component->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { component->param()[j] = param_[optOffset + j]; } - optOffset += component->size_opt(); + optOffset += component->sizeParams(); } // Reset counter @@ -560,11 +560,11 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += bus->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { bus->param()[j] = param_[optOffset + j]; } - optOffset += bus->size_opt(); + optOffset += bus->sizeParams(); } @@ -579,11 +579,11 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += component->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { component->param()[j] = param_[optOffset + j]; } - optOffset += component->size_opt(); + optOffset += component->sizeParams(); } @@ -645,11 +645,11 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += bus->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { bus->param()[j] = param_[optOffset + j]; } - optOffset += bus->size_opt(); + optOffset += bus->sizeParams(); } for(const auto& component : components_) @@ -663,17 +663,17 @@ class SystemModel : public ModelEvaluatorImpl } varOffset += component->size(); - for(IdxT j=0; jsize_opt(); ++j) + for(IdxT j=0; jsizeParams(); ++j) { component->param()[j] = param_[optOffset + j]; } - optOffset += component->size_opt(); + optOffset += component->sizeParams(); } // Evaluate integrand and update global vector for(const auto& component : components_) { - if(component->size_quad() == 1) + if(component->sizeQuadrature() == 1) { component->evaluateAdjointIntegrand(); for(IdxT j=0; j Date: Thu, 20 Feb 2025 17:47:42 -0500 Subject: [PATCH 2/8] Remove ModelEvaluatorImpl dependence for power electronics models. --- .../PowerElectronics/CircuitComponent.hpp | 213 +++++++++++++++++- 1 file changed, 206 insertions(+), 7 deletions(-) diff --git a/src/Model/PowerElectronics/CircuitComponent.hpp b/src/Model/PowerElectronics/CircuitComponent.hpp index fd88925c..1513c86d 100644 --- a/src/Model/PowerElectronics/CircuitComponent.hpp +++ b/src/Model/PowerElectronics/CircuitComponent.hpp @@ -1,9 +1,8 @@ -#ifndef _CIRCCOMP_HPP_ -#define _CIRCCOMP_HPP_ +#pragma once -#include +#include #include #include #include @@ -16,11 +15,13 @@ namespace GridKit * */ template - class CircuitComponent : public ModelEvaluatorImpl + class CircuitComponent : public Model::Evaluator { - public: + typedef typename Model::Evaluator::real_type real_type; + CircuitComponent() = default; + ~CircuitComponent() = default; void updateTime(ScalarT t, ScalarT a) { @@ -71,6 +72,174 @@ namespace GridKit return connection_nodes_.at(local_index); } + public: + virtual IdxT size() + { + return size_; + } + + virtual IdxT nnz() + { + return nnz_; + } + + virtual IdxT sizeQuadrature() + { + return size_quad_; + } + + virtual IdxT sizeParams() + { + return size_opt_; + } + + virtual void setTolerances(real_type& rel_tol, real_type& abs_tol) const + { + rel_tol = rel_tol_; + abs_tol = abs_tol_; + } + + virtual void setMaxSteps(IdxT& msa) const + { + msa = max_steps_; + } + + std::vector& y() + { + return y_; + } + + const std::vector& y() const + { + return y_; + } + + std::vector& yp() + { + return yp_; + } + + const std::vector& yp() const + { + return yp_; + } + + std::vector& tag() + { + return tag_; + } + + const std::vector& tag() const + { + return tag_; + } + + std::vector& yB() + { + return yB_; + } + + const std::vector& yB() const + { + return yB_; + } + + std::vector& ypB() + { + return ypB_; + } + + const std::vector& ypB() const + { + return ypB_; + } + + std::vector& param() + { + return param_; + } + + const std::vector& param() const + { + return param_; + } + + std::vector& param_up() + { + return param_up_; + } + + const std::vector& param_up() const + { + return param_up_; + } + + std::vector& param_lo() + { + return param_lo_; + } + + const std::vector& param_lo() const + { + return param_lo_; + } + + std::vector& getResidual() + { + return f_; + } + + const std::vector& getResidual() const + { + return f_; + } + + COO_Matrix& getJacobian() + { + return jac_; + } + + const COO_Matrix& getJacobian() const + { + return jac_; + } + + std::vector& getIntegrand() + { + return g_; + } + + const std::vector& getIntegrand() const + { + return g_; + } + + std::vector& getAdjointResidual() + { + return fB_; + } + + const std::vector& getAdjointResidual() const + { + return fB_; + } + + std::vector& getAdjointIntegrand() + { + return gB_; + } + + const std::vector& getAdjointIntegrand() const + { + return gB_; + } + + //@todo Fix ID naming + IdxT getIDcomponent() + { + return idc_; + } + protected: size_t n_extern_; size_t n_intern_; @@ -78,9 +247,39 @@ namespace GridKit ///@todo may want to replace the mapping of connection_nodes to Node objects instead of IdxT. Allows for container free setup std::map connection_nodes_; + protected: + IdxT size_{0}; + IdxT nnz_{0}; + IdxT size_quad_{0}; + IdxT size_opt_{0}; + + std::vector y_; + std::vector yp_; + std::vector tag_; + std::vector f_; + std::vector g_; + + std::vector yB_; + std::vector ypB_; + std::vector fB_; + std::vector gB_; + + COO_Matrix jac_; + + std::vector param_; + std::vector param_up_; + std::vector param_lo_; + + real_type time_; + real_type alpha_; + + real_type rel_tol_; + real_type abs_tol_; + + IdxT max_steps_; + + IdxT idc_; }; } - -#endif From 94d522c86f6a125962e219165f441aba876681d0 Mon Sep 17 00:00:00 2001 From: Slaven Peles Date: Thu, 20 Feb 2025 18:02:58 -0500 Subject: [PATCH 3/8] Remove ModelEvaluatorImpl includes from power electronics components. --- .../PowerElectronics/Capacitor/Capacitor.hpp | 2 +- .../DistributedGenerator.hpp | 2 +- .../InductionMotor/InductionMotor.hpp | 2 +- .../PowerElectronics/Inductor/Inductor.hpp | 2 +- .../LinearTransformer/LinearTransformer.hpp | 2 +- .../MicrogridBusDQ/MicrogridBusDQ.hpp | 2 +- .../MicrogridLine/MicrogridLine.hpp | 2 +- .../MicrogridLoad/MicrogridLoad.hpp | 2 +- .../PowerElectronics/Resistor/Resistor.hpp | 2 +- .../SynchronousMachine/SynchronousMachine.hpp | 2 +- .../SystemModelPowerElectronics.hpp | 44 +++++++------------ .../TransmissionLine/TransmissionLine.hpp | 2 +- .../VoltageSource/VoltageSource.hpp | 2 +- 13 files changed, 27 insertions(+), 41 deletions(-) diff --git a/src/Model/PowerElectronics/Capacitor/Capacitor.hpp b/src/Model/PowerElectronics/Capacitor/Capacitor.hpp index ad2f9bdb..b40e2f2a 100644 --- a/src/Model/PowerElectronics/Capacitor/Capacitor.hpp +++ b/src/Model/PowerElectronics/Capacitor/Capacitor.hpp @@ -3,7 +3,7 @@ #ifndef _CAP_HPP_ #define _CAP_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.hpp b/src/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.hpp index b5cd28ae..d9f7b005 100644 --- a/src/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.hpp +++ b/src/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.hpp @@ -3,7 +3,7 @@ #ifndef _CAP_HPP_ #define _CAP_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/InductionMotor/InductionMotor.hpp b/src/Model/PowerElectronics/InductionMotor/InductionMotor.hpp index 97ef633e..3223352c 100644 --- a/src/Model/PowerElectronics/InductionMotor/InductionMotor.hpp +++ b/src/Model/PowerElectronics/InductionMotor/InductionMotor.hpp @@ -3,7 +3,7 @@ #ifndef _IMOTOR_HPP_ #define _IMOTOR_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/Inductor/Inductor.hpp b/src/Model/PowerElectronics/Inductor/Inductor.hpp index 8f2a0485..fc459e90 100644 --- a/src/Model/PowerElectronics/Inductor/Inductor.hpp +++ b/src/Model/PowerElectronics/Inductor/Inductor.hpp @@ -3,7 +3,7 @@ #ifndef _IND_HPP_ #define _IND_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/LinearTransformer/LinearTransformer.hpp b/src/Model/PowerElectronics/LinearTransformer/LinearTransformer.hpp index b09373d0..3d0324ae 100644 --- a/src/Model/PowerElectronics/LinearTransformer/LinearTransformer.hpp +++ b/src/Model/PowerElectronics/LinearTransformer/LinearTransformer.hpp @@ -3,7 +3,7 @@ #ifndef _LTRANS_HPP_ #define _LTRANS_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.hpp b/src/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.hpp index 880cc494..dcbca2e9 100644 --- a/src/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.hpp +++ b/src/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.hpp @@ -3,7 +3,7 @@ #ifndef _VIRBUSDQ_HPP_ #define _VIRBUSDQ_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/MicrogridLine/MicrogridLine.hpp b/src/Model/PowerElectronics/MicrogridLine/MicrogridLine.hpp index 01affa64..21e95178 100644 --- a/src/Model/PowerElectronics/MicrogridLine/MicrogridLine.hpp +++ b/src/Model/PowerElectronics/MicrogridLine/MicrogridLine.hpp @@ -3,7 +3,7 @@ #ifndef _TRANLINE_HPP_ #define _TRANLINE_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.hpp b/src/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.hpp index e292bbf1..53d98e6d 100644 --- a/src/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.hpp +++ b/src/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.hpp @@ -3,7 +3,7 @@ #ifndef _TRANLOAD_HPP_ #define _TRANLOAD_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/Resistor/Resistor.hpp b/src/Model/PowerElectronics/Resistor/Resistor.hpp index 822f5b48..942ec24b 100644 --- a/src/Model/PowerElectronics/Resistor/Resistor.hpp +++ b/src/Model/PowerElectronics/Resistor/Resistor.hpp @@ -3,7 +3,7 @@ #ifndef _RES_HPP_ #define _RES_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.hpp b/src/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.hpp index ec1d7486..08401004 100644 --- a/src/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.hpp +++ b/src/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.hpp @@ -3,7 +3,7 @@ #ifndef _SYNMACH_HPP_ #define _SYNMACH_HPP_ -#include + #include #include #include diff --git a/src/Model/PowerElectronics/SystemModelPowerElectronics.hpp b/src/Model/PowerElectronics/SystemModelPowerElectronics.hpp index 29fa037b..66edae0d 100644 --- a/src/Model/PowerElectronics/SystemModelPowerElectronics.hpp +++ b/src/Model/PowerElectronics/SystemModelPowerElectronics.hpp @@ -7,7 +7,6 @@ #include #include -#include #include #include @@ -15,31 +14,20 @@ namespace GridKit { template - class PowerElectronicsModel : public ModelEvaluatorImpl + class PowerElectronicsModel : public CircuitComponent { typedef CircuitComponent component_type; - using ModelEvaluatorImpl::size_; - // using ModelEvaluatorImpl::size_quad_; - // using ModelEvaluatorImpl::size_opt_; - using ModelEvaluatorImpl::nnz_; - using ModelEvaluatorImpl::time_; - using ModelEvaluatorImpl::alpha_; - using ModelEvaluatorImpl::y_; - using ModelEvaluatorImpl::yp_; - // using ModelEvaluatorImpl::yB_; - // using ModelEvaluatorImpl::ypB_; - // using ModelEvaluatorImpl::tag_; - using ModelEvaluatorImpl::f_; - // using ModelEvaluatorImpl::fB_; - // using ModelEvaluatorImpl::g_; - // using ModelEvaluatorImpl::gB_; - using ModelEvaluatorImpl::jac_; - using ModelEvaluatorImpl::rel_tol_; - using ModelEvaluatorImpl::abs_tol_; - // using ModelEvaluatorImpl::param_; - // using ModelEvaluatorImpl::param_up_; - // using ModelEvaluatorImpl::param_lo_; + using CircuitComponent::size_; + using CircuitComponent::nnz_; + using CircuitComponent::time_; + using CircuitComponent::alpha_; + using CircuitComponent::y_; + using CircuitComponent::yp_; + using CircuitComponent::f_; + using CircuitComponent::jac_; + using CircuitComponent::rel_tol_; + using CircuitComponent::abs_tol_; public: /** @@ -47,7 +35,7 @@ namespace GridKit * * @post System model parameters set as default */ - PowerElectronicsModel() : ModelEvaluatorImpl(0, 0, 0) + PowerElectronicsModel() { // Set system model parameters as default rel_tol_ = 1e-4; @@ -70,7 +58,7 @@ namespace GridKit PowerElectronicsModel(double rel_tol = 1e-4, double abs_tol = 1e-4, bool use_jac = false, - IdxT max_steps = 2000) : ModelEvaluatorImpl(0, 0, 0) + IdxT max_steps = 2000) { // Set system model tolerances from input rel_tol_ = rel_tol; @@ -97,13 +85,13 @@ namespace GridKit /** * @brief allocator default * - * @todo this should throw an exception as no allocation without a graph is allowed. Or needs to be removed from the base class + * @todo this should throw an exception as no allocation without a graph is allowed. + * Or needs to be removed from the base class * * @return int */ int allocate() { - return 1; } @@ -141,7 +129,6 @@ namespace GridKit */ int allocate(IdxT s) { - // Allocate all components size_ = s; for (const auto &component : components_) @@ -164,7 +151,6 @@ namespace GridKit */ int initialize() { - // Initialize components for (const auto &component : components_) { diff --git a/src/Model/PowerElectronics/TransmissionLine/TransmissionLine.hpp b/src/Model/PowerElectronics/TransmissionLine/TransmissionLine.hpp index 0ab0b340..2dd6f33a 100644 --- a/src/Model/PowerElectronics/TransmissionLine/TransmissionLine.hpp +++ b/src/Model/PowerElectronics/TransmissionLine/TransmissionLine.hpp @@ -3,7 +3,7 @@ #ifndef _TRANLINE_HPP_ #define _TRANLINE_HPP_ -#include + #include #include diff --git a/src/Model/PowerElectronics/VoltageSource/VoltageSource.hpp b/src/Model/PowerElectronics/VoltageSource/VoltageSource.hpp index 591865fa..972f4c46 100644 --- a/src/Model/PowerElectronics/VoltageSource/VoltageSource.hpp +++ b/src/Model/PowerElectronics/VoltageSource/VoltageSource.hpp @@ -3,7 +3,7 @@ #ifndef _VOSO_HPP_ #define _VOSO_HPP_ -#include + #include #include From ef096d91df6f6c6d88f0f3e6be598f8beafbd362 Mon Sep 17 00:00:00 2001 From: Slaven Peles Date: Thu, 20 Feb 2025 18:12:33 -0500 Subject: [PATCH 4/8] Move power flow specific files to Model/PowerFlow dir. --- examples/AdjointSensitivity/AdjointSensitivity.cpp | 2 +- examples/DynamicConstrainedOpt/DynamicConstrainedOpt.cpp | 2 +- examples/GenConstLoad/GenConstLoad.cpp | 2 +- examples/GenInfiniteBus/GenInfiniteBus.cpp | 2 +- examples/ParameterEstimation/ParameterEstimation.cpp | 2 +- src/Model/PowerFlow/Branch/Branch.hpp | 2 +- src/Model/PowerFlow/Bus/BaseBus.hpp | 2 +- src/Model/PowerFlow/Generator/GeneratorBase.hpp | 2 +- src/Model/PowerFlow/Generator/GeneratorPQ.hpp | 2 +- src/Model/PowerFlow/Generator/GeneratorPV.hpp | 2 +- src/Model/PowerFlow/Generator/GeneratorSlack.hpp | 2 +- src/Model/PowerFlow/Generator2/Generator2.hpp | 2 +- src/Model/PowerFlow/Generator4/Generator4.hpp | 2 +- src/Model/PowerFlow/Generator4Governor/Generator4Governor.hpp | 2 +- src/Model/PowerFlow/Generator4Param/Generator4Param.hpp | 2 +- src/Model/PowerFlow/Load/Load.hpp | 2 +- src/Model/PowerFlow/MiniGrid/MiniGrid.hpp | 2 +- src/{ => Model/PowerFlow}/ModelEvaluatorImpl.hpp | 0 src/{ => Model/PowerFlow}/SystemModel.hpp | 2 +- src/Model/PowerFlow/SystemModelPowerFlow.hpp | 2 +- 20 files changed, 19 insertions(+), 19 deletions(-) rename src/{ => Model/PowerFlow}/ModelEvaluatorImpl.hpp (100%) rename src/{ => Model/PowerFlow}/SystemModel.hpp (99%) diff --git a/examples/AdjointSensitivity/AdjointSensitivity.cpp b/examples/AdjointSensitivity/AdjointSensitivity.cpp index 1516b255..07af5351 100644 --- a/examples/AdjointSensitivity/AdjointSensitivity.cpp +++ b/examples/AdjointSensitivity/AdjointSensitivity.cpp @@ -63,7 +63,7 @@ #include #include -#include +#include #include #include diff --git a/examples/DynamicConstrainedOpt/DynamicConstrainedOpt.cpp b/examples/DynamicConstrainedOpt/DynamicConstrainedOpt.cpp index 9f3fc632..3d07650c 100644 --- a/examples/DynamicConstrainedOpt/DynamicConstrainedOpt.cpp +++ b/examples/DynamicConstrainedOpt/DynamicConstrainedOpt.cpp @@ -63,7 +63,7 @@ #include #include -#include +#include #include #include diff --git a/examples/GenConstLoad/GenConstLoad.cpp b/examples/GenConstLoad/GenConstLoad.cpp index 5f8a9e83..2a52677d 100644 --- a/examples/GenConstLoad/GenConstLoad.cpp +++ b/examples/GenConstLoad/GenConstLoad.cpp @@ -65,7 +65,7 @@ #include #include #include -#include +#include #include #include diff --git a/examples/GenInfiniteBus/GenInfiniteBus.cpp b/examples/GenInfiniteBus/GenInfiniteBus.cpp index e60ac888..a5f0e67e 100644 --- a/examples/GenInfiniteBus/GenInfiniteBus.cpp +++ b/examples/GenInfiniteBus/GenInfiniteBus.cpp @@ -63,7 +63,7 @@ #include #include -#include +#include #include #include diff --git a/examples/ParameterEstimation/ParameterEstimation.cpp b/examples/ParameterEstimation/ParameterEstimation.cpp index 599e7ff7..d8a683c1 100644 --- a/examples/ParameterEstimation/ParameterEstimation.cpp +++ b/examples/ParameterEstimation/ParameterEstimation.cpp @@ -69,7 +69,7 @@ #include #include #include -#include +#include #include #include diff --git a/src/Model/PowerFlow/Branch/Branch.hpp b/src/Model/PowerFlow/Branch/Branch.hpp index 5a8309c0..1b977881 100644 --- a/src/Model/PowerFlow/Branch/Branch.hpp +++ b/src/Model/PowerFlow/Branch/Branch.hpp @@ -60,7 +60,7 @@ #ifndef _BRANCH_H_ #define _BRANCH_H_ -#include +#include // Forward declarations. namespace GridKit diff --git a/src/Model/PowerFlow/Bus/BaseBus.hpp b/src/Model/PowerFlow/Bus/BaseBus.hpp index b3b2ed5c..e63347aa 100644 --- a/src/Model/PowerFlow/Bus/BaseBus.hpp +++ b/src/Model/PowerFlow/Bus/BaseBus.hpp @@ -60,7 +60,7 @@ #ifndef _BASE_BUS_HPP_ #define _BASE_BUS_HPP_ -#include +#include namespace GridKit { diff --git a/src/Model/PowerFlow/Generator/GeneratorBase.hpp b/src/Model/PowerFlow/Generator/GeneratorBase.hpp index 3044d939..29609a6f 100644 --- a/src/Model/PowerFlow/Generator/GeneratorBase.hpp +++ b/src/Model/PowerFlow/Generator/GeneratorBase.hpp @@ -59,7 +59,7 @@ #pragma once -#include +#include #include namespace GridKit diff --git a/src/Model/PowerFlow/Generator/GeneratorPQ.hpp b/src/Model/PowerFlow/Generator/GeneratorPQ.hpp index c0ecfcf9..2b3b8c65 100644 --- a/src/Model/PowerFlow/Generator/GeneratorPQ.hpp +++ b/src/Model/PowerFlow/Generator/GeneratorPQ.hpp @@ -60,7 +60,7 @@ #pragma once #include -#include +#include #include #include "GeneratorBase.hpp" diff --git a/src/Model/PowerFlow/Generator/GeneratorPV.hpp b/src/Model/PowerFlow/Generator/GeneratorPV.hpp index 88c1f8aa..68671a08 100644 --- a/src/Model/PowerFlow/Generator/GeneratorPV.hpp +++ b/src/Model/PowerFlow/Generator/GeneratorPV.hpp @@ -60,7 +60,7 @@ #pragma once #include -#include +#include #include #include "GeneratorBase.hpp" diff --git a/src/Model/PowerFlow/Generator/GeneratorSlack.hpp b/src/Model/PowerFlow/Generator/GeneratorSlack.hpp index 38cac14e..2b7fdea6 100644 --- a/src/Model/PowerFlow/Generator/GeneratorSlack.hpp +++ b/src/Model/PowerFlow/Generator/GeneratorSlack.hpp @@ -59,7 +59,7 @@ #pragma once -// #include +// #include #include "GeneratorBase.hpp" #include #include diff --git a/src/Model/PowerFlow/Generator2/Generator2.hpp b/src/Model/PowerFlow/Generator2/Generator2.hpp index 8f352339..51aaf28f 100644 --- a/src/Model/PowerFlow/Generator2/Generator2.hpp +++ b/src/Model/PowerFlow/Generator2/Generator2.hpp @@ -57,7 +57,7 @@ * */ -#include +#include namespace GridKit { diff --git a/src/Model/PowerFlow/Generator4/Generator4.hpp b/src/Model/PowerFlow/Generator4/Generator4.hpp index 9af37c6e..1ceee9f8 100644 --- a/src/Model/PowerFlow/Generator4/Generator4.hpp +++ b/src/Model/PowerFlow/Generator4/Generator4.hpp @@ -60,7 +60,7 @@ #ifndef _GENERATOR_4_H_ #define _GENERATOR_4_H_ -#include +#include namespace GridKit { diff --git a/src/Model/PowerFlow/Generator4Governor/Generator4Governor.hpp b/src/Model/PowerFlow/Generator4Governor/Generator4Governor.hpp index f21be693..68ac1198 100644 --- a/src/Model/PowerFlow/Generator4Governor/Generator4Governor.hpp +++ b/src/Model/PowerFlow/Generator4Governor/Generator4Governor.hpp @@ -60,7 +60,7 @@ #ifndef _GENERATOR_4_GOVERNOR_B_HPP_ #define _GENERATOR_4_GOVERNOR_B_HPP_ -#include +#include namespace GridKit { diff --git a/src/Model/PowerFlow/Generator4Param/Generator4Param.hpp b/src/Model/PowerFlow/Generator4Param/Generator4Param.hpp index 1a2ba327..430a601f 100644 --- a/src/Model/PowerFlow/Generator4Param/Generator4Param.hpp +++ b/src/Model/PowerFlow/Generator4Param/Generator4Param.hpp @@ -60,7 +60,7 @@ #ifndef _GENERATOR_4_H_ #define _GENERATOR_4_H_ -#include +#include namespace GridKit { diff --git a/src/Model/PowerFlow/Load/Load.hpp b/src/Model/PowerFlow/Load/Load.hpp index 4586cb46..55989c11 100644 --- a/src/Model/PowerFlow/Load/Load.hpp +++ b/src/Model/PowerFlow/Load/Load.hpp @@ -60,7 +60,7 @@ #ifndef _LOAD_HPP_ #define _LOAD_HPP_ -#include +#include #include namespace GridKit diff --git a/src/Model/PowerFlow/MiniGrid/MiniGrid.hpp b/src/Model/PowerFlow/MiniGrid/MiniGrid.hpp index 5cd33bd6..1a46d008 100644 --- a/src/Model/PowerFlow/MiniGrid/MiniGrid.hpp +++ b/src/Model/PowerFlow/MiniGrid/MiniGrid.hpp @@ -58,7 +58,7 @@ */ #pragma once -#include +#include #include namespace GridKit diff --git a/src/ModelEvaluatorImpl.hpp b/src/Model/PowerFlow/ModelEvaluatorImpl.hpp similarity index 100% rename from src/ModelEvaluatorImpl.hpp rename to src/Model/PowerFlow/ModelEvaluatorImpl.hpp diff --git a/src/SystemModel.hpp b/src/Model/PowerFlow/SystemModel.hpp similarity index 99% rename from src/SystemModel.hpp rename to src/Model/PowerFlow/SystemModel.hpp index b9139d25..90588a21 100644 --- a/src/SystemModel.hpp +++ b/src/Model/PowerFlow/SystemModel.hpp @@ -65,7 +65,7 @@ #include #include -#include +#include namespace GridKit { diff --git a/src/Model/PowerFlow/SystemModelPowerFlow.hpp b/src/Model/PowerFlow/SystemModelPowerFlow.hpp index 1a7bd214..f380c4cf 100644 --- a/src/Model/PowerFlow/SystemModelPowerFlow.hpp +++ b/src/Model/PowerFlow/SystemModelPowerFlow.hpp @@ -71,7 +71,7 @@ #include #include -#include +#include namespace GridKit { From 42f0be8a1e07cc525bf16bf1644deb9eb488644d Mon Sep 17 00:00:00 2001 From: Slaven Peles Date: Thu, 20 Feb 2025 19:07:33 -0500 Subject: [PATCH 5/8] Add phasor dynamics system composer. --- src/Model/PhasorDynamics/Component.hpp | 8 +- src/Model/PhasorDynamics/SystemModel.hpp | 659 ++++++++++++++++++ tests/UnitTests/PhasorDynamics/CMakeLists.txt | 7 +- .../UnitTests/PhasorDynamics/SystemTests.hpp | 50 ++ .../PhasorDynamics/runSystemTests.cpp | 14 + 5 files changed, 733 insertions(+), 5 deletions(-) create mode 100644 src/Model/PhasorDynamics/SystemModel.hpp create mode 100644 tests/UnitTests/PhasorDynamics/SystemTests.hpp create mode 100644 tests/UnitTests/PhasorDynamics/runSystemTests.cpp diff --git a/src/Model/PhasorDynamics/Component.hpp b/src/Model/PhasorDynamics/Component.hpp index 3a724a02..25e328ef 100644 --- a/src/Model/PhasorDynamics/Component.hpp +++ b/src/Model/PhasorDynamics/Component.hpp @@ -77,8 +77,8 @@ namespace PhasorDynamics virtual void setTolerances(real_type& rtol, real_type& atol) const { - rtol = rtol_; - atol = atol_; + rtol = rel_tol_; + atol = abs_tol_; } virtual void setMaxSteps(IdxT& msa) const @@ -250,8 +250,8 @@ namespace PhasorDynamics real_type time_; real_type alpha_; - real_type rtol_; - real_type atol_; + real_type rel_tol_; + real_type abs_tol_; IdxT max_steps_; diff --git a/src/Model/PhasorDynamics/SystemModel.hpp b/src/Model/PhasorDynamics/SystemModel.hpp new file mode 100644 index 00000000..b6a4f3af --- /dev/null +++ b/src/Model/PhasorDynamics/SystemModel.hpp @@ -0,0 +1,659 @@ +#pragma once + +#include +#include +#include + +#include +#include +#include + +namespace GridKit +{ +namespace PhasorDynamics +{ + +/** + * @brief Prototype for a system model class + * + * This class maps component data to system data and implements + * Model::Evaluator for the system model. This is still work in + * progress and code is not optimized. + * + * @todo Address thread safety for the system model methods. + * + */ +template +class SystemModel : public PhasorDynamics::Component +{ + using bus_type = PhasorDynamics::BusBase; + using component_type = PhasorDynamics::Component; + using real_type = typename Model::Evaluator::real_type; + + using PhasorDynamics::Component::size_; + using PhasorDynamics::Component::size_quad_; + using PhasorDynamics::Component::size_param_; + using PhasorDynamics::Component::nnz_; + using PhasorDynamics::Component::time_; + using PhasorDynamics::Component::alpha_; + using PhasorDynamics::Component::y_; + using PhasorDynamics::Component::yp_; + using PhasorDynamics::Component::yB_; + using PhasorDynamics::Component::ypB_; + using PhasorDynamics::Component::tag_; + using PhasorDynamics::Component::f_; + using PhasorDynamics::Component::fB_; + using PhasorDynamics::Component::g_; + using PhasorDynamics::Component::gB_; + using PhasorDynamics::Component::rel_tol_; + using PhasorDynamics::Component::abs_tol_; + using PhasorDynamics::Component::param_; + using PhasorDynamics::Component::param_up_; + using PhasorDynamics::Component::param_lo_; + +public: + /** + * @brief Constructor for the system model + */ + SystemModel() + { + // Set system model tolerances + rel_tol_ = 1e-7; + abs_tol_ = 1e-9; + this->max_steps_=2000; + } + + /** + * @brief Destructor for the system model + */ + virtual ~SystemModel() + { + } + + /** + * @brief Allocate buses, components, and system objects. + * + * This method first allocates bus objects, then component objects, + * and computes system size (number of unknowns). Once the size is + * computed, system global objects are allocated. + * + * @post size_quad_ == 0 or 1 + * @post size_ >= 1 + * @post size_param_ >= 0 + * + */ + int allocate() + { + size_ = 0; + size_quad_ = 0; + size_param_ = 0; + + // Allocate all buses + for(const auto& bus: buses_) + { + bus->allocate(); + size_ += bus->size(); + size_quad_ += bus->sizeQuadrature(); + size_param_ += bus->sizeParams(); + } + + // Allocate all components + for(const auto& component : components_) + { + component->allocate(); + size_ += component->size(); + size_quad_ += component->sizeQuadrature(); + size_param_ += component->sizeParams(); + } + + // Allocate global vectors + y_.resize(size_); + yp_.resize(size_); + yB_.resize(size_); + ypB_.resize(size_); + f_.resize(size_); + fB_.resize(size_); + tag_.resize(size_); + + g_.resize(size_quad_); + gB_.resize(size_quad_*size_param_); + + param_.resize(size_param_); + param_lo_.resize(size_param_); + param_up_.resize(size_param_); + + assert(size_quad_ == 1 or size_quad_ == 0); + + return 0; + } + + /** + * @brief Assume that jacobian is not avalible + * + * @return true + * @return false + */ + bool hasJacobian() + { + return false; + } + + /** + * @brief Initialize buses first, then all the other components. + * + * @pre All buses and components must be allocated at this point. + * @pre Bus variables are written before component variables in the + * system variable vector. + * + * Buses must be initialized before other components, because other + * components may write to buses during the initialization. + * + * Also, generators may write to control devices (e.g. governors, + * exciters, etc.) during the initialization. + * + * @todo Implement writting to system vectors in a thread-safe way. + */ + int initialize() + { + // Set initial values for global solution vectors + IdxT varOffset = 0; + IdxT optOffset = 0; + + for(const auto& bus: buses_) + { + bus->initialize(); + } + + for(const auto& bus: buses_) + { + for(IdxT j=0; jsize(); ++j) + { + y_[varOffset + j] = bus->y()[j]; + yp_[varOffset + j] = bus->yp()[j]; + } + varOffset += bus->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + param_[optOffset + j] = bus->param()[j]; + param_lo_[optOffset + j] = bus->param_lo()[j]; + param_up_[optOffset + j] = bus->param_up()[j]; + } + optOffset += bus->sizeParams(); + } + + // Initialize components + for(const auto& component : components_) + { + component->initialize(); + } + + for(const auto& component : components_) + { + for(IdxT j=0; jsize(); ++j) + { + y_[varOffset + j] = component->y()[j]; + yp_[varOffset + j] = component->yp()[j]; + } + varOffset += component->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + param_[optOffset + j] = component->param()[j]; + param_lo_[optOffset + j] = component->param_lo()[j]; + param_up_[optOffset + j] = component->param_up()[j]; + } + optOffset += component->sizeParams(); + } + + return 0; + } + + /** + * @todo Tagging differential variables + * + * Identify what variables in the system of differential-algebraic + * equations are differential variables, i.e. their derivatives + * appear in the equations. + */ + int tagDifferentiable() + { + // Set initial values for global solution vectors + IdxT offset = 0; + for(const auto& bus: buses_) + { + bus->tagDifferentiable(); + for(IdxT j=0; jsize(); ++j) + { + tag_[offset + j] = bus->tag()[j]; + } + offset += bus->size(); + } + + for(const auto& component: components_) + { + component->tagDifferentiable(); + for(IdxT j=0; jsize(); ++j) + { + tag_[offset + j] = component->tag()[j]; + } + offset += component->size(); + } + + return 0; + } + + /** + * @brief Compute system residual vector + * + * First, update bus and component variables from the system solution + * vector. Next, evaluate residuals in buses and components, and + * then copy values to the global residual vector. + * + * @warning Residuals must be computed for buses, before component + * residuals are computed. Buses own residuals for active and + * power P and Q, but the contributions to these residuals come + * from components. Buses assign their residual values, while components + * add to those values by in-place adition. This is why bus residuals + * need to be computed first. + * + * @todo Here, components write to local values, which are then copied + * to global system vectors. Make components write to the system + * vectors directly. + */ + int evaluateResidual() + { + // Update variables + IdxT varOffset = 0; + IdxT optOffset = 0; + for(const auto& bus: buses_) + { + for(IdxT j=0; jsize(); ++j) + { + bus->y()[j] = y_[varOffset + j]; + bus->yp()[j] = yp_[varOffset + j]; + } + varOffset += bus->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + bus->param()[j] = param_[optOffset + j]; + } + optOffset += bus->sizeParams(); + + bus->evaluateResidual(); + } + + for(const auto& component : components_) + { + for(IdxT j=0; jsize(); ++j) + { + component->y()[j] = y_[varOffset + j]; + component->yp()[j] = yp_[varOffset + j]; + } + varOffset += component->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + component->param()[j] = param_[optOffset + j]; + } + optOffset += component->sizeParams(); + + component->evaluateResidual(); + } + + // Update residual vector + IdxT resOffset = 0; + for(const auto& bus: buses_) + { + for(IdxT j=0; jsize(); ++j) + { + f_[resOffset + j] = bus->getResidual()[j]; + } + resOffset += bus->size(); + } + + for(const auto& component : components_) + { + for(IdxT j=0; jsize(); ++j) + { + f_[resOffset + j] = component->getResidual()[j]; + } + resOffset += component->size(); + } + + return 0; + } + + /** + * @brief Evaluate system Jacobian. + * + * @todo Need to implement Jacobian. For now, using finite difference + * approximation provided by IDA. This works for dense Jacobian matrix + * only. + * + */ + int evaluateJacobian(){return 0;} + + /** + * @brief Evaluate integrands for the system quadratures. + */ + int evaluateIntegrand() + { + // Update variables + IdxT varOffset = 0; + IdxT optOffset = 0; + for(const auto& bus: buses_) + { + for(IdxT j=0; jsize(); ++j) + { + bus->y()[j] = y_[varOffset + j]; + bus->yp()[j] = yp_[varOffset + j]; + } + varOffset += bus->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + bus->param()[j] = param_[optOffset + j]; + } + optOffset += bus->sizeParams(); + + bus->evaluateIntegrand(); + } + + for(const auto& component : components_) + { + for(IdxT j=0; jsize(); ++j) + { + component->y()[j] = y_[varOffset + j]; + component->yp()[j] = yp_[varOffset + j]; + } + varOffset += component->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + component->param()[j] = param_[optOffset + j]; + } + optOffset += component->sizeParams(); + + component->evaluateIntegrand(); + } + + // Update integrand vector + IdxT intOffset = 0; + for(const auto& bus: buses_) + { + for(IdxT j=0; jsizeQuadrature(); ++j) + { + g_[intOffset + j] = bus->getIntegrand()[j]; + } + intOffset += bus->sizeQuadrature(); + } + + for(const auto& component : components_) + { + for(IdxT j=0; jsizeQuadrature(); ++j) + { + g_[intOffset + j] = component->getIntegrand()[j]; + } + intOffset += component->sizeQuadrature(); + } + + return 0; + } + + /** + * @brief Initialize system adjoint. + * + * Updates variables and optimization parameters, then initializes + * adjoints locally and copies them to the system adjoint vector. + */ + int initializeAdjoint() + { + IdxT offset = 0; + IdxT optOffset = 0; + + // Update bus variables and optimization parameters + for(const auto& bus: buses_) + { + for(IdxT j=0; jsize(); ++j) + { + bus->y()[j] = y_[offset + j]; + bus->yp()[j] = yp_[offset + j]; + } + offset += bus->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + bus->param()[j] = param_[optOffset + j]; + } + optOffset += bus->sizeParams(); + } + + // Update component variables and optimization parameters + for(const auto& component: components_) + { + for(IdxT j=0; jsize(); ++j) + { + component->y()[j] = y_[offset + j]; + component->yp()[j] = yp_[offset + j]; + } + offset += component->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + component->param()[j] = param_[optOffset + j]; + } + optOffset += component->sizeParams(); + } + + // Reset counter + offset = 0; + + // Initialize bus adjoints + for(const auto& bus: buses_) + { + bus->initializeAdjoint(); + + for(IdxT j=0; jsize(); ++j) + { + yB_[offset + j] = bus->yB()[j]; + ypB_[offset + j] = bus->ypB()[j]; + } + offset += bus->size(); + } + + // Initialize component adjoints + for(const auto& component: components_) + { + component->initializeAdjoint(); + + for(IdxT j=0; jsize(); ++j) + { + yB_[offset + j] = component->yB()[j]; + ypB_[offset + j] = component->ypB()[j]; + } + offset += component->size(); + } + + return 0; + } + + /** + * @brief Compute adjoint residual for the system model. + * + * @warning Components write to bus residuals. Do not copy bus residuals + * to system vectors before components computed their residuals. + * + */ + int evaluateAdjointResidual() + { + IdxT varOffset = 0; + IdxT optOffset = 0; + + // Update variables in component models + for(const auto& bus: buses_) + { + for(IdxT j=0; jsize(); ++j) + { + bus->y()[j] = y_[varOffset + j]; + bus->yp()[j] = yp_[varOffset + j]; + bus->yB()[j] = yB_[varOffset + j]; + bus->ypB()[j] = ypB_[varOffset + j]; + } + varOffset += bus->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + bus->param()[j] = param_[optOffset + j]; + } + optOffset += bus->sizeParams(); + + } + + for(const auto& component : components_) + { + for(IdxT j=0; jsize(); ++j) + { + component->y()[j] = y_[varOffset + j]; + component->yp()[j] = yp_[varOffset + j]; + component->yB()[j] = yB_[varOffset + j]; + component->ypB()[j] = ypB_[varOffset + j]; + } + varOffset += component->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + component->param()[j] = param_[optOffset + j]; + } + optOffset += component->sizeParams(); + + } + + for(const auto& bus: buses_) + { + bus->evaluateAdjointResidual(); + } + + for(const auto& component : components_) + { + component->evaluateAdjointResidual(); + } + + // Update residual vector + IdxT resOffset = 0; + for(const auto& bus: buses_) + { + for(IdxT j=0; jsize(); ++j) + { + fB_[resOffset + j] = bus->getAdjointResidual()[j]; + } + resOffset += bus->size(); + } + + for(const auto& component : components_) + { + for(IdxT j=0; jsize(); ++j) + { + fB_[resOffset + j] = component->getAdjointResidual()[j]; + } + resOffset += component->size(); + } + + return 0; + } + + //int evaluateAdjointJacobian(){return 0;} + + /** + * @brief Evaluate adjoint integrand for the system model. + * + * @pre Assumes there are no integrands in bus models. + * @pre Assumes integrand is implemented in only _one_ component. + * + */ + int evaluateAdjointIntegrand() + { + // First, update variables + IdxT varOffset = 0; + IdxT optOffset = 0; + for(const auto& bus: buses_) + { + for(IdxT j=0; jsize(); ++j) + { + bus->y()[j] = y_[varOffset + j]; + bus->yp()[j] = yp_[varOffset + j]; + bus->yB()[j] = yB_[varOffset + j]; + bus->ypB()[j] = ypB_[varOffset + j]; + } + varOffset += bus->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + bus->param()[j] = param_[optOffset + j]; + } + optOffset += bus->sizeParams(); + } + + for(const auto& component : components_) + { + for(IdxT j=0; jsize(); ++j) + { + component->y()[j] = y_[varOffset + j]; + component->yp()[j] = yp_[varOffset + j]; + component->yB()[j] = yB_[varOffset + j]; + component->ypB()[j] = ypB_[varOffset + j]; + } + varOffset += component->size(); + + for(IdxT j=0; jsizeParams(); ++j) + { + component->param()[j] = param_[optOffset + j]; + } + optOffset += component->sizeParams(); + } + + // Evaluate integrand and update global vector + for(const auto& component : components_) + { + if(component->sizeQuadrature() == 1) + { + component->evaluateAdjointIntegrand(); + for(IdxT j=0; jgetAdjointIntegrand()[j]; + } + break; + } + } + return 0; + } + + void updateTime(real_type t, real_type a) + { + for(const auto& component : components_) + { + component->updateTime(t, a); + } + } + + void addBus(bus_type* bus) + { + buses_.push_back(bus); + } + + void addComponent(component_type* component) + { + components_.push_back(component); + } + +private: + std::vector buses_; + std::vector components_; + + + +}; // class SystemModel + +} // namespace PhasorDynamics +} // namespace GridKit diff --git a/tests/UnitTests/PhasorDynamics/CMakeLists.txt b/tests/UnitTests/PhasorDynamics/CMakeLists.txt index 30cd2bbb..5e6be06f 100644 --- a/tests/UnitTests/PhasorDynamics/CMakeLists.txt +++ b/tests/UnitTests/PhasorDynamics/CMakeLists.txt @@ -15,8 +15,13 @@ add_executable(test_phasor_load runLoadTests.cpp) target_link_libraries(test_phasor_load GRIDKIT::phasor_dynamics_load GRIDKIT::phasor_dynamics_bus) +add_executable(test_phasor_system runSystemTests.cpp) +target_link_libraries(test_phasor_system GRIDKIT::phasor_dynamics_load + GRIDKIT::phasor_dynamics_bus) + add_test(NAME PhasorDynamicsBusTest COMMAND $) add_test(NAME PhasorDynamicsBranchTest COMMAND $) add_test(NAME PhasorDynamicsLoadTest COMMAND $) +add_test(NAME PhasorDynamicsSystemTest COMMAND $) -install(TARGETS test_phasor_bus test_phasor_branch test_phasor_load RUNTIME DESTINATION bin) +install(TARGETS test_phasor_bus test_phasor_branch test_phasor_load test_phasor_system RUNTIME DESTINATION bin) diff --git a/tests/UnitTests/PhasorDynamics/SystemTests.hpp b/tests/UnitTests/PhasorDynamics/SystemTests.hpp new file mode 100644 index 00000000..829155fd --- /dev/null +++ b/tests/UnitTests/PhasorDynamics/SystemTests.hpp @@ -0,0 +1,50 @@ +#pragma once + +#include +#include + +#include +#include +#include +#include +#include + + +namespace GridKit +{ +namespace Testing +{ + template + class SystemTests + { + public: + SystemTests() = default; + ~SystemTests() = default; + + /// Constructor, allocation, and initialization checks + TestOutcome constructor() + { + TestStatus success = true; + + // ScalarT Vr{1.0}; + // ScalarT Vi{2.0}; + + PhasorDynamics::SystemModel* system = nullptr; + + // Create an empty system + system = new PhasorDynamics::SystemModel(); + + success *= (system != nullptr); + + if (system) + { + delete system; + } + + return success.report(__func__); + } + + }; + +} // namespace Testing +} // namespace GridKit diff --git a/tests/UnitTests/PhasorDynamics/runSystemTests.cpp b/tests/UnitTests/PhasorDynamics/runSystemTests.cpp new file mode 100644 index 00000000..35ac6a73 --- /dev/null +++ b/tests/UnitTests/PhasorDynamics/runSystemTests.cpp @@ -0,0 +1,14 @@ +#include "SystemTests.hpp" + +int main() +{ + using namespace GridKit; + using namespace GridKit::Testing; + + GridKit::Testing::TestingResults result; + GridKit::Testing::SystemTests test; + + result += test.constructor(); + + return result.summary(); +} \ No newline at end of file From 429b5c2d2de5d38803170e434360386c64d2c5f3 Mon Sep 17 00:00:00 2001 From: Slaven Peles Date: Fri, 21 Feb 2025 14:32:29 -0500 Subject: [PATCH 6/8] Add test for system composer. --- tests/UnitTests/PhasorDynamics/CMakeLists.txt | 1 + .../UnitTests/PhasorDynamics/SystemTests.hpp | 49 +++++++++++++++++++ .../PhasorDynamics/runSystemTests.cpp | 1 + 3 files changed, 51 insertions(+) diff --git a/tests/UnitTests/PhasorDynamics/CMakeLists.txt b/tests/UnitTests/PhasorDynamics/CMakeLists.txt index 5e6be06f..0a33667a 100644 --- a/tests/UnitTests/PhasorDynamics/CMakeLists.txt +++ b/tests/UnitTests/PhasorDynamics/CMakeLists.txt @@ -17,6 +17,7 @@ target_link_libraries(test_phasor_load GRIDKIT::phasor_dynamics_load add_executable(test_phasor_system runSystemTests.cpp) target_link_libraries(test_phasor_system GRIDKIT::phasor_dynamics_load + GRIDKIT::phasor_dynamics_branch GRIDKIT::phasor_dynamics_bus) add_test(NAME PhasorDynamicsBusTest COMMAND $) diff --git a/tests/UnitTests/PhasorDynamics/SystemTests.hpp b/tests/UnitTests/PhasorDynamics/SystemTests.hpp index 829155fd..01deb34d 100644 --- a/tests/UnitTests/PhasorDynamics/SystemTests.hpp +++ b/tests/UnitTests/PhasorDynamics/SystemTests.hpp @@ -5,6 +5,7 @@ #include #include +#include #include #include #include @@ -17,6 +18,9 @@ namespace Testing template class SystemTests { + private: + using real_type = typename PhasorDynamics::Component::real_type; + public: SystemTests() = default; ~SystemTests() = default; @@ -44,6 +48,51 @@ namespace Testing return success.report(__func__); } + TestOutcome composer() + { + TestStatus success = true; + + real_type R{2.0}; ///< Branch series resistance + real_type X{4.0}; ///< Branch series reactance + real_type G{0.2}; ///< Branch shunt conductance + real_type B{1.2}; ///< Branch shunt charging + + ScalarT Vr1{10.0}; ///< Bus-1 real voltage + ScalarT Vi1{20.0}; ///< Bus-1 imaginary voltage + ScalarT Vr2{30.0}; ///< Bus-2 real voltage + ScalarT Vi2{40.0}; ///< Bus-2 imaginary voltage + + const ScalarT Ir1{17.0}; ///< Solution: real current entering bus-1 + const ScalarT Ii1{-10.0}; ///< Solution: imaginary current entering bus-1 + const ScalarT Ir2{15.0}; ///< Solution: real current entering bus-2 + const ScalarT Ii2{-20.0}; ///< Solution: imaginary current entering bus-2 + + // Create an empty system model + PhasorDynamics::SystemModel system; + + // Add a bus + PhasorDynamics::BusInfinite bus1(Vr1, Vi1); + system.addBus(&bus1); + + // Add a bus + PhasorDynamics::BusInfinite bus2(Vr2, Vi2); + system.addBus(&bus1); + + PhasorDynamics::Branch branch(&bus1, &bus2, R, X, G, B); + system.addComponent(&branch); + + system.allocate(); + system.initialize(); + system.evaluateResidual(); + + success *= isEqual(bus1.Ir(), Ir1); + success *= isEqual(bus1.Ii(), Ii1); + success *= isEqual(bus2.Ir(), Ir2); + success *= isEqual(bus2.Ii(), Ii2); + + return success.report(__func__); + } + }; } // namespace Testing diff --git a/tests/UnitTests/PhasorDynamics/runSystemTests.cpp b/tests/UnitTests/PhasorDynamics/runSystemTests.cpp index 35ac6a73..be0335de 100644 --- a/tests/UnitTests/PhasorDynamics/runSystemTests.cpp +++ b/tests/UnitTests/PhasorDynamics/runSystemTests.cpp @@ -9,6 +9,7 @@ int main() GridKit::Testing::SystemTests test; result += test.constructor(); + result += test.composer(); return result.summary(); } \ No newline at end of file From d942897384cf2983144e272a8a4654f372399340 Mon Sep 17 00:00:00 2001 From: Slaven Peles Date: Fri, 21 Feb 2025 15:05:48 -0500 Subject: [PATCH 7/8] Add template for implementing synchronous machine model. --- .../PhasorDynamics/Branch/CMakeLists.txt | 1 + src/Model/PhasorDynamics/CMakeLists.txt | 1 + .../SynchronousMachine/CMakeLists.txt | 13 ++ .../SynchronousMachine/SynchronousMachine.cpp | 173 ++++++++++++++++++ .../SynchronousMachine/SynchronousMachine.hpp | 127 +++++++++++++ 5 files changed, 315 insertions(+) create mode 100644 src/Model/PhasorDynamics/SynchronousMachine/CMakeLists.txt create mode 100644 src/Model/PhasorDynamics/SynchronousMachine/SynchronousMachine.cpp create mode 100644 src/Model/PhasorDynamics/SynchronousMachine/SynchronousMachine.hpp diff --git a/src/Model/PhasorDynamics/Branch/CMakeLists.txt b/src/Model/PhasorDynamics/Branch/CMakeLists.txt index b3d400e2..71478c50 100644 --- a/src/Model/PhasorDynamics/Branch/CMakeLists.txt +++ b/src/Model/PhasorDynamics/Branch/CMakeLists.txt @@ -2,6 +2,7 @@ # [[ # Author(s): # - Cameron Rutherford +# - Slaven Peles # ]] gridkit_add_library(phasor_dynamics_branch diff --git a/src/Model/PhasorDynamics/CMakeLists.txt b/src/Model/PhasorDynamics/CMakeLists.txt index 08f9613e..77dc2d3e 100644 --- a/src/Model/PhasorDynamics/CMakeLists.txt +++ b/src/Model/PhasorDynamics/CMakeLists.txt @@ -6,3 +6,4 @@ add_subdirectory(Branch) add_subdirectory(Bus) add_subdirectory(Load) +add_subdirectory(SynchronousMachine) diff --git a/src/Model/PhasorDynamics/SynchronousMachine/CMakeLists.txt b/src/Model/PhasorDynamics/SynchronousMachine/CMakeLists.txt new file mode 100644 index 00000000..49216d88 --- /dev/null +++ b/src/Model/PhasorDynamics/SynchronousMachine/CMakeLists.txt @@ -0,0 +1,13 @@ + +# [[ +# Author(s): +# - Cameron Rutherford +# - Slaven Peles +# ]] + +gridkit_add_library(phasor_dynamics_synchronous_machine + SOURCES + SynchronousMachine.cpp + OUTPUT_NAME + gridkit_phasor_dynamics_synchronous_machine) + diff --git a/src/Model/PhasorDynamics/SynchronousMachine/SynchronousMachine.cpp b/src/Model/PhasorDynamics/SynchronousMachine/SynchronousMachine.cpp new file mode 100644 index 00000000..e30adb28 --- /dev/null +++ b/src/Model/PhasorDynamics/SynchronousMachine/SynchronousMachine.cpp @@ -0,0 +1,173 @@ +/** + * @file SynchronousMachine.hpp + * @author Slaven Peles (peless@ornl.gov) + * @brief Definition of a phasor dynamics branch model. + * + * The model uses Cartesian coordinates. + * + */ + +#include +#include +#include +#include + +#include "SynchronousMachine.hpp" + +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 + SynchronousMachine::SynchronousMachine(bus_type* bus) + : bus_(bus), + R_(0.0), + X_(0.01), + G_(0.0), + B_(0.0) + { + size_ = 0; + } + + /** + * @brief Destroy the SynchronousMachine + * + * @tparam ScalarT + * @tparam IdxT + */ + template + SynchronousMachine::~SynchronousMachine() + { + //std::cout << "Destroy SynchronousMachine..." << std::endl; + } + + /*! + * @brief allocate method computes sparsity pattern of the Jacobian. + */ + template + int SynchronousMachine::allocate() + { + //std::cout << "Allocate SynchronousMachine..." << std::endl; + return 0; + } + + /** + * Initialization of the branch model + * + */ + template + int SynchronousMachine::initialize() + { + return 0; + } + + /** + * \brief Identify differential variables. + */ + template + int SynchronousMachine::tagDifferentiable() + { + return 0; + } + + /** + * \brief Residual contribution of the branch is pushed to the + * two terminal buses. + * + */ + template + int SynchronousMachine::evaluateResidual() + { + // std::cout << "Evaluating branch residual ...\n"; + // real_type b = -X_/(R_*R_ + X_*X_); + // real_type g = R_/(R_*R_ + X_*X_); + + 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 SynchronousMachine::evaluateJacobian() + { + std::cout << "Evaluate Jacobian for SynchronousMachine..." << 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 SynchronousMachine::evaluateIntegrand() + { + // std::cout << "Evaluate Integrand for SynchronousMachine..." << 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 SynchronousMachine::initializeAdjoint() + { + //std::cout << "Initialize adjoint for SynchronousMachine..." << 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 SynchronousMachine::evaluateAdjointResidual() + { + // std::cout << "Evaluate adjoint residual for SynchronousMachine..." << 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 SynchronousMachine::evaluateAdjointIntegrand() + { + // std::cout << "Evaluate adjoint Integrand for SynchronousMachine..." << std::endl; + return 0; + } + + // Available template instantiations + template class SynchronousMachine; + template class SynchronousMachine; + +} //namespace PhasorDynamics +} //namespace GridKit diff --git a/src/Model/PhasorDynamics/SynchronousMachine/SynchronousMachine.hpp b/src/Model/PhasorDynamics/SynchronousMachine/SynchronousMachine.hpp new file mode 100644 index 00000000..c3818cc5 --- /dev/null +++ b/src/Model/PhasorDynamics/SynchronousMachine/SynchronousMachine.hpp @@ -0,0 +1,127 @@ +/** + * @file SynchronousMachine.hpp + * @author Slaven Peles (peless@ornl.gov) + * @brief Declaration of a phasor dynamics branch model. + * + * The model uses Cartesian coordinates. + * + */ +#pragma once + +#include + +// Forward declarations. +namespace GridKit +{ +namespace PhasorDynamics +{ + template class BusBase; +} +} + +namespace GridKit +{ +namespace PhasorDynamics +{ + /** + * @brief Implementation of a pi-model branch between two buses. + * + * The model is implemented in Cartesian coordinates. Positive current + * direction is into the busses. + * + */ + template + class SynchronousMachine : public Component + { + using Component::size_; + using Component::nnz_; + using Component::time_; + using Component::alpha_; + using Component::y_; + using Component::yp_; + using Component::tag_; + using Component::f_; + using Component::g_; + using Component::yB_; + using Component::ypB_; + using Component::fB_; + using Component::gB_; + using Component::param_; + + using bus_type = BusBase; + using real_type = typename Component::real_type; + + public: + SynchronousMachine(bus_type* bus); + virtual ~SynchronousMachine(); + + virtual int allocate() override; + virtual int initialize() override; + virtual int tagDifferentiable() override; + virtual int evaluateResidual() override; + virtual int evaluateJacobian() override; + virtual int evaluateIntegrand() override; + + virtual int initializeAdjoint() override; + virtual int evaluateAdjointResidual() override; + // virtual int evaluateAdjointJacobian() override; + virtual int evaluateAdjointIntegrand() override; + + virtual void updateTime(real_type /* t */, real_type /* a */) override + { + } + + public: + void setR(real_type R) + { + R_ = R; + } + + void setX(real_type X) + { + // std::cout << "Setting X ...\n"; + X_ = X; + } + + void setG(real_type G) + { + G_ = G; + } + + void setB(real_type B) + { + B_ = B; + } + + private: + ScalarT& Vr() + { + return bus_->Vr(); + } + + ScalarT& Vi() + { + return bus_->Vi(); + } + + ScalarT& Ir() + { + return bus_->Ir(); + } + + ScalarT& Ii() + { + return bus_->Ii(); + } + + + private: + bus_type* bus_; + real_type R_; + real_type X_; + real_type G_; + real_type B_; + }; + +} // namespace PhasorDynamics +} // namespace GridKit From 3ca3b90381bfa38fcaee644aac2e3344820fc76f Mon Sep 17 00:00:00 2001 From: Slaven Peles Date: Fri, 21 Feb 2025 15:46:32 -0500 Subject: [PATCH 8/8] Add testing framework for the synchronous machine model --- tests/UnitTests/PhasorDynamics/CMakeLists.txt | 14 +++- .../SynchronousMachineTests.hpp | 67 +++++++++++++++++++ .../runSynchronousMachineTests.cpp | 16 +++++ 3 files changed, 95 insertions(+), 2 deletions(-) create mode 100644 tests/UnitTests/PhasorDynamics/SynchronousMachineTests.hpp create mode 100644 tests/UnitTests/PhasorDynamics/runSynchronousMachineTests.cpp diff --git a/tests/UnitTests/PhasorDynamics/CMakeLists.txt b/tests/UnitTests/PhasorDynamics/CMakeLists.txt index 0a33667a..7077c129 100644 --- a/tests/UnitTests/PhasorDynamics/CMakeLists.txt +++ b/tests/UnitTests/PhasorDynamics/CMakeLists.txt @@ -14,7 +14,12 @@ target_link_libraries(test_phasor_branch GRIDKIT::phasor_dynamics_branch add_executable(test_phasor_load runLoadTests.cpp) target_link_libraries(test_phasor_load GRIDKIT::phasor_dynamics_load GRIDKIT::phasor_dynamics_bus) - + +add_executable(test_phasor_synchronous_machine runSynchronousMachineTests.cpp) +target_link_libraries(test_phasor_synchronous_machine GRIDKIT::phasor_dynamics_synchronous_machine + GRIDKIT::phasor_dynamics_bus) + + add_executable(test_phasor_system runSystemTests.cpp) target_link_libraries(test_phasor_system GRIDKIT::phasor_dynamics_load GRIDKIT::phasor_dynamics_branch @@ -22,7 +27,12 @@ target_link_libraries(test_phasor_system GRIDKIT::phasor_dynamics_load add_test(NAME PhasorDynamicsBusTest COMMAND $) add_test(NAME PhasorDynamicsBranchTest COMMAND $) +add_test(NAME PhasorDynamicsSynchronousMachineTest COMMAND $) add_test(NAME PhasorDynamicsLoadTest COMMAND $) add_test(NAME PhasorDynamicsSystemTest COMMAND $) -install(TARGETS test_phasor_bus test_phasor_branch test_phasor_load test_phasor_system RUNTIME DESTINATION bin) +install(TARGETS test_phasor_bus + test_phasor_branch + test_phasor_load + test_phasor_synchronous_machine + test_phasor_system RUNTIME DESTINATION bin) diff --git a/tests/UnitTests/PhasorDynamics/SynchronousMachineTests.hpp b/tests/UnitTests/PhasorDynamics/SynchronousMachineTests.hpp new file mode 100644 index 00000000..b266b182 --- /dev/null +++ b/tests/UnitTests/PhasorDynamics/SynchronousMachineTests.hpp @@ -0,0 +1,67 @@ +#include +#include + +#include +#include +#include +#include +#include + +namespace GridKit +{ +namespace Testing +{ + + template + class SynchronousMachineTests + { + private: + using real_type = typename PhasorDynamics::Component::real_type; + + public: + SynchronousMachineTests() = default; + ~SynchronousMachineTests() = default; + + TestOutcome constructor() + { + TestStatus success = true; + + auto* bus = new PhasorDynamics::Bus(1.0, 0.0); + + PhasorDynamics::Component* machine = + new PhasorDynamics::SynchronousMachine(bus); + + success *= (machine != nullptr); + + if (machine) + { + delete machine; + } + delete bus; + + return success.report(__func__); + } + + TestOutcome residual() + { + TestStatus success = true; + success.skipTest(); + + + return success.report(__func__); + } + + TestOutcome accessors() + { + TestStatus success = true; + success.skipTest(); + + + return success.report(__func__); + } + }; // class SynchronousMachineTest + +} // namespace Testing +} // namespace GridKit + + diff --git a/tests/UnitTests/PhasorDynamics/runSynchronousMachineTests.cpp b/tests/UnitTests/PhasorDynamics/runSynchronousMachineTests.cpp new file mode 100644 index 00000000..e2da3b4f --- /dev/null +++ b/tests/UnitTests/PhasorDynamics/runSynchronousMachineTests.cpp @@ -0,0 +1,16 @@ +#include "SynchronousMachineTests.hpp" + +int main() +{ + using namespace GridKit; + using namespace GridKit::Testing; + + GridKit::Testing::TestingResults result; + GridKit::Testing::SynchronousMachineTests test; + + result += test.constructor(); + result += test.accessors(); + result += test.residual(); + + return result.summary(); +} \ No newline at end of file