From 775ed9cde97bee97c926ee02054b1d48760d3e2b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20LAIGRE?= Date: Fri, 29 Apr 2022 14:40:33 +0200 Subject: [PATCH 1/3] Add DanglingLineData (and unit tests) and update DanglingLineBoundary implementation to use it MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: Sébastien LAIGRE --- .../iidm/util/DanglingLineBoundary.hpp | 3 + .../powsybl/iidm/util/DanglingLineData.hpp | 63 ++++++++++++ include/powsybl/math/ComplexUtils.hpp | 27 +++++ src/CMakeLists.txt | 2 + src/iidm/util/DanglingLineBoundary.cpp | 23 +++++ src/iidm/util/DanglingLineData.cpp | 98 +++++++++++++++++++ src/math/ComplexUtils.cpp | 30 ++++++ test/iidm/CMakeLists.txt | 1 + test/iidm/DanglingLineTest.cpp | 17 +++- test/iidm/util/DanglingLineDataTest.cpp | 86 ++++++++++++++++ 10 files changed, 345 insertions(+), 5 deletions(-) create mode 100644 include/powsybl/iidm/util/DanglingLineData.hpp create mode 100644 include/powsybl/math/ComplexUtils.hpp create mode 100644 src/iidm/util/DanglingLineData.cpp create mode 100644 src/math/ComplexUtils.cpp create mode 100644 test/iidm/util/DanglingLineDataTest.cpp diff --git a/include/powsybl/iidm/util/DanglingLineBoundary.hpp b/include/powsybl/iidm/util/DanglingLineBoundary.hpp index 6061d009c..898cc88b4 100644 --- a/include/powsybl/iidm/util/DanglingLineBoundary.hpp +++ b/include/powsybl/iidm/util/DanglingLineBoundary.hpp @@ -59,6 +59,9 @@ class Boundary : public iidm::Boundary { Boundary& operator=(Boundary&&) noexcept = delete; +private: + static bool valid(double p0, double q0); + private: DanglingLine& m_parent; }; diff --git a/include/powsybl/iidm/util/DanglingLineData.hpp b/include/powsybl/iidm/util/DanglingLineData.hpp new file mode 100644 index 000000000..add4204f7 --- /dev/null +++ b/include/powsybl/iidm/util/DanglingLineData.hpp @@ -0,0 +1,63 @@ +/** + * Copyright (c) 2022, RTE (http://www.rte-france.com) + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at http://mozilla.org/MPL/2.0/. + */ + +#ifndef POWSYBL_IIDM_UTIL_DANGLINGLINEDATA_HPP +#define POWSYBL_IIDM_UTIL_DANGLINGLINEDATA_HPP + +#include + +#include + +namespace powsybl { + +namespace iidm { + +class DanglingLine; + +class DanglingLineData { +public: + explicit DanglingLineData(const DanglingLine& danglingLine); + + DanglingLineData(const DanglingLine& danglingLine, bool splitShuntAdmittance); + + double getBoundaryBusTheta() const; + + double getBoundaryBusU() const; + + const std::string& getId() const; + +private: + static double getTheta(const DanglingLine& danglingLine); + + static double getV(const DanglingLine& danglingLine); + + static bool valid(double v, double theta); + +private: + std::string m_id; + + double m_r; + double m_x; + double m_g1; + double m_g2; + double m_b1; + double m_b2; + + double m_u1; + double m_theta1; + double m_p0; + double m_q0; + + double m_boundaryBusU = stdcxx::nan(); + double m_boundaryBusTheta = stdcxx::nan(); +}; + +} // namespace iidm + +} // namespace powsybl + +#endif // POWSYBL_IIDM_UTIL_DANGLINGLINEDATA_HPP diff --git a/include/powsybl/math/ComplexUtils.hpp b/include/powsybl/math/ComplexUtils.hpp new file mode 100644 index 000000000..291ec1636 --- /dev/null +++ b/include/powsybl/math/ComplexUtils.hpp @@ -0,0 +1,27 @@ +/** + * Copyright (c) 2022, RTE (http://www.rte-france.com) + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at http://mozilla.org/MPL/2.0/. + */ + +#ifndef POWSYBL_MATH_COMPLEXUTILS_HPP +#define POWSYBL_MATH_COMPLEXUTILS_HPP + +#include + +namespace powsybl { + +namespace math { + +namespace ComplexUtils { + +std::complex polar2Complex(double r, double theta); + +} // namespace ComplexUtils + +} // namespace math + +} // namespace powsybl + +#endif // POWSYBL_MATH_COMPLEXUTILS_HPP diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 49d078892..cb9da53b8 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -200,6 +200,7 @@ set(IIDM_SOURCES iidm/util/AbstractHalfLineBoundary.cpp iidm/util/ConnectedComponents.cpp iidm/util/DanglingLineBoundary.cpp + iidm/util/DanglingLineData.cpp iidm/util/DistinctPredicate.cpp iidm/util/Identifiables.cpp iidm/util/LimitViolationUtils.cpp @@ -217,6 +218,7 @@ set(IIDM_SOURCES logging/LogMessage.cpp logging/NoopLogger.cpp + math/ComplexUtils.cpp math/ConnectedComponentsComputationResult.cpp math/GraphUtil.cpp diff --git a/src/iidm/util/DanglingLineBoundary.cpp b/src/iidm/util/DanglingLineBoundary.cpp index 8d97542c9..74ecac705 100644 --- a/src/iidm/util/DanglingLineBoundary.cpp +++ b/src/iidm/util/DanglingLineBoundary.cpp @@ -10,6 +10,7 @@ #include #include #include +#include #include #include @@ -26,6 +27,11 @@ Boundary::Boundary(DanglingLine& parent) : } double Boundary::getAngle() const { + if (valid(m_parent.getP0(), m_parent.getQ0())) { + DanglingLineData danglingLineData(m_parent, true); + return stdcxx::toDegrees * danglingLineData.getBoundaryBusTheta(); + } + const Terminal& t = m_parent.getTerminal(); const stdcxx::CReference& b = t.getBusView().getBus(); return SV(t.getP(), t.getQ(), iidm::Boundary::getV(b), iidm::Boundary::getAngle(b), Branch::Side::ONE).otherSideA(m_parent, true); @@ -40,12 +46,20 @@ Connectable& Boundary::getConnectable() { } double Boundary::getP() const { + if (valid(m_parent.getP0(), m_parent.getQ0())) { + return -m_parent.getP0(); + } + const Terminal& t = m_parent.getTerminal(); const auto& b = t.getBusView().getBus(); return SV(t.getP(), t.getQ(), iidm::Boundary::getV(b), iidm::Boundary::getAngle(b), Branch::Side::ONE).otherSideP(m_parent, true); } double Boundary::getQ() const { + if (valid(m_parent.getP0(), m_parent.getQ0())) { + return -m_parent.getQ0(); + } + const Terminal& t = m_parent.getTerminal(); const auto& b = t.getBusView().getBus(); return SV(t.getP(), t.getQ(), iidm::Boundary::getV(b), iidm::Boundary::getAngle(b), Branch::Side::ONE).otherSideQ(m_parent, true); @@ -56,6 +70,11 @@ stdcxx::optional Boundary::getSide() const { } double Boundary::getV() const { + if (valid(m_parent.getP0(), m_parent.getQ0())) { + DanglingLineData danglingLineData(m_parent, true); + return danglingLineData.getBoundaryBusU(); + } + const Terminal& t = m_parent.getTerminal(); const auto& b = t.getBusView().getBus(); return SV(t.getP(), t.getQ(), iidm::Boundary::getV(b), iidm::Boundary::getAngle(b), Branch::Side::ONE).otherSideU(m_parent, true); @@ -69,6 +88,10 @@ VoltageLevel& Boundary::getVoltageLevel() { return m_parent.getTerminal().getVoltageLevel(); } +bool Boundary::valid(double p0, double q0) { + return !std::isnan(p0) && !std::isnan(q0); +} + } // namespace dangling_line } // namespace util diff --git a/src/iidm/util/DanglingLineData.cpp b/src/iidm/util/DanglingLineData.cpp new file mode 100644 index 000000000..28dc52fd9 --- /dev/null +++ b/src/iidm/util/DanglingLineData.cpp @@ -0,0 +1,98 @@ +/** + * Copyright (c) 2022, RTE (http://www.rte-france.com) + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at http://mozilla.org/MPL/2.0/. + */ + +#include + +#include + +#include +#include +#include +#include + +namespace powsybl { + +namespace iidm { + +DanglingLineData::DanglingLineData(const DanglingLine& danglingLine) : + DanglingLineData(danglingLine, true) { +} + +DanglingLineData::DanglingLineData(const DanglingLine& danglingLine, bool splitShuntAdmittance) : + m_id(danglingLine.getId()), + m_r(danglingLine.getR()), + m_x(danglingLine.getX()), + m_g1(splitShuntAdmittance ? danglingLine.getG() * 0.5 : danglingLine.getG()), + m_g2(splitShuntAdmittance ? danglingLine.getG() * 0.5 : 0.0), + m_b1(splitShuntAdmittance ? danglingLine.getB() * 0.5 : danglingLine.getB()), + m_b2(splitShuntAdmittance ? danglingLine.getB() * 0.5 : 0.0), + m_u1(getV(danglingLine)), + m_theta1(getTheta(danglingLine)), + m_p0(danglingLine.getP0()), + m_q0(danglingLine.getQ0()) { + + if (!valid(m_u1, m_theta1)) { + return; + } + + const std::complex& v1 = math::ComplexUtils::polar2Complex(m_u1, m_theta1); + + std::complex vBoundaryBus(stdcxx::nan(), stdcxx::nan()); + if (m_p0 == 0.0 && m_q0 == 0.0) { + LinkData::BranchAdmittanceMatrix adm = LinkData::calculateBranchAdmittance(m_r, m_x, 1.0, 0.0, 1.0, 0.0, std::complex(m_g1, m_b1), std::complex(m_g2, m_b2)); + vBoundaryBus = (- adm.y21 * v1) / adm.y22; + } else { + // Two buses Loadflow + std::complex sBoundary(-m_p0, -m_q0); + std::complex ytr = 1.0 / std::complex(m_r, m_x); + std::complex ysh2(m_g2, m_b2); + std::complex zt = 1.0 / (ytr + ysh2); + std::complex v0 = ytr * v1 / (ytr + ysh2); + double v02 = std::abs(v0) * std::abs(v0); + + std::complex sigma = zt * std::conj(sBoundary) / v02; + double d = 0.25 + std::real(sigma) - std::imag(sigma) * std::imag(sigma); + // d < 0 Collapsed network + if (d >= 0) { + vBoundaryBus = std::complex(0.5 + std::sqrt(d), std::imag(sigma)) * v0; + } + } + + m_boundaryBusU = std::abs(vBoundaryBus); + m_boundaryBusTheta = std::arg(vBoundaryBus); +} + +double DanglingLineData::getBoundaryBusTheta() const { + return m_boundaryBusTheta; +} + +double DanglingLineData::getBoundaryBusU() const { + return m_boundaryBusU; +} + +const std::string& DanglingLineData::getId() const { + return m_id; +} + +double DanglingLineData::getTheta(const DanglingLine& danglingLine) { + return danglingLine.getTerminal().isConnected() ? stdcxx::toRadians * danglingLine.getTerminal().getBusView().getBus().get().getAngle() : stdcxx::nan(); +} + +double DanglingLineData::getV(const DanglingLine& danglingLine) { + return danglingLine.getTerminal().isConnected() ? danglingLine.getTerminal().getBusView().getBus().get().getV() : stdcxx::nan(); +} + +bool DanglingLineData::valid(double v, double theta) { + if (std::isnan(v) || v <= 0.0) { + return false; + } + return !std::isnan(theta); +} + +} // namespace iidm + +} // namespace powsybl diff --git a/src/math/ComplexUtils.cpp b/src/math/ComplexUtils.cpp new file mode 100644 index 000000000..bb686bf31 --- /dev/null +++ b/src/math/ComplexUtils.cpp @@ -0,0 +1,30 @@ +/** + * Copyright (c) 2022, RTE (http://www.rte-france.com) + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at http://mozilla.org/MPL/2.0/. + */ + +#include + +#include +#include + +namespace powsybl { + +namespace math { + +namespace ComplexUtils { + +std::complex polar2Complex(double r, double theta) { + if (r < 0.0) { + throw PowsyblException(stdcxx::format("negative complex module %1%", r)); + } + return std::complex(r * std::cos(theta), r * std::sin(theta)); +} + +} // namespace ComplexUtils + +} // namespace math + +} // namespace powsybl diff --git a/test/iidm/CMakeLists.txt b/test/iidm/CMakeLists.txt index 38774e8d2..150dc3d2a 100644 --- a/test/iidm/CMakeLists.txt +++ b/test/iidm/CMakeLists.txt @@ -58,6 +58,7 @@ set(UNIT_TEST_SOURCES extensions/LoadDetailTest.cpp extensions/SlackTerminalTest.cpp + util/DanglingLineDataTest.cpp util/SVTest.cpp util/TerminalFinderTest.cpp ) diff --git a/test/iidm/DanglingLineTest.cpp b/test/iidm/DanglingLineTest.cpp index bc6d639b2..3742aa47e 100644 --- a/test/iidm/DanglingLineTest.cpp +++ b/test/iidm/DanglingLineTest.cpp @@ -445,19 +445,26 @@ BOOST_AUTO_TEST_CASE(getBoundary) { DanglingLine& danglingLine = network.getDanglingLine("DL1"); const DanglingLine& cDanglingLine = network.getDanglingLine("DL1"); + danglingLine.setB(0.0); + danglingLine.setG(2.0); + danglingLine.setP0(-200.0); + danglingLine.setQ0(0); + danglingLine.setR(5.0); + danglingLine.setX(0.0); + danglingLine.getTerminal().getBusView().getBus().get().setAngle(2); danglingLine.getTerminal().setP(3); - danglingLine.getTerminal().setQ(4); + danglingLine.getTerminal().setQ(0); danglingLine.getTerminal().getBusView().getBus().get().setV(5); BOOST_CHECK(stdcxx::areSame(cDanglingLine.getBoundary(), danglingLine.getBoundary())); const Boundary& cBoundary = danglingLine.getBoundary(); Boundary& boundary = danglingLine.getBoundary(); constexpr double ACCEPTABLE_THRESHOLD = 1e-6; - BOOST_CHECK_CLOSE(82.47271661854765, boundary.getAngle(), ACCEPTABLE_THRESHOLD); - BOOST_CHECK_CLOSE(2065.500000000001, boundary.getP(), ACCEPTABLE_THRESHOLD); - BOOST_CHECK_CLOSE(-781.1250000000001, boundary.getQ(), ACCEPTABLE_THRESHOLD); - BOOST_CHECK_CLOSE(43.5, boundary.getV(), ACCEPTABLE_THRESHOLD); + BOOST_CHECK_CLOSE(2.0, boundary.getAngle(), ACCEPTABLE_THRESHOLD); + BOOST_CHECK_CLOSE(200, boundary.getP(), ACCEPTABLE_THRESHOLD); + BOOST_CHECK_CLOSE(0, boundary.getQ(), ACCEPTABLE_THRESHOLD); + BOOST_CHECK_CLOSE(13.333333333, boundary.getV(), ACCEPTABLE_THRESHOLD); BOOST_CHECK(stdcxx::areSame(cDanglingLine, cBoundary.getConnectable())); BOOST_CHECK(stdcxx::areSame(cDanglingLine, boundary.getConnectable())); BOOST_CHECK(!boundary.getSide()); diff --git a/test/iidm/util/DanglingLineDataTest.cpp b/test/iidm/util/DanglingLineDataTest.cpp new file mode 100644 index 000000000..e165ab661 --- /dev/null +++ b/test/iidm/util/DanglingLineDataTest.cpp @@ -0,0 +1,86 @@ +/** + * Copyright (c) 2022, RTE (http://www.rte-france.com) + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at http://mozilla.org/MPL/2.0/. + */ + +#include + +#include +#include +#include +#include +#include +#include +#include +#include + +namespace powsybl { + +namespace iidm { + +Network createDanglingDataNetwork() { + Network network("test", "test"); + network.setCaseDate(stdcxx::DateTime::parse("2016-06-27T12:27:58.535+02:00")); + Substation& s = network.newSubstation() + .setId("S") + .setCountry(Country::FR) + .add(); + VoltageLevel& vl = s.newVoltageLevel() + .setId("VL") + .setNominalV(400) + .setTopologyKind(TopologyKind::BUS_BREAKER) + .add(); + vl.getBusBreakerView().newBus() + .setId("B") + .add(); + vl.getBusBreakerView().getBus("B").get().setAngle(-8.60); + vl.getBusBreakerView().getBus("B").get().setV(406.62); + vl.newDanglingLine() + .setId("DL") + .setBus("B") + .setR(0.05) + .setX(0.2) + .setG(0.0) + .setB(0.000001) + .setP0(-367.40) + .setQ0(63.73) + .add(); + + return network; +} + +bool dlCompareBoundaryBusVoltage(const DanglingLineData& dlData, double boundaryBusU, double boundaryBusAngle) { + double tol = 0.00001; + if (std::abs(dlData.getBoundaryBusU() - boundaryBusU) > tol || std::abs(stdcxx::toDegrees * dlData.getBoundaryBusTheta() - boundaryBusAngle) > tol) { + return false; + } + return true; +} + +BOOST_AUTO_TEST_SUITE(DanglingLineDataTestSuite) + +BOOST_AUTO_TEST_CASE(test) { + Network network = createDanglingDataNetwork(); + const DanglingLine& danglingLine = network.getDanglingLine("DL"); + DanglingLineData dlData(danglingLine); + + BOOST_CHECK(dlCompareBoundaryBusVoltage(dlData, 406.63382758266334, -8.573434828294932)); +} + +BOOST_AUTO_TEST_CASE(testP0Q0zero) { + Network network = createDanglingDataNetwork(); + DanglingLine& danglingLine = network.getDanglingLine("DL"); + danglingLine.setP0(0.0); + danglingLine.setQ0(0.0); + DanglingLineData dlData(danglingLine); + + BOOST_CHECK(dlCompareBoundaryBusVoltage(dlData, 406.6200406620039, -8.60000143239463)); +} + +BOOST_AUTO_TEST_SUITE_END() + +} // namespace iidm + +} // namespace powsybl From feadf14b49b97103a7c60f2a1a0654a7d6bffd17 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20LAIGRE?= Date: Fri, 29 Apr 2022 18:54:03 +0200 Subject: [PATCH 2/3] Add Twtdata and unit tests MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: Sébastien LAIGRE --- include/powsybl/iidm/util/LinkData.hpp | 28 +- include/powsybl/iidm/util/TwtData.hpp | 210 ++++++++ src/CMakeLists.txt | 1 + src/iidm/util/LinkData.cpp | 113 +++++ src/iidm/util/TwtData.cpp | 610 ++++++++++++++++++++++++ test/iidm/CMakeLists.txt | 1 + test/iidm/util/DanglingLineDataTest.cpp | 5 +- test/iidm/util/TwtDataTest.cpp | 496 +++++++++++++++++++ 8 files changed, 1458 insertions(+), 6 deletions(-) create mode 100644 include/powsybl/iidm/util/TwtData.hpp create mode 100644 src/iidm/util/TwtData.cpp create mode 100644 test/iidm/util/TwtDataTest.cpp diff --git a/include/powsybl/iidm/util/LinkData.hpp b/include/powsybl/iidm/util/LinkData.hpp index 400bf35bc..df23aef05 100644 --- a/include/powsybl/iidm/util/LinkData.hpp +++ b/include/powsybl/iidm/util/LinkData.hpp @@ -10,6 +10,8 @@ #include +#include + namespace powsybl { namespace iidm { @@ -25,6 +27,10 @@ struct Flow { struct BranchAdmittanceMatrix { public: + BranchAdmittanceMatrix() = default; + + BranchAdmittanceMatrix(const std::complex& y11, const std::complex& y12, const std::complex& y21, const std::complex& y22); + std::complex y11; std::complex y12; @@ -34,8 +40,26 @@ struct BranchAdmittanceMatrix { std::complex y22; }; -BranchAdmittanceMatrix calculateBranchAdmittance(double r, double x, double ratio1, double angle1, - double ratio2, double angle2, const std::complex& ysh1, const std::complex& ysh2); +BranchAdmittanceMatrix calculateBranchAdmittance(double r, double x, double ratio1, double angle1, double ratio2, double angle2, const std::complex& ysh1, + const std::complex& ysh2); + +Flow flowBothEnds(const std::complex& y11, const std::complex& y12, const std::complex& y21, const std::complex& y22, double u1, + double theta1, double u2, double theta2); + +Flow flowBothEnds(const std::complex& y11, const std::complex& y12, const std::complex& y21, const std::complex& y22, + const std::complex& v1, const std::complex& v2); + +std::complex flowYshunt(const std::complex& ysh, double u, double theta); + +double getFixedX(double x, double epsilonX, bool applyReactanceCorrection); + +double getPhaseAngleClockDegrees(int phaseAngleClock); + +std::complex kronAntenna(const std::complex& y11, const std::complex& y12, const std::complex& y21, const std::complex& y22, + bool isOpenFrom); + +BranchAdmittanceMatrix kronChain(const BranchAdmittanceMatrix& firstAdm, const Branch::Side& firstChainNodeSide, const BranchAdmittanceMatrix& secondAdm, + const Branch::Side& secondChainNodeSide); } // namespace LinkData diff --git a/include/powsybl/iidm/util/TwtData.hpp b/include/powsybl/iidm/util/TwtData.hpp new file mode 100644 index 000000000..c85bb9fa9 --- /dev/null +++ b/include/powsybl/iidm/util/TwtData.hpp @@ -0,0 +1,210 @@ +/** + * Copyright (c) 2022, RTE (http://www.rte-france.com) + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at http://mozilla.org/MPL/2.0/. + */ + +#ifndef POWSYBL_IIDM_UTIL_TWTDATA_HPP +#define POWSYBL_IIDM_UTIL_TWTDATA_HPP + +#include +#include + +#include +#include + +namespace powsybl { + +namespace iidm { + +class ThreeWindingsTransformer; + +class TwtData { +public: + TwtData(const ThreeWindingsTransformer& twt, double epsilonX, bool applyReactanceCorrection); + + TwtData(const ThreeWindingsTransformer& twt, double epsilonX, bool applyReactanceCorrection, bool twtSplitShuntAdmittance); + + TwtData(const ThreeWindingsTransformer& twt, int phaseAngleClock2, int phaseAngleClock3, double epsilonX, bool applyReactanceCorrection, + bool twtSplitShuntAdmittance); + + double getB1(const ThreeWindingsTransformer::Side& side) const; + + double getB2(const ThreeWindingsTransformer::Side& side) const; + + double getComputedP(const ThreeWindingsTransformer::Side& side) const; + + double getComputedQ(const ThreeWindingsTransformer::Side& side) const; + + double getG1(const ThreeWindingsTransformer::Side& side) const; + + double getG2(const ThreeWindingsTransformer::Side& side) const; + + const std::string& getId() const; + + double getP(const ThreeWindingsTransformer::Side& side) const; + + double getQ(const ThreeWindingsTransformer::Side& side) const; + + double getR(const ThreeWindingsTransformer::Side& side) const; + + double getRatedU(const ThreeWindingsTransformer::Side& side) const; + + double getStarTheta() const; + + double getStarU() const; + + double getTheta(const ThreeWindingsTransformer::Side& side) const; + + double getU(const ThreeWindingsTransformer::Side& side) const; + + double getX(const ThreeWindingsTransformer::Side& side) const; + + bool isConnected(const ThreeWindingsTransformer::Side& side) const; + + bool isMainComponent(const ThreeWindingsTransformer::Side& side) const; + + int getPhaseAngleClock2() const; + + int getPhaseAngleClock3() const; + + double getRatedU0() const; + +private: + static double alpha(const ThreeWindingsTransformer::Leg& leg); + + static double getB1(const ThreeWindingsTransformer::Leg& leg, bool twtSplitShuntAdmittance); + + static double getB2(const ThreeWindingsTransformer::Leg& leg, bool twtSplitShuntAdmittance); + + static double getG1(const ThreeWindingsTransformer::Leg& leg, bool twtSplitShuntAdmittance); + + static double getG2(const ThreeWindingsTransformer::Leg& leg, bool twtSplitShuntAdmittance); + + static double getR(const ThreeWindingsTransformer::Leg& leg); + + static double getTheta(const ThreeWindingsTransformer::Leg& leg); + + static double getV(const ThreeWindingsTransformer::Leg& leg); + + static double getValue(double initialValue, double rtcStepValue, double ptcStepValue); + + static double getX(const ThreeWindingsTransformer::Leg& leg); + + static bool isMainComponent(const ThreeWindingsTransformer::Leg& leg); + + static double rho(const ThreeWindingsTransformer::Leg& leg, double ratedU0); + + static bool valid(double voltage, double theta); + +private: + std::complex calculateOneConnectedLegFlow(double u, double theta, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixFirstOpenLeg, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixSecondOpenLeg) const; + + std::complex calculateOneConnectedLegShunt(const LinkData::BranchAdmittanceMatrix& closeLeg, + const LinkData::BranchAdmittanceMatrix& firstOpenLeg, + const LinkData::BranchAdmittanceMatrix& secondOpenLeg) const; + + std::complex calculateOneConnectedLegStarBusVoltage(double u, double theta, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixFirstOpenLeg, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixSecondOpenLeg) const; + + void calculateThreeConnectedLegsFlowAndStarBusVoltage(double u1, double theta1, double u2, double theta2, + double u3, double theta3, + const LinkData::BranchAdmittanceMatrix& branchAdmittanceLeg1, + const LinkData::BranchAdmittanceMatrix& branchAdmittanceLeg2, + const LinkData::BranchAdmittanceMatrix& branchAdmittanceLeg3); + + LinkData::BranchAdmittanceMatrix calculateTwoConnectedLegsAdmittance(const LinkData::BranchAdmittanceMatrix& firstCloseLeg, + const LinkData::BranchAdmittanceMatrix& secondCloseLeg, + const LinkData::BranchAdmittanceMatrix& openLeg) const; + + LinkData::Flow calculateTwoConnectedLegsFlow(double u1, double theta1, double u2, double theta2, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg1, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg2, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixOpenLeg) const; + + std::complex calculateTwoConnectedLegsStarBusVoltage(double u1, double theta1, double u2, double theta2, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg1, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg2, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixOpenLeg) const; + +private: + std::string m_id; + + double m_p1; + double m_q1; + double m_p2; + double m_q2; + double m_p3; + double m_q3; + + double m_u1; + double m_theta1; + double m_u2; + double m_theta2; + double m_u3; + double m_theta3; + + double m_r1; + double m_x1; + double m_r2; + double m_x2; + double m_r3; + double m_x3; + + double m_g11; + double m_b11; + double m_g12; + double m_b12; + double m_g21; + double m_b21; + double m_g22; + double m_b22; + double m_g31; + double m_b31; + double m_g32; + double m_b32; + + double m_rho1; + double m_alpha1; + double m_rho2; + double m_alpha2; + double m_rho3; + double m_alpha3; + + double m_ratedU1; + double m_ratedU2; + double m_ratedU3; + + bool m_connected1; + bool m_connected2; + bool m_connected3; + bool m_mainComponent1; + bool m_mainComponent2; + bool m_mainComponent3; + + double m_computedP1; + double m_computedQ1; + double m_computedP2; + double m_computedQ2; + double m_computedP3; + double m_computedQ3; + + double m_starU; + double m_starTheta; + + int m_phaseAngleClock2; + int m_phaseAngleClock3; + double m_ratedU0; +}; + +} // namespace iidm + +} // namespace powsybl + +#endif // POWSYBL_IIDM_UTIL_TWTDATA_HPP diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index cb9da53b8..6e93ced1a 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -209,6 +209,7 @@ set(IIDM_SOURCES iidm/util/Substations.cpp iidm/util/SV.cpp iidm/util/TerminalFinder.cpp + iidm/util/TwtData.cpp iidm/util/VoltageLevels.cpp logging/ConsoleLogger.cpp diff --git a/src/iidm/util/LinkData.cpp b/src/iidm/util/LinkData.cpp index 0a7478d3a..0ac9a2d35 100644 --- a/src/iidm/util/LinkData.cpp +++ b/src/iidm/util/LinkData.cpp @@ -7,12 +7,22 @@ #include +#include + namespace powsybl { namespace iidm { namespace LinkData { +BranchAdmittanceMatrix::BranchAdmittanceMatrix(const std::complex& y11, const std::complex& y12, const std::complex& y21, + const std::complex& y22) : + y11(y11), + y12(y12), + y21(y21), + y22(y22) { +} + BranchAdmittanceMatrix calculateBranchAdmittance(double r, double x, double ratio1, double angle1, double ratio2, double angle2, const std::complex& ysh1, const std::complex& ysh2) { std::complex a1 = std::polar(ratio1, angle1); @@ -33,6 +43,109 @@ BranchAdmittanceMatrix calculateBranchAdmittance(double r, double x, double rati return branchAdmittance; } +Flow flowBothEnds(const std::complex& y11, const std::complex& y12, const std::complex& y21, const std::complex& y22, double u1, + double theta1, double u2, double theta2) { + + const std::complex& v1 = math::ComplexUtils::polar2Complex(u1, theta1); + const std::complex& v2 = math::ComplexUtils::polar2Complex(u2, theta2); + + return flowBothEnds(y11, y12, y21, y22, v1, v2); +} + +Flow flowBothEnds(const std::complex& y11, const std::complex& y12, const std::complex& y21, const std::complex& y22, + const std::complex& v1, const std::complex& v2) { + + Flow flow; + const std::complex& ift = y12 * v2 + y11 * v1; + flow.fromTo = std::conj(ift) * v1; + + const std::complex& itf = y21 * v1 + y22 * v2; + flow.toFrom = std::conj(itf) * v2; + + return flow; +} + +std::complex flowYshunt(const std::complex& ysh, double u, double theta) { + std::complex v = math::ComplexUtils::polar2Complex(u, theta); + + return std::conj(ysh) * std::conj(v) *v; +} + +double getFixedX(double x, double epsilonX, bool applyReactanceCorrection) { + return std::abs(x) < epsilonX && applyReactanceCorrection ? epsilonX : x; +} + +double getPhaseAngleClockDegrees(int phaseAngleClock) { + double phaseAngleClockDegree = 0.0; + phaseAngleClockDegree += phaseAngleClock * 30.0; + phaseAngleClockDegree = std::remainder(phaseAngleClockDegree, 360.0); + if (phaseAngleClockDegree > 180.0) { + phaseAngleClockDegree -= 360.0; + } + return phaseAngleClockDegree; +} + +std::complex kronAntenna(const std::complex& y11, const std::complex& y12, const std::complex& y21, const std::complex& y22, + bool isOpenFrom) { + std::complex ysh(0.0, 0.0); + + if (isOpenFrom) { + if (y11 != std::complex(0.0, 0.0)) { + ysh = y22 - (y21 * y12 / y11); + } + } + else { + if (y22 != std::complex(0.0, 0.0)) { + ysh = y11 - (y12 * y21 / y22); + } + } + return ysh; +} + +BranchAdmittanceMatrix kronChain(const BranchAdmittanceMatrix& firstAdm, const Branch::Side& firstChainNodeSide, const BranchAdmittanceMatrix& secondAdm, + const Branch::Side& secondChainNodeSide) { + BranchAdmittanceMatrix admittance; + + std::complex yFirst11; + std::complex yFirst1C; + std::complex yFirstC1; + std::complex yFirstCC; + if (firstChainNodeSide == Branch::Side::TWO) { + yFirst11 = firstAdm.y11; + yFirst1C = firstAdm.y12; + yFirstC1 = firstAdm.y21; + yFirstCC = firstAdm.y22; + } else { + yFirst11 = firstAdm.y22; + yFirst1C = firstAdm.y21; + yFirstC1 = firstAdm.y12; + yFirstCC = firstAdm.y11; + } + + std::complex ySecond22; + std::complex ySecond2C; + std::complex ySecondC2; + std::complex ySecondCC; + if (secondChainNodeSide == Branch::Side::TWO) { + ySecond22 = secondAdm.y11; + ySecond2C = secondAdm.y12; + ySecondC2 = secondAdm.y21; + ySecondCC = secondAdm.y22; + } else { + ySecond22 = secondAdm.y22; + ySecond2C = secondAdm.y21; + ySecondC2 = secondAdm.y12; + ySecondCC = secondAdm.y11; + } + + admittance.y11 = yFirst11 - (yFirst1C * yFirstC1 / (yFirstCC + ySecondCC)); + admittance.y12 = yFirst1C * ySecondC2 / -(yFirstCC + ySecondCC); + admittance.y21 = ySecond2C * yFirstC1 / -(yFirstCC + ySecondCC); + admittance.y22 = ySecond22 - (ySecond2C * ySecondC2 / (yFirstCC + ySecondCC)); + + return admittance; +} + } // namespace LinkData } // namespace iidm diff --git a/src/iidm/util/TwtData.cpp b/src/iidm/util/TwtData.cpp new file mode 100644 index 000000000..5e1254731 --- /dev/null +++ b/src/iidm/util/TwtData.cpp @@ -0,0 +1,610 @@ +/** + * Copyright (c) 2022, RTE (http://www.rte-france.com) + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at http://mozilla.org/MPL/2.0/. + */ + +#include + +#include +#include +#include +#include +#include +#include +#include + +namespace powsybl { + +namespace iidm { + +const std::string UNEXPECTED_SIDE = "Unexpected side"; + +TwtData::TwtData(const ThreeWindingsTransformer& twt, double epsilonX, bool applyReactanceCorrection) : + TwtData(twt, 0, 0, epsilonX, applyReactanceCorrection, false) { +} + +TwtData::TwtData(const ThreeWindingsTransformer& twt, double epsilonX, bool applyReactanceCorrection, bool twtSplitShuntAdmittance) : + TwtData(twt, 0, 0, epsilonX, applyReactanceCorrection, twtSplitShuntAdmittance) { +} + +TwtData::TwtData(const ThreeWindingsTransformer& twt, int phaseAngleClock2, int phaseAngleClock3, double epsilonX, bool applyReactanceCorrection, + bool twtSplitShuntAdmittance) : + m_id(twt.getId()), + m_p1(twt.getLeg1().getTerminal().getP()), + m_q1(twt.getLeg1().getTerminal().getQ()), + m_p2(twt.getLeg2().getTerminal().getP()), + m_q2(twt.getLeg2().getTerminal().getQ()), + m_p3(twt.getLeg3().getTerminal().getP()), + m_q3(twt.getLeg3().getTerminal().getQ()), + m_u1(getV(twt.getLeg1())), + m_theta1(getTheta(twt.getLeg1())), + m_u2(getV(twt.getLeg2())), + m_theta2(getTheta(twt.getLeg2())), + m_u3(getV(twt.getLeg3())), + m_theta3(getTheta(twt.getLeg3())), + m_r1(getR(twt.getLeg1())), + m_x1(LinkData::getFixedX(getX(twt.getLeg1()), epsilonX, applyReactanceCorrection)), + m_r2(getR(twt.getLeg2())), + m_x2(LinkData::getFixedX(getX(twt.getLeg2()), epsilonX, applyReactanceCorrection)), + m_r3(getR(twt.getLeg3())), + m_x3(LinkData::getFixedX(getX(twt.getLeg3()), epsilonX, applyReactanceCorrection)), + m_g11(getG1(twt.getLeg1(), twtSplitShuntAdmittance)), + m_b11(getB1(twt.getLeg1(), twtSplitShuntAdmittance)), + m_g12(getG2(twt.getLeg1(), twtSplitShuntAdmittance)), + m_b12(getB2(twt.getLeg1(), twtSplitShuntAdmittance)), + m_g21(getG1(twt.getLeg2(), twtSplitShuntAdmittance)), + m_b21(getB1(twt.getLeg2(), twtSplitShuntAdmittance)), + m_g22(getG2(twt.getLeg2(), twtSplitShuntAdmittance)), + m_b22(getB2(twt.getLeg2(), twtSplitShuntAdmittance)), + m_g31(getG1(twt.getLeg3(), twtSplitShuntAdmittance)), + m_b31(getB1(twt.getLeg3(), twtSplitShuntAdmittance)), + m_g32(getG2(twt.getLeg3(), twtSplitShuntAdmittance)), + m_b32(getB2(twt.getLeg3(), twtSplitShuntAdmittance)), + m_alpha1(alpha(twt.getLeg1())), + m_alpha2(alpha(twt.getLeg2())), + m_alpha3(alpha(twt.getLeg3())), + m_ratedU1(twt.getLeg1().getRatedU()), + m_ratedU2(twt.getLeg2().getRatedU()), + m_ratedU3(twt.getLeg3().getRatedU()), + m_connected1(twt.getLeg1().getTerminal().isConnected()), + m_connected2(twt.getLeg2().getTerminal().isConnected()), + m_connected3(twt.getLeg3().getTerminal().isConnected()), + m_mainComponent1(isMainComponent(twt.getLeg1())), + m_mainComponent2(isMainComponent(twt.getLeg2())), + m_mainComponent3(isMainComponent(twt.getLeg3())), + m_phaseAngleClock2(phaseAngleClock2), + m_phaseAngleClock3(phaseAngleClock3), + m_ratedU0(twt.getRatedU0()) { + m_rho1 = rho(twt.getLeg1(), m_ratedU0); + m_rho2 = rho(twt.getLeg2(), m_ratedU0); + m_rho3 = rho(twt.getLeg3(), m_ratedU0); + + double rhof = 1.0; + double alphaf = 0.0; + + double angle1 = -m_alpha1; + double angle2 = -m_alpha2 - stdcxx::toRadians * LinkData::getPhaseAngleClockDegrees(phaseAngleClock2); + double angle3 = -m_alpha3 - stdcxx::toRadians * LinkData::getPhaseAngleClockDegrees(phaseAngleClock3); + double anglef = -alphaf; + + LinkData::BranchAdmittanceMatrix branchAdmittanceLeg1 = LinkData::calculateBranchAdmittance(m_r1, m_x1, 1 / m_rho1, angle1, 1 / rhof, anglef, std::complex(m_g11, m_b11), std::complex(m_g12, m_b12)); + + LinkData::BranchAdmittanceMatrix branchAdmittanceLeg2 = LinkData::calculateBranchAdmittance(m_r2, m_x2, 1 / m_rho2, angle2, 1 / rhof, anglef, std::complex(m_g21, m_b21), std::complex(m_g22, m_b22)); + + LinkData::BranchAdmittanceMatrix branchAdmittanceLeg3 = LinkData::calculateBranchAdmittance(m_r3, m_x3, 1 / m_rho3, angle3, 1 / rhof, anglef, std::complex(m_g31, m_b31), std::complex(m_g32, m_b32)); + + // Assume the ratedU at the star bus is equal to ratedU of Leg1 + + if (m_connected1 && m_connected2 && m_connected3 && valid(m_u1, m_theta1) && valid(m_u2, m_theta2) && valid(m_u3, m_theta3)) { + + calculateThreeConnectedLegsFlowAndStarBusVoltage(m_u1, m_theta1, m_u2, m_theta2, m_u3, m_theta3, branchAdmittanceLeg1, + branchAdmittanceLeg2, branchAdmittanceLeg3); + } else if (m_connected1 && m_connected2 && valid(m_u1, m_theta1) && valid(m_u2, m_theta2)) { + + LinkData::Flow flow = calculateTwoConnectedLegsFlow(m_u1, m_theta1, m_u2, m_theta2, + branchAdmittanceLeg1, branchAdmittanceLeg2, branchAdmittanceLeg3); + m_computedP1 = std::real(flow.fromTo); + m_computedQ1 = std::imag(flow.fromTo); + m_computedP2 = std::real(flow.toFrom); + m_computedQ2 = std::imag(flow.toFrom); + m_computedP3 = 0.0; + m_computedQ3 = 0.0; + + std::complex v0 = calculateTwoConnectedLegsStarBusVoltage(m_u1, m_theta1, m_u2, m_theta2, + branchAdmittanceLeg1, branchAdmittanceLeg2, branchAdmittanceLeg3); + m_starU = std::abs(v0); + m_starTheta = std::arg(v0); + } else if (m_connected1 && m_connected3 && valid(m_u1, m_theta1) && valid(m_u3, m_theta3)) { + + LinkData::Flow flow = calculateTwoConnectedLegsFlow(m_u1, m_theta1, m_u3, m_theta3, + branchAdmittanceLeg1, branchAdmittanceLeg3, branchAdmittanceLeg2); + m_computedP1 = std::real(flow.fromTo); + m_computedQ1 = std::imag(flow.fromTo); + m_computedP2 = 0.0; + m_computedQ2 = 0.0; + m_computedP3 = std::real(flow.toFrom); + m_computedQ3 = std::imag(flow.toFrom); + + std::complex v0 = calculateTwoConnectedLegsStarBusVoltage(m_u1, m_theta1, m_u3, m_theta3, + branchAdmittanceLeg1, branchAdmittanceLeg3, branchAdmittanceLeg2); + + m_starU = std::abs(v0); + m_starTheta = std::arg(v0); + } else if (m_connected2 && m_connected3 && valid(m_u2, m_theta2) && valid(m_u3, m_theta3)) { + + LinkData::Flow flow = calculateTwoConnectedLegsFlow(m_u2, m_theta2, m_u3, m_theta3, + branchAdmittanceLeg2, branchAdmittanceLeg3, branchAdmittanceLeg1); + m_computedP1 = 0.0; + m_computedQ1 = 0.0; + m_computedP2 = std::real(flow.fromTo); + m_computedQ2 = std::imag(flow.fromTo); + m_computedP3 = std::real(flow.toFrom); + m_computedQ3 = std::imag(flow.toFrom); + + std::complex v0 = calculateTwoConnectedLegsStarBusVoltage(m_u2, m_theta2, m_u3, m_theta3, + branchAdmittanceLeg2, branchAdmittanceLeg3, branchAdmittanceLeg1); + m_starU = std::abs(v0); + m_starTheta = std::arg(v0); + } else if (m_connected1 && valid(m_u1, m_theta1)) { + + std::complex flow = calculateOneConnectedLegFlow(m_u1, m_theta1, branchAdmittanceLeg1, + branchAdmittanceLeg2, branchAdmittanceLeg3); + m_computedP1 = std::real(flow); + m_computedQ1 = std::imag(flow); + m_computedP2 = 0.0; + m_computedQ2 = 0.0; + m_computedP3 = 0.0; + m_computedQ3 = 0.0; + + std::complex v0 = calculateOneConnectedLegStarBusVoltage(m_u1, m_theta1, branchAdmittanceLeg1, + branchAdmittanceLeg2, branchAdmittanceLeg3); + + m_starU = std::abs(v0); + m_starTheta = std::arg(v0); + } else if (m_connected2 && valid(m_u2, m_theta2)) { + + std::complex flow = calculateOneConnectedLegFlow(m_u2, m_theta2, branchAdmittanceLeg2, + branchAdmittanceLeg1, branchAdmittanceLeg3); + + m_computedP1 = 0.0; + m_computedQ1 = 0.0; + m_computedP2 = std::real(flow); + m_computedQ2 = std::imag(flow); + m_computedP3 = 0.0; + m_computedQ3 = 0.0; + + std::complex v0 = calculateOneConnectedLegStarBusVoltage(m_u2, m_theta2, branchAdmittanceLeg2, + branchAdmittanceLeg1, branchAdmittanceLeg3); + m_starU = std::abs(v0); + m_starTheta = std::arg(v0); + } else if (m_connected3 && valid(m_u3, m_theta3)) { + + std::complex flow = calculateOneConnectedLegFlow(m_u3, m_theta3, branchAdmittanceLeg3, + branchAdmittanceLeg1, branchAdmittanceLeg2); + + m_computedP1 = 0.0; + m_computedQ1 = 0.0; + m_computedP2 = 0.0; + m_computedQ2 = 0.0; + m_computedP3 = std::real(flow); + m_computedQ3 = std::imag(flow); + + std::complex v0 = calculateOneConnectedLegStarBusVoltage(m_u3, m_theta3, branchAdmittanceLeg3, + branchAdmittanceLeg1, branchAdmittanceLeg2); + m_starU = std::abs(v0); + m_starTheta = std::arg(v0); + } else { + + m_computedP1 = stdcxx::nan(); + m_computedQ1 = stdcxx::nan(); + m_computedP2 = stdcxx::nan(); + m_computedQ2 = stdcxx::nan(); + m_computedP3 = stdcxx::nan(); + m_computedQ3 = stdcxx::nan(); + + m_starU = stdcxx::nan(); + m_starTheta = stdcxx::nan(); + } +} + +double TwtData::alpha(const ThreeWindingsTransformer::Leg& leg) { + return leg.getOptionalPhaseTapChanger() ? leg.getPhaseTapChanger().getCurrentStep().getAlpha() : 0.0; +} + +std::complex TwtData::calculateOneConnectedLegFlow(double u, double theta, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixFirstOpenLeg, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixSecondOpenLeg) const { + + std::complex ysh = calculateOneConnectedLegShunt(admittanceMatrixLeg, + admittanceMatrixFirstOpenLeg, admittanceMatrixSecondOpenLeg); + + return LinkData::flowYshunt(ysh, u, theta); +} + +std::complex TwtData::calculateOneConnectedLegShunt(const LinkData::BranchAdmittanceMatrix& closeLeg, + const LinkData::BranchAdmittanceMatrix& firstOpenLeg, + const LinkData::BranchAdmittanceMatrix& secondOpenLeg) const { + std::complex ysh1 = LinkData::kronAntenna(firstOpenLeg.y11, firstOpenLeg.y12, firstOpenLeg.y21, firstOpenLeg.y22, true); + std::complex ysh2 = LinkData::kronAntenna(secondOpenLeg.y11, secondOpenLeg.y12, secondOpenLeg.y21, secondOpenLeg.y22, true); + std::complex y22 = closeLeg.y22 + ysh1 + ysh2; + + return LinkData::kronAntenna(closeLeg.y11, closeLeg.y12, closeLeg.y21, y22, false); +} + +std::complex TwtData::calculateOneConnectedLegStarBusVoltage(double u, double theta, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixFirstOpenLeg, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixSecondOpenLeg) const { + + std::complex v = math::ComplexUtils::polar2Complex(u, theta); + + std::complex ysh1O = LinkData::kronAntenna(admittanceMatrixFirstOpenLeg.y11, admittanceMatrixFirstOpenLeg.y12, + admittanceMatrixFirstOpenLeg.y21, admittanceMatrixFirstOpenLeg.y22, true); + std::complex ysh2O = LinkData::kronAntenna(admittanceMatrixSecondOpenLeg.y11, admittanceMatrixSecondOpenLeg.y12, + admittanceMatrixSecondOpenLeg.y21, admittanceMatrixSecondOpenLeg.y22, true); + + return (-admittanceMatrixLeg.y21 * v) / (admittanceMatrixLeg.y22 + ysh1O + ysh2O); +} + +void TwtData::calculateThreeConnectedLegsFlowAndStarBusVoltage(double u1, double theta1, double u2, double theta2, + double u3, double theta3, + const LinkData::BranchAdmittanceMatrix& branchAdmittanceLeg1, + const LinkData::BranchAdmittanceMatrix& branchAdmittanceLeg2, + const LinkData::BranchAdmittanceMatrix& branchAdmittanceLeg3) { + std::complex v1 = math::ComplexUtils::polar2Complex(u1, theta1); + std::complex v2 = math::ComplexUtils::polar2Complex(u2, theta2); + std::complex v3 = math::ComplexUtils::polar2Complex(u3, theta3); + + std::complex v0 = -(branchAdmittanceLeg1.y21 * v1 + branchAdmittanceLeg2.y21 * v2 + branchAdmittanceLeg3.y21 * v3) + / (branchAdmittanceLeg1.y22 + branchAdmittanceLeg2.y22 + branchAdmittanceLeg3.y22); + + LinkData::Flow flowLeg1 = LinkData::flowBothEnds(branchAdmittanceLeg1.y11, branchAdmittanceLeg1.y12, + branchAdmittanceLeg1.y21, branchAdmittanceLeg1.y22, v1, v0); + + LinkData::Flow flowLeg2 = LinkData::flowBothEnds(branchAdmittanceLeg2.y11, branchAdmittanceLeg2.y12, + branchAdmittanceLeg2.y21, branchAdmittanceLeg2.y22, v2, v0); + + LinkData::Flow flowLeg3 = LinkData::flowBothEnds(branchAdmittanceLeg3.y11, branchAdmittanceLeg3.y12, + branchAdmittanceLeg3.y21, branchAdmittanceLeg3.y22, v3, v0); + + m_computedP1 = std::real(flowLeg1.fromTo); + m_computedQ1 = std::imag(flowLeg1.fromTo); + m_computedP2 = std::real(flowLeg2.fromTo); + m_computedQ2 = std::imag(flowLeg2.fromTo); + m_computedP3 = std::real(flowLeg3.fromTo); + m_computedQ3 = std::imag(flowLeg3.fromTo); + + m_starU = std::abs(v0); + m_starTheta = std::arg(v0); +} + +LinkData::BranchAdmittanceMatrix TwtData::calculateTwoConnectedLegsAdmittance(const LinkData::BranchAdmittanceMatrix& firstCloseLeg, + const LinkData::BranchAdmittanceMatrix& secondCloseLeg, + const LinkData::BranchAdmittanceMatrix& openLeg) const { + + std::complex ysh = LinkData::kronAntenna(openLeg.y11, openLeg.y12, openLeg.y21, openLeg.y22, true); + LinkData::BranchAdmittanceMatrix secondCloseLegMod(secondCloseLeg.y11, secondCloseLeg.y12, secondCloseLeg.y21, secondCloseLeg.y22 + ysh); + return LinkData::kronChain(firstCloseLeg, Branch::Side::TWO, secondCloseLegMod, Branch::Side::TWO); +} + +LinkData::Flow TwtData::calculateTwoConnectedLegsFlow(double u1, double theta1, double u2, double theta2, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg1, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg2, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixOpenLeg) const { + + std::complex v1 = math::ComplexUtils::polar2Complex(u1, theta1); + std::complex v2 = math::ComplexUtils::polar2Complex(u2, theta2); + + LinkData::BranchAdmittanceMatrix admittance = calculateTwoConnectedLegsAdmittance(admittanceMatrixLeg1, + admittanceMatrixLeg2, admittanceMatrixOpenLeg); + + return LinkData::flowBothEnds(admittance.y11, admittance.y12, admittance.y21, admittance.y22, v1, v2); +} + +std::complex TwtData::calculateTwoConnectedLegsStarBusVoltage(double u1, double theta1, double u2, double theta2, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg1, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixLeg2, + const LinkData::BranchAdmittanceMatrix& admittanceMatrixOpenLeg) const { + + std::complex v1 = math::ComplexUtils::polar2Complex(u1, theta1); + std::complex v2 = math::ComplexUtils::polar2Complex(u2, theta2); + + std::complex yshO = LinkData::kronAntenna(admittanceMatrixOpenLeg.y11, admittanceMatrixOpenLeg.y12, admittanceMatrixOpenLeg.y21, admittanceMatrixOpenLeg.y22, true); + return -(admittanceMatrixLeg1.y21 * v1 + (admittanceMatrixLeg2.y21 * v2)) + / (admittanceMatrixLeg1.y22 + admittanceMatrixLeg2.y22 + yshO); +} + +double TwtData::getB1(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_b11; + case ThreeWindingsTransformer::Side::TWO: + return m_b21; + case ThreeWindingsTransformer::Side::THREE: + return m_b31; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getB1(const ThreeWindingsTransformer::Leg& leg, bool twtSplitShuntAdmittance) { + return getValue(twtSplitShuntAdmittance ? leg.getB() / 2 : leg.getB(), + leg.getOptionalRatioTapChanger() ? leg.getRatioTapChanger().getCurrentStep().getB() : 0.0, + leg.getOptionalPhaseTapChanger() ? leg.getPhaseTapChanger().getCurrentStep().getB() : 0.0); +} + +double TwtData::getB2(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_b12; + case ThreeWindingsTransformer::Side::TWO: + return m_b22; + case ThreeWindingsTransformer::Side::THREE: + return m_b32; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getB2(const ThreeWindingsTransformer::Leg& leg, bool twtSplitShuntAdmittance) { + return getValue(twtSplitShuntAdmittance ? leg.getB() / 2 : 0.0, + leg.getOptionalRatioTapChanger() ? leg.getRatioTapChanger().getCurrentStep().getB() : 0.0, + leg.getOptionalPhaseTapChanger() ? leg.getPhaseTapChanger().getCurrentStep().getB() : 0.0); +} + +double TwtData::getComputedP(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_computedP1; + case ThreeWindingsTransformer::Side::TWO: + return m_computedP2; + case ThreeWindingsTransformer::Side::THREE: + return m_computedP3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getComputedQ(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_computedQ1; + case ThreeWindingsTransformer::Side::TWO: + return m_computedQ2; + case ThreeWindingsTransformer::Side::THREE: + return m_computedQ3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getG1(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_g11; + case ThreeWindingsTransformer::Side::TWO: + return m_g21; + case ThreeWindingsTransformer::Side::THREE: + return m_g31; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getG1(const ThreeWindingsTransformer::Leg& leg, bool twtSplitShuntAdmittance) { + return getValue(twtSplitShuntAdmittance ? leg.getG() / 2 : leg.getG(), + leg.getOptionalRatioTapChanger() ? leg.getRatioTapChanger().getCurrentStep().getG() : 0.0, + leg.getOptionalPhaseTapChanger() ? leg.getPhaseTapChanger().getCurrentStep().getG() : 0.0); +} + +double TwtData::getG2(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_g12; + case ThreeWindingsTransformer::Side::TWO: + return m_g22; + case ThreeWindingsTransformer::Side::THREE: + return m_g32; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getG2(const ThreeWindingsTransformer::Leg& leg, bool twtSplitShuntAdmittance) { + return getValue(twtSplitShuntAdmittance ? leg.getG() / 2 : 0.0, + leg.getOptionalRatioTapChanger() ? leg.getRatioTapChanger().getCurrentStep().getG() : 0.0, + leg.getOptionalPhaseTapChanger() ? leg.getPhaseTapChanger().getCurrentStep().getG() : 0.0); +} + +const std::string& TwtData::getId() const { + return m_id; +} + +double TwtData::getP(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_p1; + case ThreeWindingsTransformer::Side::TWO: + return m_p2; + case ThreeWindingsTransformer::Side::THREE: + return m_p3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +int TwtData::getPhaseAngleClock2() const { + return m_phaseAngleClock2; +} + +int TwtData::getPhaseAngleClock3() const { + return m_phaseAngleClock3; +} + +double TwtData::getQ(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_q1; + case ThreeWindingsTransformer::Side::TWO: + return m_q2; + case ThreeWindingsTransformer::Side::THREE: + return m_q3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getR(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_r1; + case ThreeWindingsTransformer::Side::TWO: + return m_r2; + case ThreeWindingsTransformer::Side::THREE: + return m_r3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getR(const ThreeWindingsTransformer::Leg& leg) { + return getValue(leg.getR(), + leg.getOptionalRatioTapChanger() ? leg.getRatioTapChanger().getCurrentStep().getR() : 0.0, + leg.getOptionalPhaseTapChanger() ? leg.getPhaseTapChanger().getCurrentStep().getR() : 0.0); +} + +double TwtData::getRatedU(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_ratedU1; + case ThreeWindingsTransformer::Side::TWO: + return m_ratedU2; + case ThreeWindingsTransformer::Side::THREE: + return m_ratedU3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getRatedU0() const { + return m_ratedU0; +} + +double TwtData::getStarTheta() const { + return m_starTheta; +} + +double TwtData::getStarU() const { + return m_starU; +} + +double TwtData::getTheta(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_theta1; + case ThreeWindingsTransformer::Side::TWO: + return m_theta2; + case ThreeWindingsTransformer::Side::THREE: + return m_theta3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getTheta(const ThreeWindingsTransformer::Leg& leg) { + return leg.getTerminal().isConnected() ? stdcxx::toRadians * leg.getTerminal().getBusView().getBus().get().getAngle() : stdcxx::nan(); +} + +double TwtData::getU(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_u1; + case ThreeWindingsTransformer::Side::TWO: + return m_u2; + case ThreeWindingsTransformer::Side::THREE: + return m_u3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getV(const ThreeWindingsTransformer::Leg& leg) { + return leg.getTerminal().isConnected() ? leg.getTerminal().getBusView().getBus().get().getV() : stdcxx::nan(); +} + +double TwtData::getValue(double initialValue, double rtcStepValue, double ptcStepValue) { + return initialValue * (1 + rtcStepValue / 100) * (1 + ptcStepValue / 100); +} + +double TwtData::getX(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_x1; + case ThreeWindingsTransformer::Side::TWO: + return m_x2; + case ThreeWindingsTransformer::Side::THREE: + return m_x3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +double TwtData::getX(const ThreeWindingsTransformer::Leg& leg) { + return getValue(leg.getX(), + leg.getOptionalRatioTapChanger() ? leg.getRatioTapChanger().getCurrentStep().getX() : 0.0, + leg.getOptionalPhaseTapChanger() ? leg.getPhaseTapChanger().getCurrentStep().getX() : 0.0); +} + +bool TwtData::isConnected(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_connected1; + case ThreeWindingsTransformer::Side::TWO: + return m_connected2; + case ThreeWindingsTransformer::Side::THREE: + return m_connected3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +bool TwtData::isMainComponent(const ThreeWindingsTransformer::Side& side) const { + switch (side) { + case ThreeWindingsTransformer::Side::ONE: + return m_mainComponent1; + case ThreeWindingsTransformer::Side::TWO: + return m_mainComponent2; + case ThreeWindingsTransformer::Side::THREE: + return m_mainComponent3; + default: + throw AssertionError(UNEXPECTED_SIDE + ": " + Enum::toString(side)); + } +} + +bool TwtData::isMainComponent(const ThreeWindingsTransformer::Leg& leg) { + const auto& bus = leg.getTerminal().getBusView().getBus(); + const auto& connectableBus = leg.getTerminal().getBusView().getConnectableBus(); + bool connectableMainComponent = connectableBus && connectableBus.get().isInMainConnectedComponent(); + return bus ? bus.get().isInMainConnectedComponent() : connectableMainComponent; +} + +double TwtData::rho(const ThreeWindingsTransformer::Leg& leg, double ratedU0) { + double rho = ratedU0 / leg.getRatedU(); + rho *= leg.getOptionalRatioTapChanger() ? leg.getRatioTapChanger().getCurrentStep().getRho() : 1.0; + rho *= leg.getOptionalPhaseTapChanger() ? leg.getPhaseTapChanger().getCurrentStep().getRho() : 1.0; + return rho; +} + +bool TwtData::valid(double voltage, double theta) { + if (std::isnan(voltage) || voltage <= 0.0) { + return false; + } + return !std::isnan(theta); +} + +} // namespace iidm + +} // namespace powsybl diff --git a/test/iidm/CMakeLists.txt b/test/iidm/CMakeLists.txt index 150dc3d2a..f858e3425 100644 --- a/test/iidm/CMakeLists.txt +++ b/test/iidm/CMakeLists.txt @@ -61,6 +61,7 @@ set(UNIT_TEST_SOURCES util/DanglingLineDataTest.cpp util/SVTest.cpp util/TerminalFinderTest.cpp + util/TwtDataTest.cpp ) add_executable(unit-tests-iidm ${UNIT_TEST_SOURCES}) diff --git a/test/iidm/util/DanglingLineDataTest.cpp b/test/iidm/util/DanglingLineDataTest.cpp index e165ab661..a2a1316a9 100644 --- a/test/iidm/util/DanglingLineDataTest.cpp +++ b/test/iidm/util/DanglingLineDataTest.cpp @@ -53,10 +53,7 @@ Network createDanglingDataNetwork() { bool dlCompareBoundaryBusVoltage(const DanglingLineData& dlData, double boundaryBusU, double boundaryBusAngle) { double tol = 0.00001; - if (std::abs(dlData.getBoundaryBusU() - boundaryBusU) > tol || std::abs(stdcxx::toDegrees * dlData.getBoundaryBusTheta() - boundaryBusAngle) > tol) { - return false; - } - return true; + return !(std::abs(dlData.getBoundaryBusU() - boundaryBusU) > tol || std::abs(stdcxx::toDegrees * dlData.getBoundaryBusTheta() - boundaryBusAngle) > tol); } BOOST_AUTO_TEST_SUITE(DanglingLineDataTestSuite) diff --git a/test/iidm/util/TwtDataTest.cpp b/test/iidm/util/TwtDataTest.cpp new file mode 100644 index 000000000..2f8f9e050 --- /dev/null +++ b/test/iidm/util/TwtDataTest.cpp @@ -0,0 +1,496 @@ +/** + * Copyright (c) 2022, RTE (http://www.rte-france.com) + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at http://mozilla.org/MPL/2.0/. + */ + +#include + +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace powsybl { + +namespace iidm { + +double P1 = 99.218431; +double Q1 = 2.7304328; +double P2 = -216.19819; +double Q2 = -85.368180; +double P3 = 118; +double Q3 = 92.612077; + +double U1 = 412.989001; +double ANGLE1 = -6.78071; +double U2 = 224.315268; +double ANGLE2 = -8.77012; +double U3 = 21.987; +double ANGLE3 = -6.6508; + +double STAR_U = 412.662007016922757; +double STAR_ANGLE = -7.353686938578365; + +double R1 = 0.898462; +double X1 = 17.204128; +double G11 = 0; +double B11 = 2.4375E-6; +double G12 = 0; +double B12 = 0; +double RATED_U1 = 400; +double R2 = 1.070770247933884; +double X2 = 19.6664; +double G21 = 0; +double B21 = 0; +double G22 = 0; +double B22 = 0; +double RATED_U2 = 220; +double R3 = 4.837006802721089; +double X3 = 21.76072562358277; +double G31 = 0; +double B31 = 0; +double G32 = 0; +double B32 = 0; +double RATED_U3 = 21; +int PHASE_ANGLE_CLOCK_2 = 0; +int PHASE_ANGLE_CLOCK_3 = 0; +double RATED_U0 = RATED_U1; + +Network createTwtDataTestNetwork() { + Network network("test", "test"); + Substation& substation = network.newSubstation() + .setId("S1") + .setName("S1_NAME") + .setCountry(Country::FR) + .setTso("TSO") + .add(); + + VoltageLevel& vl1 = substation.newVoltageLevel() + .setId("VL1") + .setName("VL1_NAME") + .setTopologyKind(TopologyKind::BUS_BREAKER) + .setNominalV(380.0) + .setLowVoltageLimit(340.0) + .setHighVoltageLimit(420.0) + .add(); + + Bus& vl1Bus1 = vl1.getBusBreakerView().newBus() + .setId("VL1_BUS1") + .add(); + + vl1.newLoad() + .setId("LOAD1") + .setName("LOAD1_NAME") + .setBus(vl1Bus1.getId()) + .setConnectableBus(vl1Bus1.getId()) + .setLoadType(LoadType::UNDEFINED) + .setP0(50.0) + .setQ0(40.0) + .add(); + + VoltageLevel& vl2 = substation.newVoltageLevel() + .setId("VL2") + .setName("VL2_NAME") + .setTopologyKind(TopologyKind::BUS_BREAKER) + .setNominalV(225.0) + .setLowVoltageLimit(200.0) + .setHighVoltageLimit(260.0) + .add(); + + Bus& vl2Bus1 = vl2.getBusBreakerView().newBus() + .setId("VL2_BUS1") + .add(); + + vl2.newLoad() + .setId("LOAD2") + .setName("LOAD2_NAME") + .setBus(vl2Bus1.getId()) + .setConnectableBus(vl2Bus1.getId()) + .setLoadType(LoadType::UNDEFINED) + .setP0(60.0) + .setQ0(70.0) + .add(); + + VoltageLevel& vl3 = substation.newVoltageLevel() + .setId("VL3") + .setName("VL3_NAME") + .setTopologyKind(TopologyKind::BUS_BREAKER) + .setNominalV(380.0) + .setLowVoltageLimit(340.0) + .setHighVoltageLimit(420.0) + .add(); + + Bus& vl3Bus1 = vl3.getBusBreakerView().newBus() + .setId("VL3_BUS1") + .add(); + + Substation& substation2 = network.newSubstation() + .setId("S2") + .setName("S2_NAME") + .setCountry(Country::FR) + .setTso("TSO") + .add(); + + VoltageLevel& vl4 = substation2.newVoltageLevel() + .setId("VL4") + .setName("VL4_NAME") + .setTopologyKind(TopologyKind::BUS_BREAKER) + .setNominalV(225.0) + .setLowVoltageLimit(200.0) + .setHighVoltageLimit(260.0) + .add(); + + vl4.getBusBreakerView().newBus() + .setId("VL4_BUS1") + .add(); + + substation.newThreeWindingsTransformer() + .setId("3WT_VL1_VL2_VL3") + .setName("3WT_VL1_VL2_VL3_NAME") + .setRatedU0(RATED_U0) + .newLeg1() + .setR(1.3) + .setX(1.4) + .setG(1.6) + .setB(1.7) + .setRatedU(1.1) + .setRatedS(2.2) + .setVoltageLevel(vl1.getId()) + .setBus(vl1Bus1.getId()) + .setConnectableBus(vl1Bus1.getId()) + .add() + .newLeg2() + .setR(2.3) + .setX(2.4) + .setG(0.0) + .setB(0.0) + .setRatedU(2.1) + .setVoltageLevel(vl2.getId()) + .setBus(vl2Bus1.getId()) + .setConnectableBus(vl2Bus1.getId()) + .add() + .newLeg3() + .setR(3.3) + .setX(3.4) + .setG(0.0) + .setB(0.0) + .setRatedU(3.1) + .setVoltageLevel(vl3.getId()) + .setBus(vl3Bus1.getId()) + .setConnectableBus(vl3Bus1.getId()) + .add() + .add(); + + ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + + twt.getLeg1().getTerminal().getBusBreakerView().getBus().get().setV(U1); + twt.getLeg1().getTerminal().getBusBreakerView().getBus().get().setAngle(ANGLE1); + twt.getLeg1().getTerminal().setP(P1); + twt.getLeg1().getTerminal().setQ(Q1); + twt.getLeg1().setR(R1); + twt.getLeg1().setX(X1); + twt.getLeg1().setRatedU(RATED_U1); + twt.getLeg1().setB(B11 + B12); + twt.getLeg1().setG(G11 + G12); + + twt.getLeg2().getTerminal().getBusBreakerView().getBus().get().setV(U2); + twt.getLeg2().getTerminal().getBusBreakerView().getBus().get().setAngle(ANGLE2); + twt.getLeg2().getTerminal().setP(P2); + twt.getLeg2().getTerminal().setQ(Q2); + twt.getLeg2().setR(R2); + twt.getLeg2().setX(X2); + twt.getLeg2().setRatedU(RATED_U2); + twt.getLeg2().setB(B21 + B22); + twt.getLeg2().setG(G21 + G22); + + twt.getLeg3().getTerminal().getBusBreakerView().getBus().get().setV(U3); + twt.getLeg3().getTerminal().getBusBreakerView().getBus().get().setAngle(ANGLE3); + twt.getLeg3().getTerminal().setP(P3); + twt.getLeg3().getTerminal().setQ(Q3); + twt.getLeg3().setR(R3); + twt.getLeg3().setX(X3); + twt.getLeg3().setRatedU(RATED_U3); + twt.getLeg3().setB(B31 + B32); + twt.getLeg3().setG(G31 + G32); + + return network; +} + +struct T3xFlow { + double p1 = stdcxx::nan(); + double q1 = stdcxx::nan(); + double p2 = stdcxx::nan(); + double q2 = stdcxx::nan(); + double p3 = stdcxx::nan(); + double q3 = stdcxx::nan(); +}; + +bool sameFlow(const T3xFlow& expected, const T3xFlow& actual) { + double tol = 0.00001; + return !((!std::isnan(expected.p1) && std::isnan(actual.p1)) || + (!std::isnan(expected.q1) && std::isnan(actual.q1)) || + (!std::isnan(expected.p2) && std::isnan(actual.p2)) || + (!std::isnan(expected.q2) && std::isnan(actual.q2)) || + (!std::isnan(expected.p3) && std::isnan(actual.p3)) || + (!std::isnan(expected.q3) && std::isnan(actual.q3)) || + std::abs(expected.p1 - actual.p1) > tol || + std::abs(expected.q1 - actual.q1) > tol || + std::abs(expected.p2 - actual.p2) > tol || + std::abs(expected.q2 - actual.q2) > tol || + std::abs(expected.p3 - actual.p3) > tol || + std::abs(expected.q3 - actual.q3) > tol); +} + +bool t3xCompareFlow(const TwtData& twtData, double p1, double q1, double p2, double q2, double p3, double q3) { + T3xFlow actual; + actual.p1 = twtData.getComputedP(ThreeWindingsTransformer::Side::ONE); + actual.q1 = twtData.getComputedQ(ThreeWindingsTransformer::Side::ONE); + actual.p2 = twtData.getComputedP(ThreeWindingsTransformer::Side::TWO); + actual.q2 = twtData.getComputedQ(ThreeWindingsTransformer::Side::TWO); + actual.p3 = twtData.getComputedP(ThreeWindingsTransformer::Side::THREE); + actual.q3 = twtData.getComputedQ(ThreeWindingsTransformer::Side::THREE); + + T3xFlow expected; + expected.p1 = p1; + expected.q1 = q1; + expected.p2 = p2; + expected.q2 = q2; + expected.p3 = p3; + expected.q3 = q3; + + return sameFlow(expected, actual); +} + +bool t3xCompareStarBusVoltage(const TwtData& twtData, double starU, double starAngle) { + double tol = 0.00001; + return !(std::isnan(twtData.getStarU()) && !std::isnan(starU)) || + (std::isnan(twtData.getStarTheta()) && !std::isnan(starAngle)) || + std::abs(twtData.getStarU() - starU) > tol || + std::abs(stdcxx::toDegrees * twtData.getStarTheta() - starAngle) > tol; +} + +BOOST_AUTO_TEST_SUITE(TwtDataTestSuite) + +BOOST_AUTO_TEST_CASE(test) { + Network network = createTwtDataTestNetwork(); + const ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + TwtData twtData(twt, 0, false, false); + + BOOST_CHECK_CLOSE(P1, twtData.getComputedP(ThreeWindingsTransformer::Side::ONE), 0.3); + BOOST_CHECK_CLOSE(Q1, twtData.getComputedQ(ThreeWindingsTransformer::Side::ONE), 0.7); // 2.7304328 != 2.7471471852083678 + BOOST_CHECK_CLOSE(P2, twtData.getComputedP(ThreeWindingsTransformer::Side::TWO), 0.3); + BOOST_CHECK_CLOSE(Q2, twtData.getComputedQ(ThreeWindingsTransformer::Side::TWO), 0.3); + BOOST_CHECK_CLOSE(P3, twtData.getComputedP(ThreeWindingsTransformer::Side::THREE), 0.3); + BOOST_CHECK_CLOSE(Q3, twtData.getComputedQ(ThreeWindingsTransformer::Side::THREE), 0.3); + + BOOST_CHECK_CLOSE(STAR_U, twtData.getStarU(), 0.0001); + BOOST_CHECK_CLOSE(STAR_ANGLE, stdcxx::toDegrees * twtData.getStarTheta(), 0.0001); +} + +BOOST_AUTO_TEST_CASE(testSplitShuntAdmittance) { + Network network = createTwtDataTestNetwork(); + const ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + TwtData twtData(twt, 0, false, true); + + bool ok = t3xCompareFlow(twtData, 99.231950, 2.876479, -216.194348, -85.558437, 117.981856, 92.439531); + BOOST_CHECK(ok); +} + +BOOST_AUTO_TEST_CASE(testEnd1End2End3Connected) { + Network network = createTwtDataTestNetwork(); + const ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + TwtData twtData(twt, 0, false); + + bool ok = t3xCompareFlow(twtData, 99.227288294050, 2.747147185208, -216.195866533486, -85.490493190353, 117.988318295633, 92.500849015581); + BOOST_CHECK(ok); + + ok = t3xCompareStarBusVoltage(twtData, 412.66200701692287, -7.353686938578365); + BOOST_CHECK(ok); +} + +BOOST_AUTO_TEST_CASE(testEnd2End3Connected) { + Network network = createTwtDataTestNetwork(); + ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + twt.getLeg1().getTerminal().disconnect(); + + TwtData twtData(twt, 0, false); + bool ok = t3xCompareFlow(twtData, 0.0, 0.0, -164.099476216398, -81.835885442800, 165.291731946141, 89.787051339157); + BOOST_CHECK(ok); + + ok = t3xCompareStarBusVoltage(twtData, 412.29478568401856, -7.700275244269859); + BOOST_CHECK(ok); +} + +BOOST_AUTO_TEST_CASE(testEnd1End3Connected) { + Network network = createTwtDataTestNetwork(); + ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + twt.getLeg2().getTerminal().disconnect(); + + TwtData twtData(twt, 0, false); + bool ok = t3xCompareFlow(twtData, -18.723067158829, -59.239225729782, 0.0, 0.0, 18.851212571411, 59.694062940578); + BOOST_CHECK(ok); + + ok = t3xCompareStarBusVoltage(twtData, 415.4806896701992, -6.690799426080698); + BOOST_CHECK(ok); +} + +BOOST_AUTO_TEST_CASE(testEnd1End2Connected) { + Network network = createTwtDataTestNetwork(); + ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + twt.getLeg3().getTerminal().disconnect(); + + TwtData twtData(twt, 0, false); + bool ok = t3xCompareFlow(twtData, 161.351352526949, 51.327798049323, -161.019856627996, -45.536840365345, 0.0, 0.0); + BOOST_CHECK(ok); + + ok = t3xCompareStarBusVoltage(twtData, 410.53566804098494, -7.703116461849692); + BOOST_CHECK(ok); +} + +BOOST_AUTO_TEST_CASE(testEnd1Connected) { + Network network = createTwtDataTestNetwork(); + ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + twt.getLeg2().getTerminal().disconnect(); + twt.getLeg3().getTerminal().disconnect(); + + TwtData twtData(twt, 0, false); + bool ok = t3xCompareFlow(twtData, 0.0, -0.415739792683, 0.0, 0.0, 0.0, 0.0); + BOOST_CHECK(ok); + + ok = t3xCompareStarBusVoltage(twtData, 412.9890009999999, -6.78071000000000); + BOOST_CHECK(ok); +} + +BOOST_AUTO_TEST_CASE(testEnd2Connected) { + Network network = createTwtDataTestNetwork(); + ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + twt.getLeg1().getTerminal().disconnect(); + twt.getLeg3().getTerminal().disconnect(); + + TwtData twtData(twt, 0, false); + bool ok = t3xCompareFlow(twtData, 0.0, 0.0, 0.000001946510, -0.405486077928, 0.0, 0.0); + BOOST_CHECK(ok); + + ok = t3xCompareStarBusVoltage(twtData, 407.8654944214268, -8.77026956158324); + BOOST_CHECK(ok); +} + +BOOST_AUTO_TEST_CASE(testEnd3Connected) { + Network network = createTwtDataTestNetwork(); + ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + twt.getLeg1().getTerminal().disconnect(); + twt.getLeg2().getTerminal().disconnect(); + + TwtData twtData(twt, 0, false); + bool ok = t3xCompareFlow(twtData, 0.0, 0.0, 0.0, 0.0, 0.000005977974, -0.427562118410); + BOOST_CHECK(ok); + + ok = t3xCompareStarBusVoltage(twtData, 418.82221596280823, -6.65147559975559); + BOOST_CHECK(ok); +} + +BOOST_AUTO_TEST_CASE(testEnd1End2End3Disconnected) { + Network network = createTwtDataTestNetwork(); + ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + twt.getLeg1().getTerminal().disconnect(); + twt.getLeg2().getTerminal().disconnect(); + twt.getLeg3().getTerminal().disconnect(); + + TwtData twtData(twt, 0, false); + bool ok = t3xCompareFlow(twtData, stdcxx::nan(), stdcxx::nan(), stdcxx::nan(), stdcxx::nan(), stdcxx::nan(), stdcxx::nan()); + BOOST_CHECK(ok); + + ok = t3xCompareStarBusVoltage(twtData, stdcxx::nan(), stdcxx::nan()); + BOOST_CHECK(ok); +} + +BOOST_AUTO_TEST_CASE(checkInputs) { + Network network = createTwtDataTestNetwork(); + ThreeWindingsTransformer& twt = network.getThreeWindingsTransformer("3WT_VL1_VL2_VL3"); + + BOOST_CHECK_CLOSE(2.4375e-6, TwtData(twt, 0, false).getB1(ThreeWindingsTransformer::Side::ONE), std::numeric_limits::epsilon()); + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getB1(ThreeWindingsTransformer::Side::TWO)); + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getB1(ThreeWindingsTransformer::Side::THREE)); + + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getB2(ThreeWindingsTransformer::Side::ONE)); + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getB2(ThreeWindingsTransformer::Side::TWO)); + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getB2(ThreeWindingsTransformer::Side::THREE)); + + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getG1(ThreeWindingsTransformer::Side::ONE)); + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getG1(ThreeWindingsTransformer::Side::TWO)); + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getG1(ThreeWindingsTransformer::Side::THREE)); + + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getG2(ThreeWindingsTransformer::Side::ONE)); + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getG2(ThreeWindingsTransformer::Side::TWO)); + BOOST_CHECK_EQUAL(0.0, TwtData(twt, 0, false).getG2(ThreeWindingsTransformer::Side::THREE)); + + BOOST_CHECK_EQUAL("3WT_VL1_VL2_VL3", TwtData(twt, 0, false).getId()); + + BOOST_CHECK_EQUAL(0, TwtData(twt, 0, false).getPhaseAngleClock2()); + BOOST_CHECK_EQUAL(0, TwtData(twt, 0, false).getPhaseAngleClock3()); + + BOOST_CHECK_CLOSE(2.7304328, TwtData(twt, 0, false).getQ(ThreeWindingsTransformer::Side::ONE), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(-85.36818, TwtData(twt, 0, false).getQ(ThreeWindingsTransformer::Side::TWO), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(92.612077, TwtData(twt, 0, false).getQ(ThreeWindingsTransformer::Side::THREE), std::numeric_limits::epsilon()); + + BOOST_CHECK_CLOSE(0.898462, TwtData(twt, 0, false).getR(ThreeWindingsTransformer::Side::ONE), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(1.070770247933884, TwtData(twt, 0, false).getR(ThreeWindingsTransformer::Side::TWO), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(4.837006802721089, TwtData(twt, 0, false).getR(ThreeWindingsTransformer::Side::THREE), std::numeric_limits::epsilon()); + + BOOST_CHECK_CLOSE(400.0, TwtData(twt, 0, false).getRatedU(ThreeWindingsTransformer::Side::ONE), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(220.0, TwtData(twt, 0, false).getRatedU(ThreeWindingsTransformer::Side::TWO), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(21.0, TwtData(twt, 0, false).getRatedU(ThreeWindingsTransformer::Side::THREE), std::numeric_limits::epsilon()); + + BOOST_CHECK_CLOSE(400, TwtData(twt, 0, false).getRatedU0(), std::numeric_limits::epsilon()); + + BOOST_CHECK_CLOSE(-0.1183457151229047, TwtData(twt, 0, false).getTheta(ThreeWindingsTransformer::Side::ONE), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(-0.1530674697950051, TwtData(twt, 0, false).getTheta(ThreeWindingsTransformer::Side::TWO), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(-0.11607835789163888, TwtData(twt, 0, false).getTheta(ThreeWindingsTransformer::Side::THREE), std::numeric_limits::epsilon()); + + BOOST_CHECK_CLOSE(412.989001, TwtData(twt, 0, false).getU(ThreeWindingsTransformer::Side::ONE), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(224.315268, TwtData(twt, 0, false).getU(ThreeWindingsTransformer::Side::TWO), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(21.987, TwtData(twt, 0, false).getU(ThreeWindingsTransformer::Side::THREE), std::numeric_limits::epsilon()); + + BOOST_CHECK_CLOSE(17.204128, TwtData(twt, 0, false).getX(ThreeWindingsTransformer::Side::ONE), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(19.6664, TwtData(twt, 0, false).getX(ThreeWindingsTransformer::Side::TWO), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(21.76072562358277, TwtData(twt, 0, false).getX(ThreeWindingsTransformer::Side::THREE), std::numeric_limits::epsilon()); + + BOOST_CHECK_CLOSE(99.218431, TwtData(twt, 0, false).getP(ThreeWindingsTransformer::Side::ONE), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(-216.19819, TwtData(twt, 0, false).getP(ThreeWindingsTransformer::Side::TWO), std::numeric_limits::epsilon()); + BOOST_CHECK_CLOSE(118.0, TwtData(twt, 0, false).getP(ThreeWindingsTransformer::Side::THREE), std::numeric_limits::epsilon()); + + BOOST_CHECK(TwtData(twt, 0, false).isConnected(ThreeWindingsTransformer::Side::ONE)); + BOOST_CHECK(TwtData(twt, 0, false).isMainComponent(ThreeWindingsTransformer::Side::ONE)); + twt.getLeg1().getTerminal().disconnect(); + BOOST_CHECK(!TwtData(twt, 0, false).isConnected(ThreeWindingsTransformer::Side::ONE)); + BOOST_CHECK(!TwtData(twt, 0, false).isMainComponent(ThreeWindingsTransformer::Side::ONE)); + + BOOST_CHECK(TwtData(twt, 0, false).isConnected(ThreeWindingsTransformer::Side::TWO)); + BOOST_CHECK(TwtData(twt, 0, false).isMainComponent(ThreeWindingsTransformer::Side::TWO)); + twt.getLeg2().getTerminal().disconnect(); + BOOST_CHECK(!TwtData(twt, 0, false).isConnected(ThreeWindingsTransformer::Side::TWO)); + BOOST_CHECK(!TwtData(twt, 0, false).isMainComponent(ThreeWindingsTransformer::Side::TWO)); + + BOOST_CHECK(TwtData(twt, 0, false).isConnected(ThreeWindingsTransformer::Side::THREE)); + BOOST_CHECK(TwtData(twt, 0, false).isMainComponent(ThreeWindingsTransformer::Side::THREE)); + twt.getLeg3().getTerminal().disconnect(); + BOOST_CHECK(!TwtData(twt, 0, false).isConnected(ThreeWindingsTransformer::Side::THREE)); + BOOST_CHECK(!TwtData(twt, 0, false).isMainComponent(ThreeWindingsTransformer::Side::THREE)); +} + +BOOST_AUTO_TEST_SUITE_END() + +} // namespace iidm + +} // namespace powsybl From ffa4c163aa2a92b5e10a1ba76094e65d591cccf1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20LAIGRE?= Date: Tue, 24 May 2022 12:10:06 +0200 Subject: [PATCH 3/3] Review MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: Sébastien LAIGRE --- .../powsybl/iidm/util/DanglingLineData.hpp | 14 +------ include/powsybl/iidm/util/LinkData.hpp | 3 +- src/iidm/util/DanglingLineData.cpp | 40 +++++++++---------- src/iidm/util/LinkData.cpp | 10 ++--- src/iidm/util/TwtData.cpp | 4 +- 5 files changed, 29 insertions(+), 42 deletions(-) diff --git a/include/powsybl/iidm/util/DanglingLineData.hpp b/include/powsybl/iidm/util/DanglingLineData.hpp index add4204f7..fd06e5066 100644 --- a/include/powsybl/iidm/util/DanglingLineData.hpp +++ b/include/powsybl/iidm/util/DanglingLineData.hpp @@ -38,19 +38,7 @@ class DanglingLineData { static bool valid(double v, double theta); private: - std::string m_id; - - double m_r; - double m_x; - double m_g1; - double m_g2; - double m_b1; - double m_b2; - - double m_u1; - double m_theta1; - double m_p0; - double m_q0; + const DanglingLine& m_danglingLine; double m_boundaryBusU = stdcxx::nan(); double m_boundaryBusTheta = stdcxx::nan(); diff --git a/include/powsybl/iidm/util/LinkData.hpp b/include/powsybl/iidm/util/LinkData.hpp index df23aef05..e32c53dae 100644 --- a/include/powsybl/iidm/util/LinkData.hpp +++ b/include/powsybl/iidm/util/LinkData.hpp @@ -25,12 +25,13 @@ struct Flow { std::complex toFrom; }; -struct BranchAdmittanceMatrix { +class BranchAdmittanceMatrix { public: BranchAdmittanceMatrix() = default; BranchAdmittanceMatrix(const std::complex& y11, const std::complex& y12, const std::complex& y21, const std::complex& y22); +public: std::complex y11; std::complex y12; diff --git a/src/iidm/util/DanglingLineData.cpp b/src/iidm/util/DanglingLineData.cpp index 28dc52fd9..ba9db97ab 100644 --- a/src/iidm/util/DanglingLineData.cpp +++ b/src/iidm/util/DanglingLineData.cpp @@ -23,33 +23,33 @@ DanglingLineData::DanglingLineData(const DanglingLine& danglingLine) : } DanglingLineData::DanglingLineData(const DanglingLine& danglingLine, bool splitShuntAdmittance) : - m_id(danglingLine.getId()), - m_r(danglingLine.getR()), - m_x(danglingLine.getX()), - m_g1(splitShuntAdmittance ? danglingLine.getG() * 0.5 : danglingLine.getG()), - m_g2(splitShuntAdmittance ? danglingLine.getG() * 0.5 : 0.0), - m_b1(splitShuntAdmittance ? danglingLine.getB() * 0.5 : danglingLine.getB()), - m_b2(splitShuntAdmittance ? danglingLine.getB() * 0.5 : 0.0), - m_u1(getV(danglingLine)), - m_theta1(getTheta(danglingLine)), - m_p0(danglingLine.getP0()), - m_q0(danglingLine.getQ0()) { - - if (!valid(m_u1, m_theta1)) { + m_danglingLine(danglingLine){ + + double g1 = splitShuntAdmittance ? danglingLine.getG() * 0.5 : danglingLine.getG(); + double g2 = splitShuntAdmittance ? danglingLine.getG() * 0.5 : 0.0; + double b1 = splitShuntAdmittance ? danglingLine.getB() * 0.5 : danglingLine.getB(); + double b2 = splitShuntAdmittance ? danglingLine.getB() * 0.5 : 0.0; + + double u1 = getV(danglingLine); + double theta1 = getTheta(danglingLine); + + if (!valid(u1, theta1)) { + m_boundaryBusU = stdcxx::nan(); + m_boundaryBusTheta = stdcxx::nan(); return; } - const std::complex& v1 = math::ComplexUtils::polar2Complex(m_u1, m_theta1); + const std::complex& v1 = math::ComplexUtils::polar2Complex(u1, theta1); std::complex vBoundaryBus(stdcxx::nan(), stdcxx::nan()); - if (m_p0 == 0.0 && m_q0 == 0.0) { - LinkData::BranchAdmittanceMatrix adm = LinkData::calculateBranchAdmittance(m_r, m_x, 1.0, 0.0, 1.0, 0.0, std::complex(m_g1, m_b1), std::complex(m_g2, m_b2)); + if (danglingLine.getP0() == 0.0 && danglingLine.getQ0() == 0.0) { + LinkData::BranchAdmittanceMatrix adm = LinkData::calculateBranchAdmittance(danglingLine.getR(), danglingLine.getX(), 1.0, 0.0, 1.0, 0.0, std::complex(g1, b1), std::complex(g2, b2)); vBoundaryBus = (- adm.y21 * v1) / adm.y22; } else { // Two buses Loadflow - std::complex sBoundary(-m_p0, -m_q0); - std::complex ytr = 1.0 / std::complex(m_r, m_x); - std::complex ysh2(m_g2, m_b2); + std::complex sBoundary(-danglingLine.getP0(), -danglingLine.getQ0()); + std::complex ytr = 1.0 / std::complex(danglingLine.getR(), danglingLine.getX()); + std::complex ysh2(g2, b2); std::complex zt = 1.0 / (ytr + ysh2); std::complex v0 = ytr * v1 / (ytr + ysh2); double v02 = std::abs(v0) * std::abs(v0); @@ -75,7 +75,7 @@ double DanglingLineData::getBoundaryBusU() const { } const std::string& DanglingLineData::getId() const { - return m_id; + return m_danglingLine.getId(); } double DanglingLineData::getTheta(const DanglingLine& danglingLine) { diff --git a/src/iidm/util/LinkData.cpp b/src/iidm/util/LinkData.cpp index 0ac9a2d35..729d01b8a 100644 --- a/src/iidm/util/LinkData.cpp +++ b/src/iidm/util/LinkData.cpp @@ -68,7 +68,7 @@ Flow flowBothEnds(const std::complex& y11, const std::complex& y std::complex flowYshunt(const std::complex& ysh, double u, double theta) { std::complex v = math::ComplexUtils::polar2Complex(u, theta); - return std::conj(ysh) * std::conj(v) *v; + return std::conj(ysh) * std::conj(v) * v; } double getFixedX(double x, double epsilonX, bool applyReactanceCorrection) { @@ -76,9 +76,7 @@ double getFixedX(double x, double epsilonX, bool applyReactanceCorrection) { } double getPhaseAngleClockDegrees(int phaseAngleClock) { - double phaseAngleClockDegree = 0.0; - phaseAngleClockDegree += phaseAngleClock * 30.0; - phaseAngleClockDegree = std::remainder(phaseAngleClockDegree, 360.0); + double phaseAngleClockDegree = std::remainder(phaseAngleClock * 30.0, 360.0); if (phaseAngleClockDegree > 180.0) { phaseAngleClockDegree -= 360.0; } @@ -139,8 +137,8 @@ BranchAdmittanceMatrix kronChain(const BranchAdmittanceMatrix& firstAdm, const B } admittance.y11 = yFirst11 - (yFirst1C * yFirstC1 / (yFirstCC + ySecondCC)); - admittance.y12 = yFirst1C * ySecondC2 / -(yFirstCC + ySecondCC); - admittance.y21 = ySecond2C * yFirstC1 / -(yFirstCC + ySecondCC); + admittance.y12 = - yFirst1C * ySecondC2 / (yFirstCC + ySecondCC); + admittance.y21 = - ySecond2C * yFirstC1 / (yFirstCC + ySecondCC); admittance.y22 = ySecond22 - (ySecond2C * ySecondC2 / (yFirstCC + ySecondCC)); return admittance; diff --git a/src/iidm/util/TwtData.cpp b/src/iidm/util/TwtData.cpp index 5e1254731..8003204f5 100644 --- a/src/iidm/util/TwtData.cpp +++ b/src/iidm/util/TwtData.cpp @@ -210,7 +210,7 @@ TwtData::TwtData(const ThreeWindingsTransformer& twt, int phaseAngleClock2, int } double TwtData::alpha(const ThreeWindingsTransformer::Leg& leg) { - return leg.getOptionalPhaseTapChanger() ? leg.getPhaseTapChanger().getCurrentStep().getAlpha() : 0.0; + return leg.getOptionalPhaseTapChanger() ? stdcxx::toRadians * leg.getPhaseTapChanger().getCurrentStep().getAlpha() : 0.0; } std::complex TwtData::calculateOneConnectedLegFlow(double u, double theta, @@ -313,7 +313,7 @@ std::complex TwtData::calculateTwoConnectedLegsStarBusVoltage(double u1, std::complex v2 = math::ComplexUtils::polar2Complex(u2, theta2); std::complex yshO = LinkData::kronAntenna(admittanceMatrixOpenLeg.y11, admittanceMatrixOpenLeg.y12, admittanceMatrixOpenLeg.y21, admittanceMatrixOpenLeg.y22, true); - return -(admittanceMatrixLeg1.y21 * v1 + (admittanceMatrixLeg2.y21 * v2)) + return -(admittanceMatrixLeg1.y21 * v1 + admittanceMatrixLeg2.y21 * v2) / (admittanceMatrixLeg1.y22 + admittanceMatrixLeg2.y22 + yshO); }