add print_pf_summary; add first draft for diagnostic function

pandapipes/diagnostic.py

@@ -0,0 +1,145 @@
+# Copyright (c) 2020-2023 by Fraunhofer Institute for Energy Economics
+# and Energy System Technology (IEE), Kassel, and University of Kassel. All rights reserved.
+# Use of this source code is governed by a BSD-style license that can be found in the LICENSE file.
+
+import pandapipes as pp
+import numpy as np
+
+from pandapipes import PipeflowNotConverged
+
+try:
+    import pandaplan.core.pplog as logging
+except ImportError:
+    import logging
+
+logger = logging.getLogger(__name__)
+
+
+def check_net(net, low_length_limit_km=0.01, check_scaling_factor=1e-5):
+    """
+    Run some diagnostic checks on the net to identify potential flaws.
+    """
+    net = net.deepcopy()  # do not modify the direct input
+    try:
+        pp.pipeflow(net)
+        if net.converged:
+            logger.info("The initial, unmodified pipeflow converges.")
+        else:
+            logger.warning("The initial, unmodified pipeflow does NOT converge.")
+    except Exception as e:
+        logger.info(f"The initial, unmodified pipeflow does NOT converge.\n"
+                    f"\t\tThis exception is raised:\n\t\t{e}")
+
+    # check ext_grid
+    if net.fluid.is_gas & (not hasattr(net, "ext_grid") | net.ext_grid.empty):
+        logger.warning("The net does not have an external grid! "
+                       "An external grid is required for gas networks.")
+
+    # check zero / low length
+    zl = net.pipe.loc[net.pipe.length_km == 0]
+    ll = net.pipe.loc[net.pipe.length_km <= low_length_limit_km]
+    if not zl.empty:
+        logger.warning(f"{len(zl.index)} pipes have a length of 0.0 km. (IDs: {zl.index})")
+    if not ll.empty:
+        logger.warning(f"{len(ll.index)} pipes have a length below"
+                       f" {low_length_limit_km} km. "
+                       f"This could lead to convergence issues. The lowest length in the net is "
+                       f"{ll.length_km.min()} km.\n"
+                       f"(IDs of pipelines with low length: {ll.index})")
+
+    net2 = net.deepcopy()
+    net2.pipe.loc[net2.pipe.length_km < low_length_limit_km].length_km = low_length_limit_km
+    try:
+        pp.pipeflow(net2)
+        if net2.converged:
+            logger.info(f"If all short pipelines (< {low_length_limit_km} km) were set to "
+                        f"{low_length_limit_km} km, the pipeflow would converge.")
+        else:
+            logger.warning(f"If all short pipelines (< {low_length_limit_km} km) were set to "
+                        f"{low_length_limit_km} km, the pipeflow would still NOT converge.")
+    except Exception as e:
+        logger.info(f"Pipeflow does not converge, even if all short pipelines (< 10 m) were set "
+                    f"to {low_length_limit_km} km. \n"
+                    f"\t\tThe error message is: {e}")
+
+    # check iterations
+    iter = 200
+    try:
+        pp.pipeflow(net, iter=iter)
+        logger.info(f"The pipeflow converges after {net._internal_results['iterations']:d} "
+                    f"iterations.")
+    except PipeflowNotConverged:
+        logger.info(f"After {iter:d} iterations the pipeflow did NOT converge.")
+
+    # check with little sink and source scaling
+    logger.info("Testing with scaled-down sinks and sources.")
+    net3 = net.deepcopy()
+    if hasattr(net, "sink"):
+        net3.sink.scaling *= check_scaling_factor
+    if hasattr(net, "source"):
+        net3.source.scaling *= check_scaling_factor
+    try:
+        pp.pipeflow(net3)
+        if net3.converged:
+            logger.info(f"If sinks and sources were scaled with a factor of to "
+                        f"{check_scaling_factor}, the pipeflow would converge.")
+        else:
+            logger.warning(f"If sinks and sources were scaled with a factor of to "
+                           f"{check_scaling_factor}, the pipeflow would still NOT converge.")
+    except Exception as e:
+        logger.info(f"Pipeflow does not converge with sinks/sources scaled by"
+                    f" {check_scaling_factor}.\n"
+                    f"\t\tThe error message is: {e}")
+
+    # check k
+    if any(net.pipe.k_mm > 0.5):
+        logger.warning(f"Some pipes have a friction factor k_mm > 0.5 (extremely rough). The "
+                       f"highest value in the net is {net.pipe.k_mm.max()}. Up to "
+                       f"0.2 mm is a common value for old steel pipes."
+                       f"\nRough pipes: {net.pipe.loc[net.pipe.k_mm > 0.5]}.")
+    net4 = net.deepcopy()
+    net4.pipe.k_mm = 1e-5
+    try:
+        pp.pipeflow(net4)
+        if net4.converged:
+            logger.info(f"If the friction factor would be reduced to 1e-5 for all pipes, "
+                        f"the pipeflow would converge.")
+        else:
+            logger.warning(f"If the friction factor would be reduced to 1e-5 for all pipes, "
+                           f"the pipeflow would still NOT converge.")
+    except Exception as e:
+        logger.info(f"Pipeflow does not converge with k_mm = 1-e5 for all pipes.\n"
+                    f"\t\tThe error message is: {e}")
+
+    # check sink and source junctions:
+    node_component = ["sink", "source", "ext_grid"]
+    for nc in node_component:
+        if hasattr(net, nc):
+            missing = np.setdiff1d(net[nc].junction, net.junction.index)
+            if len(missing):
+                logger.warning(f"Some {nc}s are connected to non-existing junctions!"
+                               f"\n{nc}s:{net[nc].loc[net[nc].junction.isin(missing)]}"
+                               f"\nmissing junctions:{missing}")
+
+    # check from and to junctions
+    branch_component = ["pipe", "valve", "compressor", "pump", "heat_exchanger", "circulation_pump"]
+    for bc in branch_component:
+        if hasattr(net, bc):
+            missing_f = np.setdiff1d(net[bc].from_junction, net.junction.index)
+            missing_t = np.setdiff1d(net[bc].to_junction, net.junction.index)
+            if len(missing_t) | len(missing_t):
+                logger.warning(f"Some {bc}s are connected to non-existing junctions!")
+                logger.warning(f"missing 'from' junctions:{missing_f}")
+                logger.warning(f"missing 'to' junctions:{missing_t}")
+
+
+if __name__ == '__main__':
+    import pandapipes.networks
+    net = pandapipes.networks.schutterwald()
+    net.ext_grid.p_bar = 0.08
+    net.sink.loc[1505, "mdot_kg_per_s"] = 1000
+    try:
+        pandapipes.pipeflow(net)
+    except Exception as e:
+        print(f"pipeflow raised: \n {e}")
+    check_net(net)

pandapipes/toolbox.py

@@ -646,3 +646,40 @@ def get_internal_tables_pandas(net):
     branch_table["TABLE_IDX"].replace(branch_table_lookup["n2t"], inplace=True)
 
     return node_table, branch_table
+
+
+def print_pf_summary(net):
+    """Print some basic results of the pipeflow.
+
+    Min./max. pressure, junctions with NaN results, max. v, sum of sinks / sources.
+    """
+    if not net.converged:
+        try:
+            pandapipes.pipeflow(net)
+        except Exception as e:
+            return logger.Error(f"Could not print pipeflow summary because the pipeflow "
+                                f"calculation was not successful (Exception: {e})")
+        if not net.converged:
+            return logger.Error(f"Could not print pipeflow summary because the pipeflow "
+                                f"calculation did not converge.")
+    else:
+        if any(net.res_junction.loc[net.junction.in_service].p_bar.isna()):
+            print(f"For {sum(net.res_junction.loc[net.junction.in_service].p_bar.isna())} "
+                  f"junctions, no pressure could be calculated (NaN)!")
+        print(f"The minimum pressure is {net.res_junction.p_bar.min():.2f} bar.")
+        print(f"The maximum pressure is {net.res_junction.p_bar.max():.2f} bar.")
+        print(f"The highest velocity is {net.res_pipe.v_mean_m_per_s.abs().max():.2f} m/s.")
+        if hasattr(net, "source") & (~net.source.empty):
+            total_source_mdot = net.res_source.mdot_kg_per_s.sum()
+            print(f"The total gas infeed from sources is {total_source_mdot:.2f} kg/s "
+                  f"(i.e. {total_source_mdot * float(net.fluid.get_property('hhv'))*3.6:.2f} "
+                  f"MW_th (HHV).)")
+        else:
+            print("There are no sources connected to the net.")
+        if hasattr(net, "sink") & (~net.sink.empty):
+            total_sink_mdot = net.res_sink.mdot_kg_per_s.sum()
+            print(f"The total gas demand from sinks is {total_sink_mdot:.2f} kg/s "
+                  f"(i.e. {total_sink_mdot * float(net.fluid.get_property('hhv'))*3.6:.2f} MW_th "
+                  f"(HHV).)")
+        else:
+            print("There are no sinks connected to the net.")