Merge pull request #562 from PowerGridModel/feature/transformer-tap-r…

…anking-part-1

Transformer ranking part 1

power_grid_model_c/power_grid_model/include/power_grid_model/common/exception.hpp

@@ -118,6 +118,14 @@ class InvalidRegulatedObject : public PowerGridError {
     }
 };
 
+class AutomaticTapCalculationError : public PowerGridError {
+  public:
+    AutomaticTapCalculationError(ID id) {
+        append_msg("Automatic tap changing regulator with tap_side at LV side is not supported. Found at id" +
+                   std::to_string(id)); // NOSONAR
+    }
+};
+
 class IDWrongType : public PowerGridError {
   public:
     explicit IDWrongType(ID id) { append_msg("Wrong type for object with id " + std::to_string(id) + '\n'); }

power_grid_model_c/power_grid_model/include/power_grid_model/container.hpp

@@ -137,7 +137,9 @@ class Container<RetrievableTypes<GettableTypes...>, StorageableTypes...> {
     }
 
     // get size
-    template <class Gettable> Idx size() const {
+    template <class Gettable>
+        requires(std::same_as<Gettable, GettableTypes> || ...)
+    Idx size() const {
         assert(construction_complete_);
         return size_[get_cls_pos_v<Gettable, GettableTypes...>];
     }

power_grid_model_c/power_grid_model/include/power_grid_model/main_core/input.hpp

@@ -11,6 +11,9 @@
 
 namespace power_grid_model::main_core {
 
+constexpr std::array<Branch3Side, 3> const branch3_sides = {Branch3Side::side_1, Branch3Side::side_2,
+                                                            Branch3Side::side_3};
+
 // template to construct components
 // using forward interators
 // different selection based on component type
@@ -124,6 +127,31 @@ inline void add_component(MainModelState<ComponentContainer>& state, ForwardIter
                 }
             }();
 
+            if (regulated_object_idx.group == get_component_type_index<Transformer>(state)) {
+                auto const& regulated_object = get_component<Transformer>(state, regulated_object_idx);
+
+                auto const non_tap_side =
+                    regulated_object.tap_side() == BranchSide::from ? BranchSide::to : BranchSide::from;
+                if (get_component<Node>(state, regulated_object.node(regulated_object.tap_side())).u_rated() <
+                    get_component<Node>(state, regulated_object.node(non_tap_side)).u_rated()) {
+                    throw AutomaticTapCalculationError(id);
+                }
+            } else if (regulated_object_idx.group == get_component_type_index<ThreeWindingTransformer>(state)) {
+                auto const& regulated_object = get_component<ThreeWindingTransformer>(state, regulated_object_idx);
+                auto const tap_side_u_rated =
+                    get_component<Node>(state, regulated_object.node(regulated_object.tap_side())).u_rated();
+                for (auto const side : branch3_sides) {
+                    if (side == regulated_object.tap_side()) {
+                        continue;
+                    }
+                    if (tap_side_u_rated < get_component<Node>(state, regulated_object.node(side)).u_rated()) {
+                        throw AutomaticTapCalculationError(id);
+                    }
+                }
+            } else {
+                throw InvalidRegulatedObject(input.regulated_object, Component::name);
+            }
+
             auto const regulated_object_type = get_component<Base>(state, regulated_object_idx).math_model_type();
             double const u_rated = get_component<Node>(state, regulated_terminal).u_rated();
 

power_grid_model_c/power_grid_model/include/power_grid_model/optimizer/tap_position_optimizer.hpp

@@ -6,10 +6,14 @@
 
 #include "base_optimizer.hpp"
 
+#include "../all_components.hpp"
 #include "../auxiliary/dataset.hpp"
 #include "../common/enum.hpp"
+#include "../common/exception.hpp"
+#include "../main_core/state_queries.hpp"
 
 #include <boost/graph/compressed_sparse_row_graph.hpp>
+
 #include <functional>
 #include <queue>
 #include <vector>
@@ -21,114 +25,239 @@ namespace detail = power_grid_model::optimizer::detail;
 
 using TrafoGraphIdx = Idx;
 using EdgeWeight = int64_t;
-using WeightedTrafo = std::pair<Idx2D, EdgeWeight>;
-using WeightedTrafoList = std::vector<WeightedTrafo>;
-using RankedTransformerGroups = std::vector<std::vector<Idx2D>>;
 constexpr auto infty = std::numeric_limits<Idx>::max();
+constexpr Idx2D unregulated_idx = {-1, -1};
 
 struct TrafoGraphVertex {
-    bool is_source{}; // is_source = true if the vertex is a source
+    bool is_source{};
 };
 
 struct TrafoGraphEdge {
-    Idx2D pos{};
+    Idx2D regulated_idx{};
     EdgeWeight weight{};
+
+    bool operator==(const TrafoGraphEdge& other) const {
+        return regulated_idx == other.regulated_idx && weight == other.weight;
+    } // thanks boost
+
+    auto operator<=>(const TrafoGraphEdge& other) const {
+        if (auto cmp = weight <=> other.weight; cmp != 0) { // NOLINT(modernize-use-nullptr)
+            return cmp;
+        }
+        if (auto cmp = regulated_idx.group <=> other.regulated_idx.group; cmp != 0) { // NOLINT(modernize-use-nullptr)
+            return cmp;
+        }
+        return regulated_idx.pos <=> other.regulated_idx.pos;
+    }
+};
+
+using TrafoGraphEdges = std::vector<std::pair<TrafoGraphIdx, TrafoGraphIdx>>;
+using TrafoGraphEdgeProperties = std::vector<TrafoGraphEdge>;
+using RankedTransformerGroups = std::vector<Idx2D>;
+
+struct RegulatedObjects {
+    std::set<Idx> transformers{};
+    std::set<Idx> transformers3w{};
 };
 
-// TODO(mgovers): investigate whether this really is the correct graph structure
 using TransformerGraph = boost::compressed_sparse_row_graph<boost::directedS, TrafoGraphVertex, TrafoGraphEdge,
                                                             boost::no_property, TrafoGraphIdx, TrafoGraphIdx>;
 
+inline void add_to_edge(TrafoGraphEdges& edges, TrafoGraphEdgeProperties& edge_props, Idx const& start, Idx const& end,
+                        TrafoGraphEdge const& edge_prop) {
+    edges.emplace_back(start, end);
+    edge_props.emplace_back(edge_prop);
+}
+
+inline void process_trafo3w_edge(ThreeWindingTransformer const& transformer3w, bool const& trafo3w_is_regulated,
+                                 Idx2D const& trafo3w_idx, TrafoGraphEdges& edges,
+                                 TrafoGraphEdgeProperties& edge_props) {
+    using enum Branch3Side;
+
+    constexpr std::array<std::tuple<Branch3Side, Branch3Side>, 3> const branch3_combinations{
+        {{side_1, side_2}, {side_2, side_3}, {side_3, side_1}}};
+
+    for (auto const& [first_side, second_side] : branch3_combinations) {
+        if (!transformer3w.status(first_side) || !transformer3w.status(second_side)) {
+            continue;
+        }
+        auto const& from_node = transformer3w.node(first_side);
+        auto const& to_node = transformer3w.node(second_side);
+
+        auto const tap_at_first_side = transformer3w.tap_side() == first_side;
+        auto const single_direction_condition =
+            trafo3w_is_regulated && (tap_at_first_side || transformer3w.tap_side() == second_side);
+        // ranking
+        if (single_direction_condition) {
+            auto const& tap_side_node = tap_at_first_side ? from_node : to_node;
+            auto const& non_tap_side_node = tap_at_first_side ? to_node : from_node;
+            // add regulated idx only when the first side node is tap side node.
+            // This is done to add only one directional edge with regulated idx.
+            Idx2D const regulated_idx = from_node == tap_side_node ? unregulated_idx : trafo3w_idx;
+            add_to_edge(edges, edge_props, tap_side_node, non_tap_side_node, {regulated_idx, 1});
+        } else {
+            add_to_edge(edges, edge_props, from_node, to_node, {unregulated_idx, 1});
+            add_to_edge(edges, edge_props, to_node, from_node, {unregulated_idx, 1});
+        }
+    }
+}
+
+template <std::derived_from<ThreeWindingTransformer> Component, class ComponentContainer>
+    requires main_core::model_component_state_c<main_core::MainModelState, ComponentContainer, Component>
+constexpr void add_edge(main_core::MainModelState<ComponentContainer> const& state,
+                        RegulatedObjects const& regulated_objects, TrafoGraphEdges& edges,
+                        TrafoGraphEdgeProperties& edge_props) {
+
+    for (auto const& transformer3w : state.components.template citer<ThreeWindingTransformer>()) {
+        bool const trafo3w_is_regulated = regulated_objects.transformers3w.contains(transformer3w.id());
+        Idx2D const trafo3w_idx = main_core::get_component_idx_by_id(state, transformer3w.id());
+        process_trafo3w_edge(transformer3w, trafo3w_is_regulated, trafo3w_idx, edges, edge_props);
+    }
+}
+
+template <std::derived_from<Transformer> Component, class ComponentContainer>
+    requires main_core::model_component_state_c<main_core::MainModelState, ComponentContainer, Component>
+constexpr void add_edge(main_core::MainModelState<ComponentContainer> const& state,
+                        RegulatedObjects const& regulated_objects, TrafoGraphEdges& edges,
+                        TrafoGraphEdgeProperties& edge_props) {
+    for (auto const& transformer : state.components.template citer<Transformer>()) {
+        if (!transformer.from_status() || !transformer.to_status()) {
+            continue;
+        }
+        auto const& from_node = transformer.from_node();
+        auto const& to_node = transformer.to_node();
+
+        if (regulated_objects.transformers.contains(transformer.id())) {
+            auto const tap_at_from_side = transformer.tap_side() == BranchSide::from;
+            auto const& tap_side_node = tap_at_from_side ? from_node : to_node;
+            auto const& non_tap_side_node = tap_at_from_side ? to_node : from_node;
+            add_to_edge(edges, edge_props, tap_side_node, non_tap_side_node,
+                        {main_core::get_component_idx_by_id(state, transformer.id()), 1});
+        } else {
+            add_to_edge(edges, edge_props, from_node, to_node, {unregulated_idx, 1});
+            add_to_edge(edges, edge_props, to_node, from_node, {unregulated_idx, 1});
+        }
+    }
+}
+
+template <std::derived_from<Branch> Component, class ComponentContainer>
+    requires main_core::model_component_state_c<main_core::MainModelState, ComponentContainer, Component> &&
+             (!transformer_c<Component>)
+constexpr void add_edge(main_core::MainModelState<ComponentContainer> const& state,
+                        RegulatedObjects const& /* regulated_objects */, TrafoGraphEdges& edges,
+                        TrafoGraphEdgeProperties& edge_props) {
+    auto const& iter = state.components.template citer<Component>();
+    edges.reserve(std::distance(iter.begin(), iter.end()) * 2);
+    edge_props.reserve(std::distance(iter.begin(), iter.end()) * 2);
+    for (auto const& branch : iter) {
+        if (!branch.from_status() || !branch.to_status()) {
+            continue;
+        }
+        add_to_edge(edges, edge_props, branch.from_node(), branch.to_node(), {unregulated_idx, 0});
+        add_to_edge(edges, edge_props, branch.to_node(), branch.from_node(), {unregulated_idx, 0});
+    }
+}
+
+template <typename... ComponentTypes, main_core::main_model_state_c State>
+inline auto add_edges(State const& state, RegulatedObjects const& regulated_objects, TrafoGraphEdges& edges,
+                      TrafoGraphEdgeProperties& edge_props) {
+    (add_edge<ComponentTypes>(state, regulated_objects, edges, edge_props), ...);
+}
+
+template <main_core::main_model_state_c State>
+inline auto retrieve_regulator_info(State const& state) -> RegulatedObjects {
+    RegulatedObjects regulated_objects;
+    for (auto const& regulator : state.components.template citer<TransformerTapRegulator>()) {
+        if (!regulator.status()) {
+            continue;
+        }
+        if (regulator.regulated_object_type() == ComponentType::branch) {
+            regulated_objects.transformers.emplace(regulator.regulated_object());
+        } else {
+            regulated_objects.transformers3w.emplace(regulator.regulated_object());
+        }
+    }
+    return regulated_objects;
+}
+
 template <main_core::main_model_state_c State>
-inline auto build_transformer_graph(State const& /*state*/) -> TransformerGraph {
-    // TODO(nbharambe): implement
-    return {};
+inline auto build_transformer_graph(State const& state) -> TransformerGraph {
+    TrafoGraphEdges edges;
+    TrafoGraphEdgeProperties edge_props;
+
+    const RegulatedObjects regulated_objects = retrieve_regulator_info(state);
+
+    add_edges<Transformer, ThreeWindingTransformer, Line, Link>(state, regulated_objects, edges, edge_props);
+
+    // build graph
+    TransformerGraph trafo_graph{boost::edges_are_unsorted_multi_pass, edges.cbegin(), edges.cend(),
+                                 edge_props.cbegin(),
+                                 static_cast<TrafoGraphIdx>(state.components.template size<Node>())};
+
+    BGL_FORALL_VERTICES(v, trafo_graph, TransformerGraph) { trafo_graph[v].is_source = false; }
+
+    // Mark sources
+    for (auto const& source : state.components.template citer<Source>()) {
+        // ignore disabled sources
+        trafo_graph[source.node()].is_source = source.status();
+    }
+
+    return trafo_graph;
 }
 
-inline auto process_edges_dijkstra(Idx v, std::vector<EdgeWeight>& edge_weight, std::vector<Idx2D>& edge_pos,
-                                   TransformerGraph const& graph) -> void {
+inline void process_edges_dijkstra(Idx v, std::vector<EdgeWeight>& vertex_distances, TransformerGraph const& graph) {
     using TrafoGraphElement = std::pair<EdgeWeight, TrafoGraphIdx>;
     std::priority_queue<TrafoGraphElement, std::vector<TrafoGraphElement>, std::greater<>> pq;
-    edge_weight[v] = 0;
-    edge_pos[v] = {v, v};
+    vertex_distances[v] = 0;
     pq.push({0, v});
 
     while (!pq.empty()) {
         auto [dist, u] = pq.top();
         pq.pop();
 
-        if (dist != edge_weight[u]) {
+        if (dist != vertex_distances[u]) {
             continue;
         }
 
-        for (auto e : boost::make_iterator_range(boost::out_edges(u, graph))) {
+        BGL_FORALL_OUTEDGES(u, e, graph, TransformerGraph) {
             auto v = boost::target(e, graph);
             const EdgeWeight weight = graph[e].weight;
 
-            if (edge_weight[u] + weight < edge_weight[v]) {
-                edge_weight[v] = edge_weight[u] + weight;
-                edge_pos[v] = graph[e].pos;
-                pq.push({edge_weight[v], v});
+            if (vertex_distances[u] + weight < vertex_distances[v]) {
+                vertex_distances[v] = vertex_distances[u] + weight;
+                pq.push({vertex_distances[v], v});
             }
         }
     }
 }
 
-// Step 2: Initialize the rank of all vertices (transformer nodes) as infinite (INT_MAX)
-// Step 3: Loop all the connected sources (status == 1)
-//      a. Perform Dijkstra shortest path algorithm from the vertex with that source.
-//         This is to determine the shortest path of all vertices to this particular source.
-inline auto get_edge_weights(TransformerGraph const& graph) -> WeightedTrafoList {
-    std::vector<EdgeWeight> edge_weight(boost::num_vertices(graph), infty);
-    std::vector<Idx2D> edge_pos(boost::num_vertices(graph));
-
-    for (auto v : boost::make_iterator_range(boost::vertices(graph))) {
+inline auto get_edge_weights(TransformerGraph const& graph) -> TrafoGraphEdgeProperties {
+    std::vector<EdgeWeight> vertex_distances(boost::num_vertices(graph), infty);
+    BGL_FORALL_VERTICES(v, graph, TransformerGraph) {
         if (graph[v].is_source) {
-            process_edges_dijkstra(v, edge_weight, edge_pos, graph);
+            process_edges_dijkstra(v, vertex_distances, graph);
         }
     }
 
-    WeightedTrafoList result;
-    for (size_t i = 0; i < edge_weight.size(); ++i) {
-        result.emplace_back(edge_pos[i], edge_weight[i]);
+    TrafoGraphEdgeProperties result;
+    BGL_FORALL_EDGES(e, graph, TransformerGraph) {
+        if (graph[e].regulated_idx == unregulated_idx) {
+            continue;
+        }
+        result.push_back({graph[e].regulated_idx, vertex_distances[boost::source(e, graph)]});
     }
 
     return result;
 }
 
-// Step 4: Loop all transformers with automatic tap changers, including the transformers which are not
-//         fully connected
-//   a.Rank of the transformer <-
-//       i. Infinity(INT_MAX), if tap side of the transformer is disconnected.
-//          The transformer regulation should be ignored
-//       ii.Rank of the vertex at the tap side of the transformer, if tap side of the transformer is connected
-inline auto transformer_disconnected(Idx2D const& /*pos*/) -> bool {
-    // <TODO: jguo> waiting for the functionalities in step 1 to be implemented
-    return false;
-}
-
-inline auto rank_transformers(WeightedTrafoList const& w_trafo_list) -> RankedTransformerGroups {
+inline auto rank_transformers(TrafoGraphEdgeProperties const& w_trafo_list) -> RankedTransformerGroups {
     auto sorted_trafos = w_trafo_list;
 
-    for (auto& trafo : sorted_trafos) {
-        if (transformer_disconnected(trafo.first)) {
-            trafo.second = infty;
-        }
-    }
-
     std::sort(sorted_trafos.begin(), sorted_trafos.end(),
-              [](const WeightedTrafo& a, const WeightedTrafo& b) { return a.second < b.second; });
-
-    RankedTransformerGroups groups;
-    Idx last_weight = -1;
-    for (const auto& trafo : sorted_trafos) {
-        if (groups.empty() || last_weight != trafo.second) {
-            groups.push_back(std::vector<Idx2D>{trafo.first});
-            last_weight = trafo.second;
-        } else {
-            groups.back().push_back(trafo.first);
-        }
-    }
+              [](const TrafoGraphEdge& a, const TrafoGraphEdge& b) { return a.weight < b.weight; });
+
+    RankedTransformerGroups groups(sorted_trafos.size());
+    std::ranges::transform(sorted_trafos, groups.begin(), [](const TrafoGraphEdge& x) { return x.regulated_idx; });
     return groups;
 }