Merge pull request #265 from openego/feature/power-flow-options

Implemented troubleshooting options while performing power flow

edisgo/edisgo.py

@@ -5,8 +5,10 @@
 import pickle
 import shutil
 
+from numbers import Number
 from pathlib import PurePath
 
+import numpy as np
 import pandas as pd
 
 from edisgo.flex_opt.charging_strategies import charging_strategy
@@ -662,7 +664,16 @@ def import_generators(self, generator_scenario=None, **kwargs):
             edisgo_object=self, generator_scenario=generator_scenario, **kwargs
         )
 
-    def analyze(self, mode=None, timesteps=None, raise_not_converged=True, **kwargs):
+    def analyze(
+        self,
+        mode: str | None = None,
+        timesteps: pd.Timestamp | pd.DatetimeIndex | None = None,
+        raise_not_converged: bool = True,
+        troubleshooting_mode: str | None = None,
+        range_start: Number = 0.1,
+        range_num: int = 10,
+        **kwargs,
+    ):
         """
         Conducts a static, non-linear power flow analysis.
 
@@ -720,6 +731,30 @@ def analyze(self, mode=None, timesteps=None, raise_not_converged=True, **kwargs)
             for all time steps.
             Default: True.
 
+        troubleshooting_mode : str or None
+            Two optional troubleshooting methods in case of nonconvergence of nonlinear
+            power flow (cf. [1])
+
+            * None (default)
+                Power flow analysis is conducted using nonlinear power flow method.
+            * 'lpf'
+                Non-linear power flow initial guess is seeded with the voltage angles
+                from the linear power flow.
+            * 'iteration'
+                Power flow analysis is conducted by reducing all power values of
+                generators and loads to a fraction, e.g. 10%, solving the load flow and
+                using it as a seed for the power at 20%, iteratively up to 100%.
+
+        range_start : float, optional
+            Specifies the minimum fraction that power values are set to when using
+            troubleshooting_mode 'iteration'. Must be between 0 and 1.
+            Default: 0.1.
+
+        range_num : int, optional
+            Specifies the number of fraction samples to generate when using
+            troubleshooting_mode 'iteration'. Must be non-negative.
+            Default: 10.
+
         Other Parameters
         -----------------
         kwargs : dict
@@ -731,7 +766,37 @@ def analyze(self, mode=None, timesteps=None, raise_not_converged=True, **kwargs)
         :pandas:`pandas.DatetimeIndex<DatetimeIndex>`
             Returns the time steps for which power flow analysis did not converge.
 
+        References
+        --------
+        [1] https://pypsa.readthedocs.io/en/latest/troubleshooting.html
+
         """
+
+        def _check_convergence():
+            # get converged and not converged time steps
+            timesteps_converged = pf_results["converged"][
+                pf_results["converged"]["0"]
+            ].index
+            timesteps_not_converged = pf_results["converged"][
+                ~pf_results["converged"]["0"]
+            ].index
+
+            if raise_not_converged and len(timesteps_not_converged) > 0:
+                raise ValueError(
+                    "Power flow analysis did not converge for the "
+                    "following {} time steps: {}.".format(
+                        len(timesteps_not_converged), timesteps_not_converged
+                    )
+                )
+            elif len(timesteps_not_converged) > 0:
+                logger.warning(
+                    "Power flow analysis did not converge for the "
+                    "following {} time steps: {}.".format(
+                        len(timesteps_not_converged), timesteps_not_converged
+                    )
+                )
+            return timesteps_converged, timesteps_not_converged
+
         if timesteps is None:
             timesteps = self.timeseries.timeindex
         # check if timesteps is array-like, otherwise convert to list
@@ -740,31 +805,39 @@ def analyze(self, mode=None, timesteps=None, raise_not_converged=True, **kwargs)
 
         pypsa_network = self.to_pypsa(mode=mode, timesteps=timesteps, **kwargs)
 
-        # run power flow analysis
-        pf_results = pypsa_network.pf(timesteps, use_seed=kwargs.get("use_seed", False))
-
-        # get converged and not converged time steps
-        timesteps_converged = pf_results["converged"][
-            pf_results["converged"]["0"]
-        ].index
-        timesteps_not_converged = pf_results["converged"][
-            ~pf_results["converged"]["0"]
-        ].index
-
-        if raise_not_converged and len(timesteps_not_converged) > 0:
-            raise ValueError(
-                "Power flow analysis did not converge for the "
-                "following {} time steps: {}.".format(
-                    len(timesteps_not_converged), timesteps_not_converged
-                )
-            )
-        elif len(timesteps_not_converged) > 0:
-            logger.warning(
-                "Power flow analysis did not converge for the "
-                "following {} time steps: {}.".format(
-                    len(timesteps_not_converged), timesteps_not_converged
+        if troubleshooting_mode == "lpf":
+            # run linear power flow analysis
+            pypsa_network.lpf()
+            # run power flow analysis
+            pf_results = pypsa_network.pf(timesteps, use_seed=True)
+            # get converged and not converged time steps
+            timesteps_converged, timesteps_not_converged = _check_convergence()
+        elif troubleshooting_mode == "iteration":
+            pypsa_network_copy = pypsa_network.copy()
+            for fraction in np.linspace(range_start, 1, range_num):
+                # Reduce power values of generators, loads and storages to fraction of
+                # original value
+                for obj1, obj2 in [
+                    (pypsa_network.generators_t, pypsa_network_copy.generators_t),
+                    (pypsa_network.loads_t, pypsa_network_copy.loads_t),
+                    (pypsa_network.storage_units_t, pypsa_network_copy.storage_units_t),
+                ]:
+                    for attr in ["p_set", "q_set"]:
+                        setattr(obj1, attr, getattr(obj2, attr) * fraction)
+                # run power flow analysis
+                pf_results = pypsa_network.pf(timesteps, use_seed=True)
+                logging.warning(
+                    "Current fraction in iterative process: {}.".format(fraction)
                 )
+                # get converged and not converged time steps
+                timesteps_converged, timesteps_not_converged = _check_convergence()
+        else:
+            # run power flow analysis
+            pf_results = pypsa_network.pf(
+                timesteps, use_seed=kwargs.get("use_seed", False)
             )
+            # get converged and not converged time steps
+            timesteps_converged, timesteps_not_converged = _check_convergence()
 
         # handle converged time steps
         pypsa_io.process_pfa_results(self, pypsa_network, timesteps_converged)

edisgo/network/results.py

@@ -807,6 +807,41 @@ def reduce_memory(self, attr_to_reduce=None, to_type="float32"):
         for attr in attr_to_reduce:
             setattr(self, attr, getattr(self, attr).astype(to_type))
 
+    def equality_check(self, results_obj):
+        """
+        Checks the equality of two results objects.
+
+        Parameters
+        ----------
+        results_obj : :class:~.network.results.Results
+            Contains the results of analyze function with default settings.
+
+        Returns
+        -------
+        bool
+            True if equality check is successful, False otherwise.
+
+        """
+
+        attr_to_check = [
+            "pfa_p",
+            "pfa_q",
+            "pfa_v_ang_seed",
+            "pfa_v_mag_pu_seed",
+            "v_res",
+            "i_res",
+        ]
+        try:
+            for attr in attr_to_check:
+                pd.testing.assert_frame_equal(
+                    getattr(self, attr),
+                    getattr(results_obj, attr),
+                    check_freq=False,
+                )
+            return True
+        except AssertionError:
+            return False
+
     def to_csv(
         self, directory, parameters=None, reduce_memory=False, save_seed=False, **kwargs
     ):

examples/edisgo_simple_example.ipynb

@@ -657,6 +657,7 @@
        "    <tr style=\"text-align: right;\">\n",
        "      <th></th>\n",
        "      <th>bus</th>\n",
+       "      <th>control</th>\n",
        "      <th>p_nom</th>\n",
        "      <th>type</th>\n",
        "      <th>control</th>\n",
@@ -2068,6 +2069,34 @@
       "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.025066 seconds\n",
       "/home/birgit/virtualenvs/edisgo_emob/lib/python3.8/site-packages/pandas/core/series.py:1056: SettingWithCopyWarning:\n",
       "\n",
+      "See the documentation here:\n",
+      "https://pandas.pydata.org/pandas-docs/stable/indexing.html#deprecate-loc-reindex-listlike\n",
+      "  return self._getitem_tuple(key)\n",
+      "INFO:pypsa.pf:Performing non-linear load-flow on AC sub-network SubNetwork 0 for snapshots DatetimeIndex(['1970-01-01 00:00:00', '1970-01-01 01:00:00'], dtype='datetime64[ns]', freq='H')\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.126566 seconds\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.111184 seconds\n",
+      "INFO:pypsa.pf:Performing non-linear load-flow on AC sub-network SubNetwork 0 for snapshots DatetimeIndex(['1970-01-01 00:00:00', '1970-01-01 01:00:00'], dtype='datetime64[ns]', freq='H')\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.088311 seconds\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.112642 seconds\n",
+      "INFO:edisgo:==> Load issues were solved in 1 iteration step(s).\n",
+      "INFO:edisgo:==> Voltage issues in MV topology were solved in 0 iteration step(s).\n",
+      "INFO:pypsa.pf:Performing non-linear load-flow on AC sub-network SubNetwork 0 for snapshots DatetimeIndex(['1970-01-01 00:00:00', '1970-01-01 01:00:00'], dtype='datetime64[ns]', freq='H')\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.087276 seconds\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.105445 seconds\n",
+      "INFO:edisgo:==> Voltage issues at busbars in LV grids were solved in 1 iteration step(s).\n",
+      "INFO:pypsa.pf:Performing non-linear load-flow on AC sub-network SubNetwork 0 for snapshots DatetimeIndex(['1970-01-01 00:00:00', '1970-01-01 01:00:00'], dtype='datetime64[ns]', freq='H')\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.086914 seconds\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.139845 seconds\n",
+      "INFO:pypsa.pf:Performing non-linear load-flow on AC sub-network SubNetwork 0 for snapshots DatetimeIndex(['1970-01-01 00:00:00', '1970-01-01 01:00:00'], dtype='datetime64[ns]', freq='H')\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.100760 seconds\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.101456 seconds\n",
+      "INFO:pypsa.pf:Performing non-linear load-flow on AC sub-network SubNetwork 0 for snapshots DatetimeIndex(['1970-01-01 00:00:00', '1970-01-01 01:00:00'], dtype='datetime64[ns]', freq='H')\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.096341 seconds\n",
+      "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.106207 seconds\n",
+      "INFO:edisgo:==> Voltage issues in LV grids were solved in 3 iteration step(s).\n",
+      "INFO:edisgo:==> Load issues were rechecked and solved in 0 iteration step(s).\n",
+      "/home/birgit/virtualenvs/edisgo/lib/python3.6/site-packages/pandas/core/frame.py:6701: FutureWarning: Sorting because non-concatenation axis is not aligned. A future version\n",
+      "of pandas will change to not sort by default.\n",
       "\n",
       "A value is trying to be set on a copy of a slice from a DataFrame\n",
       "\n",
@@ -2890,6 +2919,16 @@
       "INFO:pypsa.pf:Newton-Raphson solved in 3 iterations with error of 0.000000 in 0.298490 seconds\n",
       "INFO:edisgo.flex_opt.reinforce_grid:==> Voltage issues at busbars in LV grids were solved in 7 iteration step(s).\n"
      ]
+    },
+    {
+     "data": {
+      "text/plain": [
+       "<edisgo.network.results.Results at 0x7fef449ba220>"
+      ]
+     },
+     "execution_count": 41,
+     "metadata": {},
+     "output_type": "execute_result"
     }
    ],
    "source": [

tests/network/test_electromobility.py

@@ -13,15 +13,17 @@
 
 class TestElectromobility:
     @classmethod
-    def setup_class(self):
-        self.ding0_path = pytest.ding0_test_network_3_path
-        self.simbev_path = pytest.simbev_example_scenario_path
-        self.tracbev_path = pytest.tracbev_example_scenario_path
-        self.standing_times_path = os.path.join(self.simbev_path, "simbev_run")
-        self.charging_strategies = ["dumb", "reduced", "residual"]
-
-        self.edisgo_obj = import_edisgo_from_files(
-            self.ding0_path, import_topology=True, import_timeseries=True
+    def setup_class(cls):
+        cls.ding0_path = pytest.ding0_test_network_3_path
+        cls.simbev_path = pytest.simbev_example_scenario_path
+        cls.tracbev_path = pytest.tracbev_example_scenario_path
+        cls.standing_times_path = os.path.join(cls.simbev_path, "simbev_run")
+        cls.charging_strategies = ["dumb", "reduced", "residual"]
+
+        cls.edisgo_obj = import_edisgo_from_files(
+            cls.ding0_path,
+            import_topology=True,
+            import_timeseries=True,
         )
 
     def test_import_simbev_electromobility(self):
@@ -130,7 +132,7 @@ def test_integrate_charging_parks(self):
     def test_charging_strategy(self):
         charging_demand_lst = []
 
-        for count, strategy in enumerate(self.charging_strategies):
+        for strategy in self.charging_strategies:
             charging_strategy(self.edisgo_obj, strategy=strategy)
 
             ts = self.edisgo_obj.timeseries