Merge branch 'e2nIEE:develop' into develop

CHANGELOG.rst

@@ -3,6 +3,7 @@ Change Log
 
 [upcoming release] - 2023-..-..
 -------------------------------
+- [ADDED] function to search std_types from the basic standard type library
 - [ADDED] Documentation for running powerflow using power-grid-model
 - [ADDED] exporting to :code:`GeoJSON` with all properties from :code:`bus`, :code:`res_bus` and :code:`line`, :code:`res_line`
 - [ADDED] function to run powerflow using the power-grid-model library
@@ -34,7 +35,12 @@ Change Log
 - [ADDED] option to use a second tap changer for the trafo element
 - [CHANGED] parameters of function merge_internal_net_and_equivalent_external_net()
 - [FIXED] :code:`convert_format.py`: update the attributes of the characteristic objects to match the new characteristic
+- [FIXED] fixed the wrong id numbers for pypower powerflow algorithms fdxb and fdbx
+- [FIXED] additional arguments from mpc saved to net._options: create "_options" if it does not exist
 - [CHANGED] cim2pp: extracted getting default classes, added generic setting datatypes from CGMES XMI schema
+- [ADDED] function :code:`getOTDF` to obtain Outage Transfer Distribution Factors, that can be used to analyse outages using the DC approximation of the power system
+- [ADDED] function :code:`outage_results_OTDF` to obtain the matrix of results for all outage scenarios, with rows as outage scenarios and columns as branch power flows in that scenario
+- [FIXED] add some safeguards for TDPF to avoid numerical issues in some cases
 
 
 [2.13.1] - 2023-05-12

doc/plotting/matplotlib/simple_plot.rst

@@ -2,7 +2,7 @@
 Simple Plotting
 =============================
 
-The function simple_plot() can be used for simple plotting. For advanced possibilities see the tutorials
+The function simple_plot() can be used for simple plotting. For advanced possibilities see the `tutorial <http://nbviewer.jupyter.org/github/e2nIEE/pandapower/blob/develop/tutorials/plotting_basic.ipynb>`_.
 
 .. _simple_plot:
 

pandapower/auxiliary.py

@@ -318,9 +318,9 @@ def __repr__(self):  # pragma: no cover
 
 def plural_s(number):
     if number > 1:
-        return ""
-    else:
         return "s"
+    else:
+        return ""
 
 def _preserve_dtypes(df, dtypes):
     for item, dtype in list(dtypes.items()):

pandapower/control/controller/trafo/ContinuousTapControl.py

@@ -4,6 +4,7 @@
 # and Energy System Technology (IEE), Kassel. All rights reserved.
 
 import numpy as np
+
 from pandapower.auxiliary import read_from_net, write_to_net
 from pandapower.control.controller.trafo_control import TrafoController
 
@@ -97,6 +98,8 @@ def is_converged(self, net):
             return True
 
         vm_pu = read_from_net(net, "res_bus", self.controlled_bus, "vm_pu", self._read_write_flag)
+        # this is possible in case the trafo is set out of service by the connectivity check
+        is_nan =  np.isnan(vm_pu)
         self.tap_pos = read_from_net(net, self.trafotable, self.controlled_tid, "tap_pos", self._read_write_flag)
         difference = 1 - self.vm_set_pu / vm_pu
 
@@ -110,4 +113,4 @@ def is_converged(self, net):
         else:
             converged = np.abs(difference) < self.tol
 
-        return np.all(converged)
+        return np.all(np.logical_or(converged, is_nan))

pandapower/control/controller/trafo/DiscreteTapControl.py

@@ -3,9 +3,11 @@
 # Copyright (c) 2016-2023 by University of Kassel and Fraunhofer Institute for Energy Economics
 # and Energy System Technology (IEE), Kassel. All rights reserved.
 import numpy as np
+
 from pandapower.auxiliary import read_from_net, write_to_net
 from pandapower.control.controller.trafo_control import TrafoController
 
+
 class DiscreteTapControl(TrafoController):
     """
     Trafo Controller with local tap changer voltage control.
@@ -115,6 +117,8 @@ def is_converged(self, net):
             return True
 
         vm_pu = read_from_net(net, "res_bus", self.controlled_bus, "vm_pu", self._read_write_flag)
+        # this is possible in case the trafo is set out of service by the connectivity check
+        is_nan = np.isnan(vm_pu)
         self.tap_pos = read_from_net(net, self.trafotable, self.controlled_tid, "tap_pos", self._read_write_flag)
 
         reached_limit = np.where(self.tap_side_coeff * self.tap_sign == 1,
@@ -125,5 +129,4 @@ def is_converged(self, net):
 
         converged = np.logical_or(reached_limit, np.logical_and(self.vm_lower_pu < vm_pu, vm_pu < self.vm_upper_pu))
 
-        return np.all(converged)
-
+        return np.all(np.logical_or(converged, is_nan))

pandapower/converter/matpower/from_mpc.py

@@ -66,6 +66,8 @@ def from_mpc(mpc_file, f_hz=50, casename_mpc_file='mpc', validate_conversion=Fal
         ppc = _m2ppc(mpc_file, casename_mpc_file)
     net = from_ppc(ppc, f_hz=f_hz, validate_conversion=validate_conversion, **kwargs)
     if "mpc_additional_data" in ppc:
+        if "_options" not in net:
+            net["_options"] = dict()
         net._options.update(ppc["mpc_additional_data"])
         logger.info('added fields %s in net._options' % list(ppc["mpc_additional_data"].keys()))
 

pandapower/pf/create_jacobian_tdpf.py

@@ -115,7 +115,7 @@ def calc_r_theta(t_air_pu, a0, a1, a2, i_square_pu, p_loss_pu):
     2019, pp. 1-6, doi: 10.1109/EEEIC.2019.8783234.
     """
     t_rise_pu = a0 + a1 * i_square_pu + a2 * np.square(i_square_pu) - t_air_pu
-    r_theta_pu = t_rise_pu / p_loss_pu
+    r_theta_pu = t_rise_pu / np.where(p_loss_pu == 0, 1e-6, p_loss_pu)
     return r_theta_pu
 
 

pandapower/pf/runpf_pypower.py

@@ -78,7 +78,7 @@ def _get_options(options, **kwargs):
     max_iteration = options["max_iteration"]
 
     # algorithms implemented within pypower
-    algorithm_pypower_dict = {'nr': 1, 'fdbx': 2, 'fdxb': 3, 'gs': 4}
+    algorithm_pypower_dict = {'nr': 1, 'fdxb': 2, 'fdbx': 3, 'gs': 4}
 
     ppopt = ppoption(ENFORCE_Q_LIMS=enforce_q_lims, PF_TOL=tolerance_mva,
                      PF_ALG=algorithm_pypower_dict[algorithm], **kwargs)

pandapower/plotting/generic_geodata.py

@@ -142,29 +142,32 @@ def coords_from_igraph(graph, roots, meshed=False, calculate_meshed=False):
     return list(zip(*layout.coords))
 
 
-def coords_from_nxgraph(mg=None):
+def coords_from_nxgraph(mg=None, layout_engine='neato'):
     """
     Create a list of generic coordinates from a networkx graph layout.
 
     :param mg: The networkx graph on which the coordinates shall be based
     :type mg: networkx.Graph
+    :param layout_engine: GraphViz Layout Engine for layouting a network. See https://graphviz.org/docs/layouts/
+    :type layout_engine: str
     :return: coords - list of coordinates from the graph layout
     """
     # workaround for bug in agraph
-    for u, v in mg.edges(data=False, keys=False):
+    for u, v in mg.edges(data=False):
         if 'key' in mg[int(u)][int(v)]:
             del mg[int(u)][int(v)]['key']
         if 'key' in mg[int(u)][int(v)][0]:
             del mg[int(u)][int(v)][0]['key']
     # ToDo: Insert fallback layout for nxgraph
-    return list(zip(*(list(nx.drawing.nx_agraph.graphviz_layout(mg, prog='neato').values()))))
+    return list(zip(*(list(nx.drawing.nx_agraph.graphviz_layout(mg, prog=layout_engine).values()))))
 
 
 def create_generic_coordinates(net, mg=None, library="igraph",
                                respect_switches=False,
                                geodata_table="bus_geodata",
                                buses=None,
-                               overwrite=False):
+                               overwrite=False,
+                               layout_engine='neato'):
     """
     This function will add arbitrary geo-coordinates for all buses based on an analysis of branches
     and rings. It will remove out of service buses/lines from the net. The coordinates will be
@@ -174,6 +177,8 @@ def create_generic_coordinates(net, mg=None, library="igraph",
     :type net: pandapowerNet
     :param mg: Existing networkx multigraph, if available. Convenience to save computation time.
     :type mg: networkx.Graph
+    :param respect_switches: respect switches in a network for generic coordinates
+    :type respect_switches: bool
     :param library: "igraph" to use igraph package or "networkx" to use networkx package
     :type library: str
     :param geodata_table: table to write the generic geodatas to
@@ -182,6 +187,8 @@ def create_generic_coordinates(net, mg=None, library="igraph",
     :type buses: list
     :param overwrite: overwrite existing geodata
     :type overwrite: bool
+    :param layout_engine: GraphViz Layout Engine for layouting a network. See https://graphviz.org/docs/layouts/
+    :type layout_engine: str
     :return: net - pandapower network with added geo coordinates for the buses
 
     :Example:
@@ -200,7 +207,7 @@ def create_generic_coordinates(net, mg=None, library="igraph",
                                      include_out_of_service=True)
         else:
             nxg = copy.deepcopy(mg)
-        coords = coords_from_nxgraph(nxg)
+        coords = coords_from_nxgraph(nxg, layout_engine=layout_engine)
     else:
         raise ValueError("Unknown library %s - chose 'igraph' or 'networkx'" % library)
     if len(coords):

pandapower/plotting/plotly/traces.py

@@ -1094,15 +1094,25 @@ def draw_traces(traces, on_map=False, map_style='basic', showlegend=True, figsiz
             for trace in traces:
                 xs += trace.get('x') or trace['lon']
                 ys += trace.get('y') or trace['lat']
-            x_dropna = pd.Series(xs).dropna()
-            y_dropna = pd.Series(ys).dropna()
-            xrange = x_dropna.max() - x_dropna.min()
-            yrange = y_dropna.max() - y_dropna.min()
-            ratio = xrange / yrange
-            if ratio < 1:
-                aspectratio = (ratio, 1.)
+            xs_arr = np.array(xs)
+            ys_arr = np.array(ys)
+            xrange = np.nanmax(xs_arr) - np.nanmin(xs_arr)
+            yrange = np.nanmax(ys_arr) - np.nanmin(ys_arr)
+
+            # the ratio only makes sense, if xrange and yrange != 0
+            if xrange == 0 and yrange == 0:
+                aspectratio = (1, 1)
+            elif xrange == 0:
+                aspectratio = (0.35, 1)
+            elif yrange == 0:
+                aspectratio = (1, 0.35)
+
             else:
-                aspectratio = (1., 1 / ratio)
+                ratio = xrange / yrange
+                if ratio < 1:
+                    aspectratio = (ratio, 1.)
+                else:
+                    aspectratio = (1., 1 / ratio)
 
         aspectratio = np.array(aspectratio) / max(aspectratio)
         fig['layout']['width'], fig['layout']['height'] = ([ar * figsize * 700 for ar in aspectratio])

pandapower/pypower/makeLODF.py

@@ -66,3 +66,118 @@ def makeLODF(branch, PTDF):
     update_LODF_diag(LODF)
 
     return LODF
+
+
+def makeOTDF(PTDF, LODF, outage_branches):
+    """
+    Compute the Outage Transfer Distribution Factors (OTDF) matrix.
+
+    This function creates the OTDF matrix that relates bus power injections
+    to branch flows for specified outage scenarios. It's essential that
+    outage branches do not lead to isolated nodes or disconnected islands
+    in the grid.
+
+    The grid cannot have isolated nodes or disconnected islands. Use the
+    pandapower.topology module to identify branches that, if outaged, would
+    lead to isolated nodes (determine_stubs) or islands (find_graph_characteristics).
+
+    The resulting matrix has a width equal to the number of nodes and a length
+    equal to the number of outage branches multiplied by the total number of branches.
+    The dot product of OTDF and the bus power vector in generation reference frame
+    (positive for generation, negative for consumption - the opposite of res_bus.p_mw)
+    yields an array with outage branch power flows for every outage scenario,
+    facilitating the analysis of all outage scenarios under a DC
+    power flow approximation.
+
+    Parameters
+    ----------
+    PTDF : numpy.ndarray
+        The Power Transfer Distribution Factor matrix, defining the sensitivity
+        of branch flows to bus power injections.
+    LODF : numpy.ndarray
+        The Line Outage Distribution Factor matrix, describing how branch flows
+        are affected by outages of other branches.
+    outage_branches : list or numpy.ndarray
+        Indices of branches for which outage scenarios are to be considered.
+
+    Returns
+    -------
+    OTDF : numpy.ndarray
+        The Outage Transfer Distribution Factor matrix. Rows correspond to
+        outage scenarios, and columns correspond to branch flows.
+
+    Examples
+    --------
+    >>> H = makePTDF(baseMVA, bus, branch)
+    >>> LODF = makeLODF(branch, H)
+    >>> outage_branches = [0, 2]  # Example branch indices for outage scenarios
+    >>> OTDF = makeOTDF(H, LODF, outage_branches)
+    >>> # To obtain a 2D array with the outage results:
+    >>> outage_results = (OTDF @ Pbus).reshape(len(outage_branches), -1)
+
+    Notes
+    -----
+    - The function assumes a DC power flow model.
+    - Ensure that the specified outage branches do not lead to grid
+      disconnection or isolated nodes.
+    """
+    OTDF = np.vstack([PTDF + LODF[:, [i]] @ PTDF[[i], :] for i in outage_branches])
+    return OTDF
+
+
+def outage_results_OTDF(OTDF, Pbus, outage_branches):
+    """
+    Calculate the branch power flows for each outage scenario based on the given
+    Outage Transfer Distribution Factors (OTDF), bus power injections (Pbus), and
+    specified outage branches.
+
+    This function computes how branch flows are affected under N-1 contingency
+    scenarios (i.e., for each branch outage specified). It uses the OTDF matrix and
+    the bus power vector (Pbus) to determine the branch flows in each outage scenario.
+
+    Pbus should represent the net power at each bus in the generation reference case.
+
+    Parameters
+    ----------
+    OTDF : numpy.ndarray
+        The Outage Transfer Distribution Factor matrix, which relates bus power
+        injections to branch flows under specific outage scenarios. Its shape
+        should be (num_outage_scenarios * num_branches, num_buses).
+    Pbus : numpy.ndarray
+        A vector representing the net power injections at each bus. Positive values 
+        for generation, negative for consumption. Its length should be equal to 
+        the total number of buses.
+    outage_branches : numpy.ndarray
+        An array of indices representing the branches that are outaged in each
+        scenario. Its length should be equal to the number of outage scenarios.
+
+    Returns
+    -------
+    numpy.ndarray
+        A 2D array where each row corresponds to an outage scenario and each column
+        represents the resulting power flow in a branch. The number of rows is equal
+        to the number of outage scenarios, and the number of columns is equal to the
+        number of branches.
+
+    Examples
+    --------
+    >>> OTDF = np.array([...])  # example OTDF matrix
+    >>> Pbus = np.array([...])  # example bus power vector
+    >>> outage_branches = np.array([...])  # example outage branches
+    >>> branch_flows = outage_results_OTDF(OTDF,Pbus,outage_branches)
+
+    Notes
+    -----
+    The function assumes a linear relationship between bus power injections and
+    branch flows, which is typical in DC power flow models.
+    """
+    # get branch flows as an array first:
+    nminus1_otdf = (OTDF @ Pbus.reshape(-1, 1))
+    # reshape to a 2D array with rows relating to outage scenarios and columns to
+    # the resulting branch power flows
+    nminus1_otdf = nminus1_otdf.reshape(outage_branches.shape[0], -1)
+    return nminus1_otdf
+
+
+
+

pandapower/pypower/newtonpf.py

@@ -298,7 +298,8 @@ def newtonpf(Ybus, Sbus, V0, ref, pv, pq, ppci, options, makeYbus=None):
 
         if tdpf:
             # update the R, g, b for the tdpf_lines, and the Y-matrices
-            branch[tdpf_lines, BR_R] = r = r_ref_pu * (1 + alpha_pu * (T - t_ref_pu))
+            # here: limit the change of the R to reflect a realistic range of values for T to avoid numerical issues
+            branch[tdpf_lines, BR_R] = r = r_ref_pu * (1 + alpha_pu * np.clip(np.nan_to_num(T - t_ref_pu), -50, 250 / T_base))
             # todo expansion with SSC and VSC (that are not controllable)
             Ybus, Yf, Yt = makeYbus(baseMVA, bus, branch)
             g, b = calc_g_b(r, x)