Add DanglingLineData (and unit tests) and update DanglingLineBoundary…

… implementation to use it

Signed-off-by: Sébastien LAIGRE <[email protected]>

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;
 };

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 <string>
+
+#include <powsybl/stdcxx/math.hpp>
+
+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

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 <complex>
+
+namespace powsybl {
+
+namespace math {
+
+namespace ComplexUtils {
+
+std::complex<double> polar2Complex(double r, double theta);
+
+}  // namespace ComplexUtils
+
+}  // namespace math
+
+}  // namespace powsybl
+
+#endif  // POWSYBL_MATH_COMPLEXUTILS_HPP

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
 

src/iidm/util/DanglingLineBoundary.cpp

@@ -10,6 +10,7 @@
 #include <powsybl/iidm/Bus.hpp>
 #include <powsybl/iidm/DanglingLine.hpp>
 #include <powsybl/iidm/Terminal.hpp>
+#include <powsybl/iidm/util/DanglingLineData.hpp>
 #include <powsybl/iidm/util/SV.hpp>
 #include <powsybl/stdcxx/math.hpp>
 
@@ -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<Bus>& 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<Branch::Side> 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

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 <powsybl/iidm/util/DanglingLineData.hpp>
+
+#include <complex>
+
+#include <powsybl/iidm/Bus.hpp>
+#include <powsybl/iidm/DanglingLine.hpp>
+#include <powsybl/iidm/util/LinkData.hpp>
+#include <powsybl/math/ComplexUtils.hpp>
+
+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<double>& v1 = math::ComplexUtils::polar2Complex(m_u1, m_theta1);
+
+    std::complex<double> 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<double>(m_g1, m_b1), std::complex<double>(m_g2, m_b2));
+        vBoundaryBus = (- adm.y21 * v1) / adm.y22;
+    } else {
+        // Two buses Loadflow
+        std::complex<double> sBoundary(-m_p0, -m_q0);
+        std::complex<double> ytr = 1.0 / std::complex<double>(m_r, m_x);
+        std::complex<double> ysh2(m_g2, m_b2);
+        std::complex<double> zt = 1.0 / (ytr + ysh2);
+        std::complex<double> v0 = ytr * v1 / (ytr + ysh2);
+        double v02 = std::abs(v0) * std::abs(v0);
+
+        std::complex<double> 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<double>(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

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 <powsybl/math/ComplexUtils.hpp>
+
+#include <powsybl/PowsyblException.hpp>
+#include <powsybl/stdcxx/format.hpp>
+
+namespace powsybl {
+
+namespace math {
+
+namespace ComplexUtils {
+
+std::complex<double> polar2Complex(double r, double theta) {
+    if (r < 0.0) {
+        throw PowsyblException(stdcxx::format("negative complex module %1%", r));
+    }
+    return std::complex<double>(r * std::cos(theta), r * std::sin(theta));
+}
+
+}  // namespace ComplexUtils
+
+}  // namespace math
+
+}  // namespace powsybl

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
 )

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());