diff --git a/CHANGELOG.rst b/CHANGELOG.rst
index fb2e9937..20a3f60f 100644
--- a/CHANGELOG.rst
+++ b/CHANGELOG.rst
@@ -6,6 +6,7 @@ Change Log
 - [ADDED] multiple creation of heat exchanger
 - [ADDED] support Python 3.11 (now included in test pipeline)
 - [CHANGED] dropped support for Python 3.7 (no longer included in test pipeline)
+- [CHANGED] connectivity check now separated by hydraulics and heat_transfer calculation, so that also results can differ in some rows (NaN or not)
 - [CHANGED] dynamic creation of lookups for getting pit as pandas tables
 - [FIXED] in STANET converter: bug fix for heat exchanger creation and external temperatures of pipes added
 - [REMOVED] broken travis badge removed from readme
diff --git a/pandapipes/component_models/abstract_models/base_component.py b/pandapipes/component_models/abstract_models/base_component.py
index cdfb0f09..b5ce3a60 100644
--- a/pandapipes/component_models/abstract_models/base_component.py
+++ b/pandapipes/component_models/abstract_models/base_component.py
@@ -33,7 +33,7 @@ def init_results(cls, net):
         return res_table
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         """
         Function that extracts certain results.
 
@@ -43,10 +43,8 @@ def extract_results(cls, net, options, branch_results, nodes_connected, branches
         :type options:
         :param branch_results:
         :type branch_results:
-        :param nodes_connected:
-        :type nodes_connected:
-        :param branches_connected:
-        :type branches_connected:
+        :param mode:
+        :type mode:
         :return: No Output.
         """
         raise NotImplementedError
diff --git a/pandapipes/component_models/abstract_models/branch_models.py b/pandapipes/component_models/abstract_models/branch_models.py
index 3ee41146..4eeac9cc 100644
--- a/pandapipes/component_models/abstract_models/branch_models.py
+++ b/pandapipes/component_models/abstract_models/branch_models.py
@@ -159,5 +159,5 @@ def calculate_temperature_lift(cls, net, branch_component_pit, node_pit):
         raise NotImplementedError
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         raise NotImplementedError
diff --git a/pandapipes/component_models/abstract_models/branch_w_internals_models.py b/pandapipes/component_models/abstract_models/branch_w_internals_models.py
index 6861223e..1b88f0ea 100644
--- a/pandapipes/component_models/abstract_models/branch_w_internals_models.py
+++ b/pandapipes/component_models/abstract_models/branch_w_internals_models.py
@@ -217,7 +217,7 @@ def get_internal_pipe_number(cls, net):
         return np.array(net[cls.table_name()].sections.values)
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         raise NotImplementedError
 
     @classmethod
diff --git a/pandapipes/component_models/abstract_models/branch_wo_internals_models.py b/pandapipes/component_models/abstract_models/branch_wo_internals_models.py
index 0efa325a..106895d5 100644
--- a/pandapipes/component_models/abstract_models/branch_wo_internals_models.py
+++ b/pandapipes/component_models/abstract_models/branch_wo_internals_models.py
@@ -100,5 +100,5 @@ def get_connected_node_type(cls):
         raise NotImplementedError
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         raise NotImplementedError
diff --git a/pandapipes/component_models/abstract_models/branch_wzerolength_models.py b/pandapipes/component_models/abstract_models/branch_wzerolength_models.py
index dc224f8a..b9da3028 100644
--- a/pandapipes/component_models/abstract_models/branch_wzerolength_models.py
+++ b/pandapipes/component_models/abstract_models/branch_wzerolength_models.py
@@ -61,7 +61,7 @@ def create_pit_branch_entries(cls, net, branch_pit):
         return branch_wzerolength_pit
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         raise NotImplementedError
 
     @classmethod
diff --git a/pandapipes/component_models/abstract_models/circulation_pump.py b/pandapipes/component_models/abstract_models/circulation_pump.py
index cbcbfb60..767dc336 100644
--- a/pandapipes/component_models/abstract_models/circulation_pump.py
+++ b/pandapipes/component_models/abstract_models/circulation_pump.py
@@ -97,14 +97,12 @@ def calculate_temperature_lift(cls, net, pipe_pit, node_pit):
         raise NotImplementedError
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         """
         Function that extracts certain results.
 
-        :param nodes_connected:
-        :type nodes_connected:
-        :param branches_connected:
-        :type branches_connected:
+        :param mode:
+        :type mode:
         :param branch_results:
         :type branch_results:
         :param net: The pandapipes network
diff --git a/pandapipes/component_models/abstract_models/const_flow_models.py b/pandapipes/component_models/abstract_models/const_flow_models.py
index c49a58ec..a10a4d0e 100644
--- a/pandapipes/component_models/abstract_models/const_flow_models.py
+++ b/pandapipes/component_models/abstract_models/const_flow_models.py
@@ -51,14 +51,12 @@ def create_pit_node_entries(cls, net, node_pit):
         node_pit[index, LOAD] += loads_sum
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         """
         Function that extracts certain results.
 
-        :param nodes_connected:
-        :type nodes_connected:
-        :param branches_connected:
-        :type branches_connected:
+        :param mode:
+        :type mode:
         :param branch_results:
         :type branch_results:
         :param net: The pandapipes network
@@ -74,8 +72,8 @@ def extract_results(cls, net, options, branch_results, nodes_connected, branches
         is_loads = loads.in_service.values
         fj, tj = get_lookup(net, "node", "from_to")[cls.get_connected_node_type().table_name()]
         junct_pit = net["_pit"]["node"][fj:tj, :]
-        nodes_connected = get_lookup(net, "node", "active")[fj:tj]
-        is_juncts = np.isin(loads.junction.values, junct_pit[nodes_connected, ELEMENT_IDX])
+        nodes_connected_hyd = get_lookup(net, "node", "active_hydraulics")[fj:tj]
+        is_juncts = np.isin(loads.junction.values, junct_pit[nodes_connected_hyd, ELEMENT_IDX])
 
         is_calc = is_loads & is_juncts
         res_table["mdot_kg_per_s"].values[is_calc] = loads.mdot_kg_per_s.values[is_calc] \
diff --git a/pandapipes/component_models/abstract_models/node_element_models.py b/pandapipes/component_models/abstract_models/node_element_models.py
index 04e594e2..0792fa53 100644
--- a/pandapipes/component_models/abstract_models/node_element_models.py
+++ b/pandapipes/component_models/abstract_models/node_element_models.py
@@ -51,5 +51,5 @@ def create_pit_node_entries(cls, net, node_pit):
         raise NotImplementedError
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         raise NotImplementedError
diff --git a/pandapipes/component_models/abstract_models/node_models.py b/pandapipes/component_models/abstract_models/node_models.py
index 5a917bef..0ce18e4a 100644
--- a/pandapipes/component_models/abstract_models/node_models.py
+++ b/pandapipes/component_models/abstract_models/node_models.py
@@ -66,5 +66,5 @@ def get_result_table(cls, net):
         raise NotImplementedError
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         raise NotImplementedError
diff --git a/pandapipes/component_models/component_toolbox.py b/pandapipes/component_models/component_toolbox.py
index 04f357b3..2db03d12 100644
--- a/pandapipes/component_models/component_toolbox.py
+++ b/pandapipes/component_models/component_toolbox.py
@@ -5,6 +5,7 @@
 import numpy as np
 import pandas as pd
 
+from pandapipes import get_fluid
 from pandapipes.constants import NORMAL_PRESSURE, TEMP_GRADIENT_KPM, AVG_TEMPERATURE_K, \
     HEIGHT_EXPONENT
 from pandapipes.idx_branch import LOAD_VEC_NODES, FROM_NODE, TO_NODE
@@ -177,3 +178,23 @@ def get_mass_flow_at_nodes(net, node_pit, branch_pit, eg_nodes, comp):
         raise UserWarning("In component %s: Something went wrong with the mass flow balance. "
                           "Please report this error at github." % comp.__name__)
     return sum_mass_flows, inverse_nodes, counts
+
+
+def standard_branch_wo_internals_result_lookup(net):
+    required_results_hyd = [
+        ("p_from_bar", "p_from"), ("p_to_bar", "p_to"), ("mdot_to_kg_per_s", "mf_to"),
+        ("mdot_from_kg_per_s", "mf_from"), ("vdot_norm_m3_per_s", "vf"), ("lambda", "lambda"),
+        ("reynolds", "reynolds")
+    ]
+    required_results_ht = [("t_from_k", "temp_from"), ("t_to_k", "temp_to")]
+
+    if get_fluid(net).is_gas:
+        required_results_hyd.extend([
+            ("v_from_m_per_s", "v_gas_from"), ("v_to_m_per_s", "v_gas_to"),
+            ("v_mean_m_per_s", "v_gas_mean"), ("normfactor_from", "normfactor_from"),
+            ("normfactor_to", "normfactor_to")
+        ])
+    else:
+        required_results_hyd.extend([("v_mean_m_per_s", "v_mps")])
+
+    return required_results_hyd, required_results_ht
\ No newline at end of file
diff --git a/pandapipes/component_models/ext_grid_component.py b/pandapipes/component_models/ext_grid_component.py
index 9570f51e..d7802832 100644
--- a/pandapipes/component_models/ext_grid_component.py
+++ b/pandapipes/component_models/ext_grid_component.py
@@ -62,20 +62,18 @@ def create_pit_node_entries(cls, net, node_pit):
         return ext_grids, press
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         """
         Function that extracts certain results.
 
-        :param nodes_connected:
-        :type nodes_connected:
-        :param branches_connected:
-        :type branches_connected:
         :param branch_results:
         :type branch_results:
         :param net: The pandapipes network
         :type net: pandapipesNet
         :param options:
         :type options:
+        :param mode:
+        :type mode:
         :return: No Output.
         """
         ext_grids = net[cls.table_name()]
diff --git a/pandapipes/component_models/flow_control_component.py b/pandapipes/component_models/flow_control_component.py
index d42f5815..739e5ada 100644
--- a/pandapipes/component_models/flow_control_component.py
+++ b/pandapipes/component_models/flow_control_component.py
@@ -5,8 +5,9 @@
 import numpy as np
 from numpy import dtype
 
-from pandapipes.component_models.junction_component import Junction
 from pandapipes.component_models.abstract_models import BranchWZeroLengthComponent, get_fluid
+from pandapipes.component_models.component_toolbox import standard_branch_wo_internals_result_lookup
+from pandapipes.component_models.junction_component import Junction
 from pandapipes.idx_branch import D, AREA, TL, JAC_DERIV_DP, JAC_DERIV_DP1, JAC_DERIV_DV, VINIT, \
     RHO, LOAD_VEC_BRANCHES, ELEMENT_IDX
 from pandapipes.pf.result_extraction import extract_branch_results_without_internals
@@ -83,24 +84,11 @@ def calculate_temperature_lift(cls, net, branch_component_pit, node_pit):
         branch_component_pit[:, TL] = 0
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
-        required_results = [
-            ("p_from_bar", "p_from"), ("p_to_bar", "p_to"), ("t_from_k", "temp_from"),
-            ("t_to_k", "temp_to"), ("mdot_to_kg_per_s", "mf_to"), ("mdot_from_kg_per_s", "mf_from"),
-            ("vdot_norm_m3_per_s", "vf"), ("lambda", "lambda"), ("reynolds", "reynolds")
-        ]
-
-        if get_fluid(net).is_gas:
-            required_results.extend([
-                ("v_from_m_per_s", "v_gas_from"), ("v_to_m_per_s", "v_gas_to"),
-                ("v_mean_m_per_s", "v_gas_mean"), ("normfactor_from", "normfactor_from"),
-                ("normfactor_to", "normfactor_to")
-            ])
-        else:
-            required_results.extend([("v_mean_m_per_s", "v_mps")])
+    def extract_results(cls, net, options, branch_results, mode):
+        required_results_hyd, required_results_ht = standard_branch_wo_internals_result_lookup(net)
 
-        extract_branch_results_without_internals(net, branch_results, required_results,
-                                                 cls.table_name(), branches_connected)
+        extract_branch_results_without_internals(net, branch_results, required_results_hyd,
+                                                 required_results_ht, cls.table_name(), mode)
 
     @classmethod
     def get_component_input(cls):
diff --git a/pandapipes/component_models/heat_exchanger_component.py b/pandapipes/component_models/heat_exchanger_component.py
index abbb2259..6304bb4d 100644
--- a/pandapipes/component_models/heat_exchanger_component.py
+++ b/pandapipes/component_models/heat_exchanger_component.py
@@ -5,6 +5,7 @@
 import numpy as np
 from numpy import dtype
 
+from pandapipes.component_models import standard_branch_wo_internals_result_lookup
 from pandapipes.component_models.abstract_models.branch_wzerolength_models import \
     BranchWZeroLengthComponent
 from pandapipes.component_models.junction_component import Junction
@@ -60,24 +61,25 @@ def create_pit_branch_entries(cls, net, branch_pit):
         heat_exchanger_pit[:, T_OUT] = 307
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
-        required_results = [
-            ("p_from_bar", "p_from"), ("p_to_bar", "p_to"), ("t_from_k", "temp_from"),
-            ("t_to_k", "temp_to"), ("mdot_to_kg_per_s", "mf_to"), ("mdot_from_kg_per_s", "mf_from"),
-            ("vdot_norm_m3_per_s", "vf"), ("lambda", "lambda"), ("reynolds", "reynolds")
-        ]
+    def extract_results(cls, net, options, branch_results, mode):
+        """
+        Class method to extract pipeflow results from the internal structure into the results table.
 
-        if get_fluid(net).is_gas:
-            required_results.extend([
-                ("v_from_m_per_s", "v_gas_from"), ("v_to_m_per_s", "v_gas_to"),
-                ("v_mean_m_per_s", "v_gas_mean"), ("normfactor_from", "normfactor_from"),
-                ("normfactor_to", "normfactor_to")
-            ])
-        else:
-            required_results.extend([("v_mean_m_per_s", "v_mps")])
+        :param net: The pandapipes network
+        :type net: pandapipesNet
+        :param options: pipeflow options
+        :type options: dict
+        :param branch_results: important branch results
+        :type branch_results: dict
+        :param mode: simulation mode
+        :type mode: str
+        :return: No Output.
+        :rtype: None
+        """
+        required_results_hyd, required_results_ht = standard_branch_wo_internals_result_lookup(net)
 
-        extract_branch_results_without_internals(net, branch_results, required_results,
-                                                 cls.table_name(), branches_connected)
+        extract_branch_results_without_internals(net, branch_results, required_results_hyd,
+                                                 required_results_ht, cls.table_name(), mode)
 
     @classmethod
     def calculate_temperature_lift(cls, net, branch_component_pit, node_pit):
diff --git a/pandapipes/component_models/junction_component.py b/pandapipes/component_models/junction_component.py
index 3426fae0..fdf2a22a 100644
--- a/pandapipes/component_models/junction_component.py
+++ b/pandapipes/component_models/junction_component.py
@@ -88,37 +88,44 @@ def create_pit_node_entries(cls, net, node_pit):
         junction_pit[:, ACTIVE_ND] = junctions.in_service.values
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         """
         Function that extracts certain results.
 
-        :param nodes_connected:
-        :type nodes_connected:
-        :param branches_connected:
-        :type branches_connected:
-        :param branch_results:
-        :type branch_results:
+        :param mode:
+        :type mode:
         :param net: The pandapipes network
         :type net: pandapipesNet
         :param options:
         :type options:
+        :param branch_results:
+        :type branch_results:
+        :param mode:
+        :type mode:
         :return: No Output.
         """
         res_table = net["res_" + cls.table_name()]
 
         f, t = get_lookup(net, "node", "from_to")[cls.table_name()]
-        fa, ta = get_lookup(net, "node", "from_to_active")[cls.table_name()]
-
-        junction_pit = net["_active_pit"]["node"][fa:ta, :]
-        junctions_active = get_lookup(net, "node", "active")[f:t]
-
-        if np.any(junction_pit[:, PINIT] < 0):
-            warn(UserWarning('Pipeflow converged, however, the results are physically incorrect '
-                             'as pressure is negative at nodes %s'
-                             % junction_pit[junction_pit[:, PINIT] < 0, ELEMENT_IDX]))
-
-        res_table["p_bar"].values[junctions_active] = junction_pit[:, PINIT]
-        res_table["t_k"].values[junctions_active] = junction_pit[:, TINIT]
+        junction_pit = net["_pit"]["node"][f:t, :]
+
+        if mode in ["hydraulics", "all"]:
+            junctions_connected_hydraulic = get_lookup(net, "node", "active_hydraulics")[f:t]
+
+            if np.any(junction_pit[junctions_connected_hydraulic, PINIT] < 0):
+                warn(UserWarning('Pipeflow converged, however, the results are physically incorrect '
+                                 'as pressure is negative at nodes %s'
+                                 % junction_pit[junction_pit[:, PINIT] < 0, ELEMENT_IDX]))
+
+        #     res_table["p_bar"].values[junctions_connected_hydraulic] = junction_pit[:, PINIT]
+        #     if mode == "hydraulics":
+        #         res_table["t_k"].values[junctions_connected_hydraulic] = junction_pit[:, TINIT]
+        #
+        # if mode in ["heat", "all"]:
+        #     junctions_connected_ht = get_lookup(net, "node", "active_heat_transfer")[f:t]
+        #     res_table["t_k"].values[junctions_connected_ht] = junction_pit[:, TINIT]
+        res_table["p_bar"].values[:] = junction_pit[:, PINIT]
+        res_table["t_k"].values[:] = junction_pit[:, TINIT]
 
     @classmethod
     def get_component_input(cls):
diff --git a/pandapipes/component_models/pipe_component.py b/pandapipes/component_models/pipe_component.py
index c81a41eb..af399f15 100644
--- a/pandapipes/component_models/pipe_component.py
+++ b/pandapipes/component_models/pipe_component.py
@@ -172,28 +172,35 @@ def calculate_temperature_lift(cls, net, branch_component_pit, node_pit):
         branch_component_pit[:, TL] = 0
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
-        res_nodes_from = [("p_from_bar", "p_from"), ("t_from_k", "temp_from"),
-                          ("mdot_from_kg_per_s", "mf_from")]
-        res_nodes_to = [("p_to_bar", "p_to"), ("t_to_k", "temp_to"), ("mdot_to_kg_per_s", "mf_to")]
-        res_mean = [("vdot_norm_m3_per_s", "vf"), ("lambda", "lambda"), ("reynolds", "reynolds")]
+    def extract_results(cls, net, options, branch_results, mode):
+        res_nodes_from_hyd = [("p_from_bar", "p_from"), ("mdot_from_kg_per_s", "mf_from")]
+        res_nodes_from_ht = [("t_from_k", "temp_from")]
+        res_nodes_to_hyd = [("p_to_bar", "p_to"), ("mdot_to_kg_per_s", "mf_to")]
+        res_nodes_to_ht = [("t_to_k", "temp_to")]
+        res_mean_hyd = [("vdot_norm_m3_per_s", "vf"), ("lambda", "lambda"),
+                        ("reynolds", "reynolds")]
 
         if get_fluid(net).is_gas:
-            res_nodes_from.extend(
-                [("v_from_m_per_s", "v_gas_from"), ("normfactor_from", "normfactor_from")])
-            res_nodes_to.extend([("v_to_m_per_s", "v_gas_to"), ("normfactor_to", "normfactor_to")])
-            res_mean.extend([("v_mean_m_per_s", "v_gas_mean")])
+            res_nodes_from_hyd.extend([("v_from_m_per_s", "v_gas_from"),
+                                       ("normfactor_from", "normfactor_from")])
+            res_nodes_to_hyd.extend([("v_to_m_per_s", "v_gas_to"),
+                                     ("normfactor_to", "normfactor_to")])
+            res_mean_hyd.extend([("v_mean_m_per_s", "v_gas_mean")])
         else:
-            res_mean.extend([("v_mean_m_per_s", "v_mps")])
+            res_mean_hyd.extend([("v_mean_m_per_s", "v_mps")])
 
         if np.any(cls.get_internal_pipe_number(net) > 1):
             extract_branch_results_with_internals(
-                net, branch_results, cls.table_name(), res_nodes_from, res_nodes_to, res_mean,
-                cls.get_connected_node_type().table_name(), branches_connected)
+                net, branch_results, cls.table_name(), res_nodes_from_hyd, res_nodes_from_ht,
+                res_nodes_to_hyd, res_nodes_to_ht, res_mean_hyd, [],
+                cls.get_connected_node_type().table_name(), mode)
         else:
-            required_results = res_nodes_from + res_nodes_to + res_mean
-            extract_branch_results_without_internals(net, branch_results, required_results,
-                                                     cls.table_name(), branches_connected)
+            required_results_hyd = res_nodes_from_hyd + res_nodes_to_hyd + res_mean_hyd
+            required_results_ht = res_nodes_from_ht + res_nodes_to_ht
+            extract_branch_results_without_internals(
+                net, branch_results, required_results_hyd, required_results_ht, cls.table_name(),
+                mode
+            )
 
     @classmethod
     def get_internal_results(cls, net, pipe):
diff --git a/pandapipes/component_models/pressure_control_component.py b/pandapipes/component_models/pressure_control_component.py
index abd4b247..89d2f220 100644
--- a/pandapipes/component_models/pressure_control_component.py
+++ b/pandapipes/component_models/pressure_control_component.py
@@ -95,14 +95,12 @@ def calculate_temperature_lift(cls, net, branch_component_pit, node_pit):
         branch_component_pit[:, TL] = 0
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         """
         Function that extracts certain results.
 
-        :param nodes_connected:
-        :type nodes_connected:
-        :param branches_connected:
-        :type branches_connected:
+        :param mode:
+        :type mode:
         :param branch_results:
         :type branch_results:
         :param net: The pandapipes network
@@ -111,22 +109,22 @@ def extract_results(cls, net, options, branch_results, nodes_connected, branches
         :type options:
         :return: No Output.
         """
-        required_results = [
-            ("p_from_bar", "p_from"), ("p_to_bar", "p_to"), ("t_from_k", "temp_from"),
-            ("t_to_k", "temp_to"), ("mdot_to_kg_per_s", "mf_to"), ("mdot_from_kg_per_s", "mf_from"),
-            ("vdot_norm_m3_per_s", "vf")
+        required_results_hyd = [
+            ("p_from_bar", "p_from"), ("p_to_bar", "p_to"), ("mdot_from_kg_per_s", "mf_from"),
+            ("mdot_to_kg_per_s", "mf_to"), ("vdot_norm_m3_per_s", "vf")
         ]
+        required_results_ht = [("t_from_k", "temp_from"), ("t_to_k", "temp_to")]
 
         if get_fluid(net).is_gas:
-            required_results.extend([
+            required_results_hyd.extend([
                 ("v_from_m_per_s", "v_gas_from"), ("v_to_m_per_s", "v_gas_to"),
                 ("normfactor_from", "normfactor_from"), ("normfactor_to", "normfactor_to")
             ])
         else:
-            required_results.extend([("v_mean_m_per_s", "v_mps")])
+            required_results_hyd.extend([("v_mean_m_per_s", "v_mps")])
 
-        extract_branch_results_without_internals(net, branch_results, required_results,
-                                                 cls.table_name(), branches_connected)
+        extract_branch_results_without_internals(net, branch_results, required_results_hyd,
+                                                 required_results_ht, cls.table_name(), mode)
 
         res_table = net["res_" + cls.table_name()]
         f, t = get_lookup(net, "branch", "from_to")[cls.table_name()]
diff --git a/pandapipes/component_models/pump_component.py b/pandapipes/component_models/pump_component.py
index 9cee8f91..3de66000 100644
--- a/pandapipes/component_models/pump_component.py
+++ b/pandapipes/component_models/pump_component.py
@@ -111,40 +111,38 @@ def calculate_temperature_lift(cls, net, branch_component_pit, node_pit):
         branch_component_pit[:, TL] = 0
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
+    def extract_results(cls, net, options, branch_results, mode):
         """
         Function that extracts certain results.
 
-        :param nodes_connected:
-        :type nodes_connected:
-        :param branches_connected:
-        :type branches_connected:
         :param branch_results:
         :type branch_results:
         :param net: The pandapipes network
         :type net: pandapipesNet
         :param options:
         :type options:
+        :param mode:
+        :type mode:
         :return: No Output.
         """
         calc_compr_pow = options['calc_compression_power']
 
-        required_results = [
-            ("p_from_bar", "p_from"), ("p_to_bar", "p_to"), ("t_from_k", "temp_from"),
-            ("t_to_k", "temp_to"), ("mdot_to_kg_per_s", "mf_to"), ("mdot_from_kg_per_s", "mf_from"),
-            ("vdot_norm_m3_per_s", "vf"), ("deltap_bar", "pl")
+        required_results_hyd = [
+            ("p_from_bar", "p_from"), ("p_to_bar", "p_to"), ("mdot_to_kg_per_s", "mf_to"),
+            ("mdot_from_kg_per_s", "mf_from"), ("vdot_norm_m3_per_s", "vf"), ("deltap_bar", "pl"),
         ]
+        required_results_ht = [("t_from_k", "temp_from"), ("t_to_k", "temp_to")]
 
         if get_fluid(net).is_gas:
-            required_results.extend([
+            required_results_hyd.extend([
                 ("v_from_m_per_s", "v_gas_from"), ("v_to_m_per_s", "v_gas_to"),
                 ("normfactor_from", "normfactor_from"), ("normfactor_to", "normfactor_to")
             ])
         else:
-            required_results.extend([("v_mean_m_per_s", "v_mps")])
+            required_results_hyd.extend([("v_mean_m_per_s", "v_mps")])
 
-        extract_branch_results_without_internals(net, branch_results, required_results,
-                                                 cls.table_name(), branches_connected)
+        extract_branch_results_without_internals(net, branch_results, required_results_hyd,
+                                                 required_results_ht, cls.table_name(), mode)
 
         if calc_compr_pow:
             f, t = get_lookup(net, "branch", "from_to")[cls.table_name()]
diff --git a/pandapipes/component_models/valve_component.py b/pandapipes/component_models/valve_component.py
index b63173e8..6c32676c 100644
--- a/pandapipes/component_models/valve_component.py
+++ b/pandapipes/component_models/valve_component.py
@@ -5,6 +5,7 @@
 import numpy as np
 from numpy import dtype
 
+from pandapipes.component_models.component_toolbox import standard_branch_wo_internals_result_lookup
 from pandapipes.component_models.abstract_models.branch_wzerolength_models import \
     BranchWZeroLengthComponent
 from pandapipes.component_models.junction_component import Junction
@@ -82,24 +83,11 @@ def get_component_input(cls):
                 ("type", dtype(object))]
 
     @classmethod
-    def extract_results(cls, net, options, branch_results, nodes_connected, branches_connected):
-        required_results = [
-            ("p_from_bar", "p_from"), ("p_to_bar", "p_to"), ("t_from_k", "temp_from"),
-            ("t_to_k", "temp_to"), ("mdot_to_kg_per_s", "mf_to"), ("mdot_from_kg_per_s", "mf_from"),
-            ("vdot_norm_m3_per_s", "vf"), ("lambda", "lambda"), ("reynolds", "reynolds")
-        ]
+    def extract_results(cls, net, options, branch_results, mode):
+        required_results_hyd, required_results_ht = standard_branch_wo_internals_result_lookup(net)
 
-        if get_fluid(net).is_gas:
-            required_results.extend([
-                ("v_from_m_per_s", "v_gas_from"), ("v_to_m_per_s", "v_gas_to"),
-                ("v_mean_m_per_s", "v_gas_mean"), ("normfactor_from", "normfactor_from"),
-                ("normfactor_to", "normfactor_to")
-            ])
-        else:
-            required_results.extend([("v_mean_m_per_s", "v_mps")])
-
-        extract_branch_results_without_internals(net, branch_results, required_results,
-                                                 cls.table_name(), branches_connected)
+        extract_branch_results_without_internals(net, branch_results, required_results_hyd,
+                                                 required_results_ht, cls.table_name(), mode)
 
     @classmethod
     def get_result_table(cls, net):
diff --git a/pandapipes/pf/pipeflow_setup.py b/pandapipes/pf/pipeflow_setup.py
index 4dac0133..51a71952 100644
--- a/pandapipes/pf/pipeflow_setup.py
+++ b/pandapipes/pf/pipeflow_setup.py
@@ -9,9 +9,10 @@
 from scipy.sparse import coo_matrix, csgraph
 
 from pandapipes.idx_branch import FROM_NODE, TO_NODE, branch_cols, \
-    ACTIVE as ACTIVE_BR
+    ACTIVE as ACTIVE_BR, VINIT
 from pandapipes.idx_node import NODE_TYPE, P, NODE_TYPE_T, node_cols, T, ACTIVE as ACTIVE_ND, \
     TABLE_IDX as TABLE_IDX_ND, ELEMENT_IDX as ELEMENT_IDX_ND
+from pandapipes.pf.internals_toolbox import _sum_by_group
 from pandapipes.properties.fluids import get_fluid
 
 try:
@@ -161,8 +162,9 @@ def get_lookup(net, pit_type="node", lookup_type="index"):
     """
     pit_type = pit_type.lower()
     lookup_type = lookup_type.lower()
-    all_lookup_types = ["index", "table", "from_to", "active", "length", "from_to_active",
-                        "index_active"]
+    all_lookup_types = ["index", "table", "from_to", "active_hydraulics", "active_heat_transfer",
+                        "length", "from_to_active_hydraulics", "from_to_active_heat_transfer",
+                        "index_active_hydraulics", "index_active_heat_transfer"]
     if lookup_type not in all_lookup_types:
         type_names = "', '".join(all_lookup_types)
         logger.error("No lookup type '%s' exists. Please choose one of '%s'."
@@ -187,7 +189,7 @@ def set_user_pf_options(net, reset=False, **kwargs):
     :type net: pandapipesNet
     :param reset: Specifies whether the user_pf_options is removed before setting new options
     :type reset: bool, default False
-    :param kwargs: pipeflow options that shall be set, e. g. tol_v = 1e-7
+    :param kwargs: pipeflow options that shall be set, e.g. tol_v = 1e-7
     :return: No output
     """
     if reset or 'user_pf_options' not in net.keys():
@@ -239,7 +241,7 @@ def init_options(net, local_parameters):
                  automatically with respect to the convergence behaviour.
 
         - **mode** (str): "hydraulics" - Define the calculation mode: what shall be calculated - \
-                solely hydraulics ('hydraulic'), solely heat transfer('heat') or both combined \
+                solely hydraulics ('hydraulics'), solely heat transfer('heat') or both combined \
                 ('all').
 
         - **only_update_hydraulic_matrix** (bool): False - If True, the system matrix is not \
@@ -435,7 +437,82 @@ def create_lookups(net):
                        "internal_nodes_lookup": internal_nodes_lookup}
 
 
-def check_connectivity(net, branch_pit, node_pit, check_heat):
+def identify_active_nodes_branches(net, branch_pit, node_pit, hydraulic=True):
+    """
+    Function that creates the connectivity lookup for nodes and branches. If the option \
+    "check_connectivity" is set, a full connectivity check is performed based on a sparse matrix \
+    graph search. Otherwise, only the nodes and branches are identified that are inactive, which \
+    means:\
+      - in case of hydraulics, just use the "ACTIVE" identifier of the respective components\
+      - in case of heat transfer, use the hydraulic result to check which branches are traversed \
+        by the fluid and a simple rule to make sure that active nodes are connected to at least one\
+        traversed branch\
+    The result of this connectivity search is stored in the lookups (e.g. as \
+    net["_lookups"]["node_active_hydraulics"])
+
+    :param net: the pandapipes net for which to identify the connectivity
+    :type net: pandapipes.pandapipesNet
+    :param branch_pit: Internal array with branch entries
+    :type branch_pit: np.array
+    :param node_pit: Internal array with node entries
+    :type node_pit: np.array
+    :param hydraulic: flag for the mode (if True, do the check for the hydraulic simulation, \
+        otherwise for the heat transfer simulation with other considerations)
+    :type hydraulic: bool, default True
+    :return: No output
+    """
+    if hydraulic:
+        # connectivity check for hydraulic simulation
+        if get_net_option(net, "check_connectivity"):
+            nodes_connected, branches_connected = check_connectivity(net, branch_pit, node_pit)
+        else:
+            # if connectivity check is switched off, still consider oos elements
+            nodes_connected = node_pit[:, ACTIVE_ND].astype(np.bool_)
+            branches_connected = branch_pit[:, ACTIVE_BR].astype(np.bool_)
+    else:
+        # connectivity check for heat simulation (needs to consider branches with 0 velocity as
+        # well)
+        if get_net_option(net, "check_connectivity"):
+            # full connectivity check for hydraulic simulation
+            nodes_connected, branches_connected = check_connectivity(net, branch_pit, node_pit,
+                                                                     mode="heat_transfer")
+        else:
+            # if no full connectivity check is performed, all nodes that are not connected to the
+            # rest of the network wrt. flow can be identified by a more performant sum_by_group_call
+            # check for branches that are not traversed (for temperature calculation, this means
+            # that they are "out of service")
+            branches_connected = get_lookup(net, "branch", "active_hydraulics") \
+                                 & branches_connected_flow(branch_pit)
+            fn = branch_pit[:, FROM_NODE].astype(np.int32)
+            tn = branch_pit[:, TO_NODE].astype(np.int32)
+            fn_tn, flow = _sum_by_group(
+                get_net_option(net, "use_numba"), np.concatenate([fn, tn]),
+                np.concatenate([branches_connected, branches_connected]).astype(np.int32)
+            )
+            nodes_connected = np.copy(get_lookup(net, "node", "active_hydraulics"))
+            # set nodes oos that are not connected to any branches with flow > 0 (0.1 is arbitrary
+            # here, any value between 0 and 1 should work, excluding 0 and 1)
+            nodes_connected[fn_tn] = nodes_connected[fn_tn] & (flow > 0.1)
+    mode = "hydraulics" if hydraulic else "heat_transfer"
+    net["_lookups"]["node_active_" + mode] = nodes_connected
+    net["_lookups"]["branch_active_" + mode] = branches_connected
+
+
+def branches_connected_flow(branch_pit):
+    """
+    Simple function to identify branches with flow based on the calculated velocity.
+
+    :param branch_pit: The pandapipes internal table of the network (including hydraulics results)
+    :type branch_pit: np.array
+    :return: branches_connected_flow - lookup array if branch is connected wrt. flow
+    :rtype: np.array
+    """
+    # TODO: is this formulation correct or could there be any caveats?
+    return ~np.isnan(branch_pit[:, VINIT]) \
+        & ~np.isclose(branch_pit[:, VINIT], 0, rtol=1e-10, atol=1e-10)
+
+
+def check_connectivity(net, branch_pit, node_pit, mode="hydraulics"):
     """
     Perform a connectivity check which means that network nodes are identified that don't have any
     connection to an external grid component. Quick overview over the steps of this function:
@@ -461,40 +538,30 @@ def check_connectivity(net, branch_pit, node_pit, check_heat):
     :type branch_pit: np.array
     :param node_pit: Internal array with node entries
     :type node_pit: np.array
-    :param check_heat: Flag which determines whether to also check for connectivity to heat \
-        external grids
-    :type check_heat: bool
-    :return: (nodes_connected_hyd, branches_connected) - Lookups of np.arrays stating which of the
+    :return: (nodes_connected, branches_connected) - Lookups of np.arrays stating which of the
             internal nodes and branches are reachable from any of the hyd_slacks (np mask).
     :rtype: tuple(np.array)
     """
-    active_branch_lookup = branch_pit[:, ACTIVE_BR].astype(bool)
-    active_node_lookup = node_pit[:, ACTIVE_ND].astype(bool)
-    from_nodes = branch_pit[:, FROM_NODE].astype(np.int32)
-    to_nodes = branch_pit[:, TO_NODE].astype(np.int32)
-    hyd_slacks = np.where((node_pit[:, NODE_TYPE] == P) & active_node_lookup)[0]
-    # hyd_slacks = np.where(((node_pit[:, NODE_TYPE] == P) | (node_pit[:, NODE_TYPE] == PC))
-    #                       & active_node_lookup)[0]
-
-    nodes_connected, branches_connected = perform_connectivity_search(
-        net, node_pit, hyd_slacks, from_nodes, to_nodes, active_node_lookup, active_branch_lookup,
-        mode="hydraulics")
-    if not check_heat:
-        return nodes_connected, branches_connected
-
-    heat_slacks = np.where((node_pit[:, NODE_TYPE_T] == T) & nodes_connected)[0]
-    if len(heat_slacks) == len(hyd_slacks) and np.all(heat_slacks == hyd_slacks):
-        return nodes_connected, branches_connected
-
-    nodes_connected, branches_connected = perform_connectivity_search(
-        net, node_pit, heat_slacks, from_nodes, to_nodes, nodes_connected, branches_connected,
-        mode="heat transfer")
-    return nodes_connected, branches_connected
+    if mode == "hydraulics":
+        active_branch_lookup = branch_pit[:, ACTIVE_BR].astype(np.bool_)
+        active_node_lookup = node_pit[:, ACTIVE_ND].astype(np.bool_)
+        slacks = np.where((node_pit[:, NODE_TYPE] == P) & active_node_lookup)[0]
+    else:
+        active_branch_lookup = branches_connected_flow(branch_pit) \
+                               & get_lookup(net, "branch", "active_hydraulics")
+        active_node_lookup = node_pit[:, ACTIVE_ND].astype(np.bool_)\
+                             & get_lookup(net, "node", "active_hydraulics")
+        slacks = np.where((node_pit[:, NODE_TYPE_T] == T) & active_node_lookup)[0]
+
+    return perform_connectivity_search(net, node_pit, branch_pit, slacks, active_node_lookup,
+                                       active_branch_lookup, mode=mode)
 
 
-def perform_connectivity_search(net, node_pit, slack_nodes, from_nodes, to_nodes,
+def perform_connectivity_search(net, node_pit, branch_pit, slack_nodes,
                                 active_node_lookup, active_branch_lookup, mode="hydraulics"):
     len_nodes = len(node_pit)
+    from_nodes = branch_pit[:, FROM_NODE].astype(np.int32)
+    to_nodes = branch_pit[:, TO_NODE].astype(np.int32)
     nobranch = np.sum(active_branch_lookup)
     active_from_nodes = from_nodes[active_branch_lookup]
     active_to_nodes = to_nodes[active_branch_lookup]
@@ -569,7 +636,7 @@ def get_table_index_list(net, pit_array, pit_indices, pit_type="node"):
             for tbl in tables]
 
 
-def reduce_pit(net, node_pit, branch_pit, nodes_connected, branches_connected):
+def reduce_pit(net, node_pit, branch_pit, mode="hydraulics"):
     """
     Create an internal ("active") pit with all nodes and branches that are actually in_service. This
     is also done for different lookups (e.g. the from_to indices for this pit and the node index
@@ -582,44 +649,46 @@ def reduce_pit(net, node_pit, branch_pit, nodes_connected, branches_connected):
     :type node_pit: np.array
     :param branch_pit: The internal structure branch array
     :type branch_pit: np.array
-    :param nodes_connected: A mask array stating which nodes are actually connected to the rest of\
-            the net
-    :type nodes_connected: np.array
-    :param branches_connected: A mask array stating which branches are actually connected to the \
-             rest of the net
-    :type branches_connected: np.array
+    :param mode: the mode of the calculation (either "hydraulics" or "heat_transfer") for storing /\
+        retrieving correct lookups
+    :type mode: str, default "hydraulics"
     :return: No output
     """
     active_pit = dict()
     els = dict()
     reduced_node_lookup = None
+    nodes_connected = get_lookup(net, "node", "active_" + mode)
+    branches_connected = get_lookup(net, "branch", "active_" + mode)
     if np.alltrue(nodes_connected):
-        net["_lookups"]["node_from_to_active"] = copy.deepcopy(get_lookup(net, "node", "from_to"))
-        net["_lookups"]["node_index_active"] = copy.deepcopy(get_lookup(net, "node", "index"))
+        net["_lookups"]["node_from_to_active_" + mode] = copy.deepcopy(
+            get_lookup(net, "node", "from_to"))
+        net["_lookups"]["node_index_active_" + mode] = copy.deepcopy(
+            get_lookup(net, "node", "index"))
         active_pit["node"] = np.copy(node_pit)
     else:
         active_pit["node"] = np.copy(node_pit[nodes_connected, :])
         reduced_node_lookup = np.cumsum(nodes_connected) - 1
         node_idx_lookup = get_lookup(net, "node", "index")
-        net["_lookups"]["node_index_active"] = {
+        net["_lookups"]["node_index_active_" + mode] = {
             tbl: reduced_node_lookup[idx_lookup[idx_lookup != -1]]
             for tbl, idx_lookup in node_idx_lookup.items()}
         els["node"] = nodes_connected
     if np.alltrue(branches_connected):
-        net["_lookups"]["branch_from_to_active"] = copy.deepcopy(get_lookup(net, "branch",
-                                                                            "from_to"))
+        net["_lookups"]["branch_from_to_active_" + mode] = copy.deepcopy(
+            get_lookup(net, "branch", "from_to"))
         active_pit["branch"] = np.copy(branch_pit)
-        net["_lookups"]["branch_index_active"] = copy.deepcopy(get_lookup(net, "branch", "index"))
+        net["_lookups"]["branch_index_active_" + mode] = copy.deepcopy(
+            get_lookup(net, "branch", "index"))
     else:
         active_pit["branch"] = np.copy(branch_pit[branches_connected, :])
         branch_idx_lookup = get_lookup(net, "branch", "index")
         if len(branch_idx_lookup):
             reduced_branch_lookup = np.cumsum(branches_connected) - 1
-            net["_lookups"]["branch_index_active"] = {
+            net["_lookups"]["branch_index_active_" + mode] = {
                 tbl: reduced_branch_lookup[idx_lookup[idx_lookup != -1]]
                 for tbl, idx_lookup in branch_idx_lookup.items()}
         else:
-            net["_lookups"]["branch_index_active"] = dict()
+            net["_lookups"]["branch_index_active_" + mode] = dict()
         els["branch"] = branches_connected
     if reduced_node_lookup is not None:
         active_pit["branch"][:, FROM_NODE] = reduced_node_lookup[
@@ -627,8 +696,6 @@ def reduce_pit(net, node_pit, branch_pit, nodes_connected, branches_connected):
         active_pit["branch"][:, TO_NODE] = reduced_node_lookup[
             branch_pit[branches_connected, TO_NODE].astype(np.int32)]
     net["_active_pit"] = active_pit
-    net["_lookups"]["node_active"] = nodes_connected
-    net["_lookups"]["branch_active"] = branches_connected
 
     for el, connected_els in els.items():
         ft_lookup = get_lookup(net, el, "from_to")
@@ -639,4 +706,4 @@ def reduce_pit(net, node_pit, branch_pit, nodes_connected, branches_connected):
         for table, (_, _, len_new) in sorted(aux_lookup.items(), key=lambda x: x[1][0]):
             from_to_active_lookup[table] = (count, count + len_new)
             count += len_new
-        net["_lookups"]["%s_from_to_active" % el] = from_to_active_lookup
+        net["_lookups"]["%s_from_to_active_%s" % (el, mode)] = from_to_active_lookup
diff --git a/pandapipes/pf/result_extraction.py b/pandapipes/pf/result_extraction.py
index c453c9f5..1782a9f4 100644
--- a/pandapipes/pf/result_extraction.py
+++ b/pandapipes/pf/result_extraction.py
@@ -2,7 +2,7 @@
 
 from pandapipes.constants import NORMAL_PRESSURE, NORMAL_TEMPERATURE
 from pandapipes.idx_branch import ELEMENT_IDX, FROM_NODE, TO_NODE, LOAD_VEC_NODES, VINIT, RE, \
-    LAMBDA, TINIT, FROM_NODE_T, TO_NODE_T, PL
+    LAMBDA, TINIT, FROM_NODE_T, TO_NODE_T, PL, T_OUT
 from pandapipes.idx_node import TABLE_IDX as TABLE_IDX_NODE, PINIT, PAMB, TINIT as TINIT_NODE
 from pandapipes.pf.internals_toolbox import _sum_by_group
 from pandapipes.pf.pipeflow_setup import get_table_number, get_lookup, get_net_option
@@ -14,13 +14,15 @@
     from pandapower.pf.no_numba import jit
 
 
-def extract_all_results(net, nodes_connected, branches_connected):
+def extract_all_results(net, calculation_mode):
     """
     Extract results from branch pit and node pit and write them to the different tables of the net,\
     as defined by the component models.
 
     :param net: pandapipes net for which to extract results into net.res_xy
     :type net: pandapipesNet
+    :param net: mode of the simulation (e.g. "hydraulics" or "heat" or "all")
+    :type net: str
     :return: No output
 
     """
@@ -49,8 +51,7 @@ def extract_all_results(net, nodes_connected, branches_connected):
         }
         branch_results.update(gas_branch_results)
     for comp in net['component_list']:
-        comp.extract_results(net, net["_options"], branch_results, nodes_connected,
-                             branches_connected)
+        comp.extract_results(net, net["_options"], branch_results, calculation_mode)
 
 
 def get_basic_branch_results(net, branch_pit, node_pit):
@@ -141,8 +142,11 @@ def get_gas_vel_numba(branch_pit, comp_from, comp_to, comp_mean, p_abs_from, p_a
     return v_gas_from, v_gas_to, v_gas_mean, normfactor_from, normfactor_to, normfactor_mean
 
 
-def extract_branch_results_with_internals(net, branch_results, table_name, res_nodes_from,
-                                          res_nodes_to, res_mean, node_name, branches_connected):
+def extract_branch_results_with_internals(net, branch_results, table_name,
+                                          res_nodes_from_hydraulics, res_nodes_from_heat,
+                                          res_nodes_to_hydraulics, res_nodes_to_heat,
+                                          res_mean_hydraulics, res_mean_heat, node_name,
+                                          simulation_mode):
     # the result table to write results to
     res_table = net["res_" + table_name]
 
@@ -155,101 +159,137 @@ def extract_branch_results_with_internals(net, branch_results, table_name, res_n
     # respective table), the placement of the indices mus be known to allocate the values correctly
     placement_table = np.argsort(net[table_name].index.values)
     idx_pit = branch_pit[f:t, ELEMENT_IDX]
-    comp_connected = branches_connected[f:t]
 
     node_pit = net["_pit"]["node"]
 
     # the id of the external node table inside the node_pit (mostly this is "junction": 0)
     ext_node_tbl_idx = get_table_number(get_lookup(net, "node", "table"), node_name)
 
-    if len(res_nodes_from) > 0:
-        # results that relate to the from_node --> in case of many internal nodes, only the single
-        # from_node that is the exterior node (e.g. junction vs. internal pipe_node) result has to
-        # be extracted from the node_pit
-        from_nodes = branch_results["from_nodes"][f:t]
-        from_nodes_external = node_pit[from_nodes, TABLE_IDX_NODE] == ext_node_tbl_idx
-        considered = from_nodes_external & comp_connected
-        external_active = comp_connected[from_nodes_external]
-        for res_name, entry in res_nodes_from:
-            res_table[res_name].values[external_active] = branch_results[entry][f:t][considered]
-
-    if len(res_nodes_to) > 0:
-        # results that relate to the to_node --> in case of many internal nodes, only the single
-        # to_node that is the exterior node (e.g. junction vs. internal pipe_node) result has to
-        # be extracted from the node_pit
-        to_nodes = branch_results["to_nodes"][f:t]
-        to_nodes_external = node_pit[to_nodes, TABLE_IDX_NODE] == ext_node_tbl_idx
-        considered = to_nodes_external & comp_connected
-        external_active = comp_connected[to_nodes_external]
-        for res_name, entry in res_nodes_to:
-            res_table[res_name].values[external_active] = branch_results[entry][f:t][considered]
-
-    if len(res_mean) > 0:
-        # results that relate to the whole branch and shall be averaged (by summing up all values
-        # and dividing by number of internal sections)
-        use_numba = get_net_option(net, "use_numba")
-        res = _sum_by_group(use_numba, idx_pit, np.ones_like(idx_pit),
-                            comp_connected.astype(np.int32),
-                            *[branch_results[rn[1]][f:t] for rn in res_mean])
-        connected_ind = res[2] > 0.99
-        num_internals = res[1][connected_ind]
-        # hint: idx_pit[placement_table] should result in the indices as ordered in the table
-        placement_table = placement_table[connected_ind]
-
-        for i, (res_name, entry) in enumerate(res_mean):
-            res_table[res_name].values[placement_table] = res[i + 3][connected_ind] / num_internals
-
-
-def extract_branch_results_without_internals(net, branch_results, required_results, table_name,
-                                             branches_connected):
+    for (result_mode, res_nodes_from, res_nodes_to, res_mean) in [
+        ("hydraulics", res_nodes_from_hydraulics, res_nodes_to_hydraulics, res_mean_hydraulics),
+        ("heat", res_nodes_from_heat, res_nodes_to_heat, res_mean_heat)
+    ]:
+        if result_mode == "hydraulics" and simulation_mode == "heat":
+            continue
+        lookup_name = "hydraulics"
+        if result_mode == "heat" and simulation_mode in ["heat", "all"]:
+            lookup_name = "heat_transfer"
+        comp_connected = get_lookup(net, "branch", "active_" + lookup_name)[f:t]
+        for (res_ext, node_name) in ((res_nodes_from, "from_nodes"), (res_nodes_to, "to_nodes")):
+            if len(res_ext) == 0:
+                continue
+            # results that relate to the from_node --> in case of many internal nodes, only the
+            # single from_node that is the exterior node (e.g. junction vs. internal pipe_node)
+            # result has to be extracted from the node_pit
+            end_nodes = branch_results[node_name][f:t]
+            end_nodes_external = node_pit[end_nodes, TABLE_IDX_NODE] == ext_node_tbl_idx
+            considered = end_nodes_external & comp_connected
+            external_active = comp_connected[end_nodes_external]
+            for res_name, entry in res_ext:
+                res_table[res_name].values[external_active] = branch_results[entry][f:t][considered]
+        if len(res_mean) > 0:
+            # results that relate to the whole branch and shall be averaged (by summing up all
+            # values and dividing by number of internal sections)
+            use_numba = get_net_option(net, "use_numba")
+            res = _sum_by_group(use_numba, idx_pit, np.ones_like(idx_pit),
+                                comp_connected.astype(np.int32),
+                                *[branch_results[rn[1]][f:t] for rn in res_mean])
+            connected_ind = res[2] > 0.99
+            num_internals = res[1][connected_ind]
+
+            # hint: idx_pit[placement_table] should result in the indices as ordered in the table
+            pt = placement_table[connected_ind]
+
+            for i, (res_name, entry) in enumerate(res_mean_hydraulics):
+                res_table[res_name].values[pt] = res[i + 3][connected_ind] / num_internals
+
+
+def extract_branch_results_without_internals(net, branch_results, required_results_hydraulic,
+                                             required_results_heat, table_name, simulation_mode):
+    """
+    Extract the results from the branch result array derived from the pit to the result table of the
+    net (only for branch components without internal nodes). Here, we need to consider which results
+    exist for hydraulic calculation and for heat transfer calculation (wrt. connectivity).
+
+    :param net: The pandapipes net that the internal structure belongs to
+    :type net: pandapipesNet
+    :param branch_results: Important branch results from the internal pit structure
+    :type branch_results: dict[np.ndarray]
+    :param required_results_hydraulic: The entries that should be extracted for the respective \
+        component for hydraulic calculation
+    :type required_results_hydraulic: list[tuple]
+    :param required_results_heat: The entries that should be extracted for the respective \
+        component for heat transfer calculation
+    :type required_results_heat: list[tuple]
+    :param table_name: The name of the table that the results should be written to
+    :type table_name: str
+    :param simulation_mode: simulation mode (e.g. "hydraulics", "heat", "all"); defines whether results from \
+        hydraulic or temperature calculation are transferred
+    :type simulation_mode: str
+    :return: No output
+    :rtype: None
+    """
     res_table = net["res_" + table_name]
     f, t = get_lookup(net, "branch", "from_to")[table_name]
-    comp_connected = branches_connected[f:t]
-
-    for res_name, entry in required_results:
-        res_table[res_name].values[:][comp_connected] = branch_results[entry][f:t][comp_connected]
-
 
-def extract_results_active_pit(net, node_pit, branch_pit, nodes_connected, branches_connected):
+    # extract hydraulic results
+    if simulation_mode in ["hydraulics", "all"]:
+        # lookup for connected branch elements (hydraulic results)
+        comp_connected_hyd = get_lookup(net, "branch", "active_hydraulics")[f:t]
+        for res_name, entry in required_results_hydraulic:
+            res_table[res_name].values[:][comp_connected_hyd] = \
+                branch_results[entry][f:t][comp_connected_hyd]
+        if simulation_mode == "hydraulics":
+            for res_name, entry in required_results_heat:
+                res_table[res_name].values[:][comp_connected_hyd] = \
+                    branch_results[entry][f:t][comp_connected_hyd]
+
+    # extract heat transfer results
+    if simulation_mode in ["heat", "all"]:
+        # lookup for connected branch elements (heat transfer results)
+        comp_connected_ht = get_lookup(net, "branch", "active_heat_transfer")[f:t]
+        for res_name, entry in required_results_heat:
+            res_table[res_name].values[:][comp_connected_ht] = \
+                branch_results[entry][f:t][comp_connected_ht]
+
+
+def extract_results_active_pit(net,  mode="hydraulics"):
     """
     Extract the pipeflow results from the internal pit structure ("_active_pit") to the general pit
     structure.
 
     :param net: The pandapipes net that the internal structure belongs to
     :type net: pandapipesNet
-    :param node_pit: The internal structure node array
-    :type node_pit: np.array
-    :param branch_pit: The internal structure branch array
-    :type branch_pit: np.array
-    :param nodes_connected: A mask array stating which nodes are actually connected to the rest of\
-            the net
-    :type nodes_connected: np.array
-    :param branches_connected: A mask array stating which branches are actually connected to the \
-             rest of the net
-    :type branches_connected: np.array
+    :param mode: defines whether results from hydraulic or temperature calculation are transferred
+    :type mode: str, default "hydraulics"
     :return: No output
 
     """
-    all_nodes_connected = np.alltrue(nodes_connected)
-    if not all_nodes_connected:
-        node_pit[~nodes_connected, PINIT] = np.NaN
-        node_pit[nodes_connected, :] = net["_active_pit"]["node"]
-        cols_br = np.array([i for i in range(branch_pit.shape[1])
-                            if i not in [FROM_NODE, TO_NODE, FROM_NODE_T, TO_NODE_T]])
-    else:
-        net["_pit"]["node"] = np.copy(net["_active_pit"]["node"])
-        cols_br = None
-
-    if not np.alltrue(branches_connected):
-        branch_pit[~branches_connected, VINIT] = np.NaN
-        rows_active_br = np.where(branches_connected)[0]
-        if all_nodes_connected:
-            branch_pit[rows_active_br, :] = net["_active_pit"]["branch"][:, :]
-        else:
-            branch_pit[rows_active_br[:, np.newaxis], cols_br[np.newaxis, :]] = \
-                net["_active_pit"]["branch"][:, cols_br]
-    else:
-        if all_nodes_connected:
-            net["_pit"]["branch"] = np.copy(net["_active_pit"]["branch"])
-        else:
-            net["_pit"]["branch"][:, cols_br] = net["_active_pit"]["branch"][:, cols_br]
+    nodes_connected = get_lookup(net, "node", "active_" + mode)
+    branches_connected = get_lookup(net, "branch", "active_" + mode)
+    result_node_col = PINIT if mode == "hydraulics" else TINIT_NODE
+    not_affected_node_col = TINIT_NODE if mode == "hydraulics" else PINIT
+    copied_node_cols = np.array([i for i in range(net["_pit"]["node"].shape[1])
+                                 if i not in [not_affected_node_col]])
+    rows_nodes = np.arange(net["_pit"]["node"].shape[0])[nodes_connected]
+
+    result_branch_col = VINIT if mode == "hydraulics" else T_OUT
+    not_affected_branch_col = T_OUT if mode == "hydraulics" else VINIT
+    copied_branch_cols = np.array([i for i in range(net["_pit"]["branch"].shape[1])
+                                   if i not in [FROM_NODE, TO_NODE, FROM_NODE_T, TO_NODE_T,
+                                                not_affected_branch_col]])
+    rows_branches = np.arange(net["_pit"]["branch"].shape[0])[branches_connected]
+
+    net["_pit"]["node"][~nodes_connected, result_node_col] = np.NaN
+    net["_pit"]["node"][rows_nodes[:, np.newaxis], copied_node_cols[np.newaxis, :]] = \
+        net["_active_pit"]["node"][:, copied_node_cols]
+    net["_pit"]["branch"][~branches_connected, result_branch_col] = np.NaN
+    net["_pit"]["branch"][rows_branches[:, np.newaxis], copied_branch_cols[np.newaxis, :]] = \
+        net["_active_pit"]["branch"][:, copied_branch_cols]
+
+
+def consider_heat(mode, results=None):
+    consider_ = mode in ["heat", "all"]
+    if results is None:
+        return consider_
+    return consider_ and any(r[2] for r in results)
diff --git a/pandapipes/pipeflow.py b/pandapipes/pipeflow.py
index 9f9544aa..1dc89367 100644
--- a/pandapipes/pipeflow.py
+++ b/pandapipes/pipeflow.py
@@ -7,15 +7,14 @@
 from pandapower.auxiliary import ppException
 from scipy.sparse.linalg import spsolve
 
-from pandapipes.idx_branch import ACTIVE as ACTIVE_BR, FROM_NODE, TO_NODE, FROM_NODE_T, \
-    TO_NODE_T, VINIT, T_OUT, VINIT_T
-from pandapipes.idx_node import PINIT, TINIT, ACTIVE as ACTIVE_ND
+from pandapipes.idx_branch import FROM_NODE, TO_NODE, FROM_NODE_T, TO_NODE_T, VINIT, T_OUT, VINIT_T
+from pandapipes.idx_node import PINIT, TINIT
 from pandapipes.pf.build_system_matrix import build_system_matrix
 from pandapipes.pf.derivative_calculation import calculate_derivatives_hydraulic
 from pandapipes.pf.pipeflow_setup import get_net_option, get_net_options, set_net_option, \
     init_options, create_internal_results, write_internal_results, get_lookup, create_lookups, \
-    initialize_pit, check_connectivity, reduce_pit, \
-    set_user_pf_options, init_all_result_tables
+    initialize_pit, reduce_pit, set_user_pf_options, init_all_result_tables, \
+    identify_active_nodes_branches
 from pandapipes.pf.result_extraction import extract_all_results, extract_results_active_pit
 
 try:
@@ -79,37 +78,36 @@ def pipeflow(net, sol_vec=None, **kwargs):
     calculate_hydraulics = calculation_mode in ["hydraulics", "all"]
     calculate_heat = calculation_mode in ["heat", "all"]
 
-    if get_net_option(net, "check_connectivity"):
-        nodes_connected, branches_connected = check_connectivity(
-            net, branch_pit, node_pit, check_heat=calculate_heat)
-    else:
-        nodes_connected = node_pit[:, ACTIVE_ND].astype(bool)
-        branches_connected = branch_pit[:, ACTIVE_BR].astype(bool)
-
-    reduce_pit(net, node_pit, branch_pit, nodes_connected, branches_connected)
+    identify_active_nodes_branches(net, branch_pit, node_pit)
 
-    if calculation_mode == "heat" and not net.user_pf_options["hyd_flag"]:
-        raise UserWarning("Converged flag not set. Make sure that hydraulic calculation results "
-                          "are available.")
-    elif calculation_mode == "heat" and net.user_pf_options["hyd_flag"]:
-        net["_active_pit"]["node"][:, PINIT] = sol_vec[:len(node_pit)]
-        net["_active_pit"]["branch"][:, VINIT] = sol_vec[len(node_pit):]
+    if calculation_mode == "heat":
+        if not net.user_pf_options["hyd_flag"]:
+            raise UserWarning("Converged flag not set. Make sure that hydraulic calculation "
+                              "results are available.")
+        else:
+            net["_pit"]["node"][:, PINIT] = sol_vec[:len(node_pit)]
+            net["_pit"]["branch"][:, VINIT] = sol_vec[len(node_pit):]
 
     if calculate_hydraulics:
+        reduce_pit(net, node_pit, branch_pit, mode="hydraulics")
         converged, _ = hydraulics(net)
         if not converged:
             raise PipeflowNotConverged("The hydraulic calculation did not converge to a solution.")
+        extract_results_active_pit(net, mode="hydraulics")
 
     if calculate_heat:
+        node_pit, branch_pit = net["_pit"]["node"], net["_pit"]["branch"]
+        identify_active_nodes_branches(net, branch_pit, node_pit, False)
+        reduce_pit(net, node_pit, branch_pit, mode="heat_transfer")
         converged, _ = heat_transfer(net)
         if not converged:
             raise PipeflowNotConverged("The heat transfer calculation did not converge to a "
                                        "solution.")
+        extract_results_active_pit(net, mode="heat_transfer")
     elif not calculate_hydraulics:
         raise UserWarning("No proper calculation mode chosen.")
 
-    extract_results_active_pit(net, node_pit, branch_pit, nodes_connected, branches_connected)
-    extract_all_results(net, nodes_connected, branches_connected)
+    extract_all_results(net, calculation_mode)
 
 
 def hydraulics(net):
@@ -191,7 +189,7 @@ def heat_transfer(net):
         error_t_out.append(linalg.norm(delta_t_out) / (len(delta_t_out)))
 
         finalize_iteration(net, niter, error_t, error_t_out, residual_norm, nonlinear_method, tol_t,
-                           tol_t, tol_res, t_init_old, t_out_old, hyraulic_mode=True)
+                           tol_t, tol_res, t_init_old, t_out_old, hydraulic_mode=False)
         logger.debug("F: %s" % epsilon.round(4))
         logger.debug("T_init_: %s" % t_init.round(4))
         logger.debug("T_out_: %s" % t_out.round(4))
@@ -222,7 +220,7 @@ def solve_hydraulics(net):
     branch_pit = net["_active_pit"]["branch"]
     node_pit = net["_active_pit"]["node"]
 
-    branch_lookups = get_lookup(net, "branch", "from_to_active")
+    branch_lookups = get_lookup(net, "branch", "from_to_active_hydraulics")
     for comp in net['component_list']:
         comp.adaption_before_derivatives_hydraulic(
             net, branch_pit, node_pit, branch_lookups, options)
@@ -258,7 +256,7 @@ def solve_temperature(net):
     options = net["_options"]
     branch_pit = net["_active_pit"]["branch"]
     node_pit = net["_active_pit"]["node"]
-    branch_lookups = get_lookup(net, "branch", "from_to_active")
+    branch_lookups = get_lookup(net, "branch", "from_to_active_heat_transfer")
 
     # Negative velocity values are turned to positive ones (including exchange of from_node and
     # to_node for temperature calculation
@@ -314,8 +312,8 @@ def set_damping_factor(net, niter, error):
 
 
 def finalize_iteration(net, niter, error_1, error_2, residual_norm, nonlinear_method, tol_1, tol_2,
-                       tol_res, vals_1_old, vals_2_old, hyraulic_mode=True):
-    col1, col2 = (PINIT, VINIT) if hyraulic_mode else (TINIT, T_OUT)
+                       tol_res, vals_1_old, vals_2_old, hydraulic_mode=True):
+    col1, col2 = (PINIT, VINIT) if hydraulic_mode else (TINIT, T_OUT)
 
     # Control of damping factor
     if nonlinear_method == "automatic":
@@ -335,7 +333,7 @@ def finalize_iteration(net, niter, error_1, error_2, residual_norm, nonlinear_me
         elif get_net_option(net, "alpha") == 1:
             set_net_option(net, "converged", True)
 
-    if hyraulic_mode:
+    if hydraulic_mode:
         logger.debug("errorv: %s" % error_1[niter])
         logger.debug("errorp: %s" % error_2[niter])
         logger.debug("alpha: %s" % get_net_option(net, "alpha"))
diff --git a/pandapipes/test/pipeflow_internals/test_inservice.py b/pandapipes/test/pipeflow_internals/test_inservice.py
index 56beecf2..1db30680 100644
--- a/pandapipes/test/pipeflow_internals/test_inservice.py
+++ b/pandapipes/test/pipeflow_internals/test_inservice.py
@@ -87,7 +87,7 @@ def complex_heat_connectivity_grid():
     pandapipes.create_pipe_from_parameters(net, j3, j5, 0.1, 0.1, alpha_w_per_m2k=5,
                                            in_service=False, index=7)
     pandapipes.create_pipe_from_parameters(net, j6, j7, 0.1, 0.1, alpha_w_per_m2k=5, index=9)
-    pandapipes.create_pipe_from_parameters(net, j5, 8, 0.1, 0.1, alpha_w_per_m2k=5,
+    pandapipes.create_pipe_from_parameters(net, j5, j8, 0.1, 0.1, alpha_w_per_m2k=5,
                                            in_service=False, index=8)
     pandapipes.create_pipe_from_parameters(net, j8, j10, 0.1, 0.1, alpha_w_per_m2k=5, index=1)
     pandapipes.create_pipe_from_parameters(net, j9, j10, 0.1, 0.1, alpha_w_per_m2k=5, index=2)
@@ -143,7 +143,7 @@ def test_inservice_gas(create_test_net, use_numba):
 
     assert np.all(np.isnan(net.res_pipe.loc[~net.pipe.in_service, :].values))
     assert np.all(np.isnan(net.res_valve.loc[~net.valve.opened, :].values))
-    assert np.all(np.isnan(net.res_junction.loc[~net.junction.in_service, :].values))
+    assert np.all(np.isnan(net.res_junction.p_bar.loc[~net.junction.in_service].values))
 
     oos_sinks = np.isin(net.sink.junction.values, net.junction.index[~net.junction.in_service]) \
         | ~net.sink.in_service.values
@@ -176,7 +176,7 @@ def test_inservice_water(create_test_net, use_numba):
 
     assert np.all(np.isnan(net.res_pipe.loc[~net.pipe.in_service, :].values))
     assert np.all(np.isnan(net.res_valve.loc[~net.valve.opened, :].values))
-    assert np.all(np.isnan(net.res_junction.loc[~net.junction.in_service, :].values))
+    assert np.all(np.isnan(net.res_junction.p_bar.loc[~net.junction.in_service].values))
 
     oos_sinks = np.isin(net.sink.junction.values, net.junction.index[~net.junction.in_service]) \
         | ~net.sink.in_service.values
@@ -210,7 +210,7 @@ def test_connectivity_hydraulic(create_test_net, use_numba):
     pandapipes.pipeflow(net, iter=100, tol_p=1e-7, tol_v=1e-7, friction_model="nikuradse",
                         use_numba=use_numba)
 
-    assert np.all(np.isnan(net.res_junction.loc[[2, 5, 6], :].values))
+    assert np.all(np.isnan(net.res_junction.p_bar.loc[[2, 5, 6]].values))
     assert np.all(np.isnan(net.res_pipe.loc[[1, 2, 3], :].values))
     assert not np.any(np.isnan(net.res_junction.loc[[0, 1, 3, 4], :].values))
     assert not np.any(np.isnan(net.res_pipe.loc[[0, 4],
@@ -222,8 +222,8 @@ def test_connectivity_hydraulic(create_test_net, use_numba):
     assert np.allclose(net.res_ext_grid.mdot_kg_per_s.sum(), -net.res_sink.mdot_kg_per_s.sum(),
                        rtol=1e-10, atol=0)
 
-    active_branches = get_lookup(net, "branch", "active")
-    active_nodes = get_lookup(net, "node", "active")
+    active_branches = get_lookup(net, "branch", "active_hydraulics")
+    active_nodes = get_lookup(net, "node", "active_hydraulics")
 
     assert np.all(active_nodes == np.array([True, True, False, True, True, False, False, False,
                                             False, True]))
@@ -257,16 +257,20 @@ def test_connectivity_hydraulic2(create_test_net, use_numba):
     pandapipes.pipeflow(net, iter=100, tol_p=1e-7, tol_v=1e-7, friction_model="nikuradse",
                         use_numba=use_numba)
 
-    active_branches = get_lookup(net, "branch", "active")
-    active_nodes = get_lookup(net, "node", "active")
+    active_branches = get_lookup(net, "branch", "active_hydraulics")
+    active_nodes = get_lookup(net, "node", "active_hydraulics")
 
     assert np.all(active_nodes == np.array([True, True, True, True, True, True, True, False,
                                             False, True, False, False, True, True, True]))
     assert np.all(active_branches)
 
-    assert not np.all(np.isnan(net.res_junction.loc[[0, 1, 2, 3, 4, 5, 9], :].values))
+    assert not np.all(np.isnan(net.res_junction.p_bar.loc[[0, 1, 2, 3, 4, 5, 9]].values))
     assert not np.all(np.isnan(net.res_pipe.values))
-    assert np.any(np.isnan(net.res_junction.loc[[7, 8, 10, 11], :].values))
+    assert np.all(np.isnan(net.res_junction.p_bar.loc[[7, 8, 10, 11]].values))
+
+    with pytest.raises(PipeflowNotConverged):
+        pandapipes.pipeflow(net, iter=100, tol_p=1e-7, tol_v=1e-7, friction_model="nikuradse",
+                            use_numba=use_numba, check_connectivity=False)
 
 
 @pytest.mark.parametrize("use_numba", [True, False])
@@ -274,21 +278,30 @@ def test_connectivity_heat1(complex_heat_connectivity_grid, use_numba):
     net = copy.deepcopy(complex_heat_connectivity_grid)
     pandapipes.pipeflow(net, mode="all", check_connectivity=True, use_numba=use_numba)
 
-    assert set(net.res_junction.loc[net.res_junction.p_bar.notnull()].index) == {8, 9, 10}
-    assert set(net.res_junction.loc[net.res_junction.p_bar.isnull()].index) \
-           == set(net.junction.index) - {8, 9, 10}
-    assert set(net.res_pipe.loc[net.res_pipe.v_mean_m_per_s.notnull()].index) == {1, 2}
-    assert set(net.res_pipe.loc[net.res_pipe.v_mean_m_per_s.isnull()].index) \
-           == set(net.pipe.index) - {1, 2}
+    oos_juncs_hyd = {4, 5, 6, 7}
+    oos_pipe_hyd = {5, 7, 8, 9}
+    oos_sink_hyd = {4, 5}
+    oos_source_hyd = {7}
+
+    assert set(net.res_junction.loc[net.res_junction.p_bar.notnull()].index) == \
+           set(net.junction.index) - oos_juncs_hyd
+    assert set(net.res_junction.loc[net.res_junction.p_bar.isnull()].index) == oos_juncs_hyd
+
+    assert set(net.res_pipe.loc[net.res_pipe.v_mean_m_per_s.notnull()].index) == \
+           set(net.pipe.index) - oos_pipe_hyd
+    assert set(net.res_pipe.loc[net.res_pipe.v_mean_m_per_s.isnull()].index) == oos_pipe_hyd
+
     assert set(net.res_valve.loc[net.res_valve.v_mean_m_per_s.notnull()].index) == set()
     assert set(net.res_valve.loc[net.res_valve.v_mean_m_per_s.isnull()].index) \
            == set(net.valve.index)
-    assert set(net.res_sink.loc[net.res_sink.mdot_kg_per_s.isnull()].index) == {3, 4, 5}
+
+    assert set(net.res_sink.loc[net.res_sink.mdot_kg_per_s.isnull()].index) == oos_sink_hyd
     assert set(net.res_sink.loc[net.res_sink.mdot_kg_per_s.notnull()].index) == \
-           set(net.sink.index) - {3, 4, 5}
-    assert set(net.res_source.loc[net.res_source.mdot_kg_per_s.isnull()].index) == \
-           set(net.source.index)
-    assert set(net.res_source.loc[net.res_source.mdot_kg_per_s.notnull()].index) == set()
+           set(net.sink.index) - oos_sink_hyd
+
+    assert set(net.res_source.loc[net.res_source.mdot_kg_per_s.isnull()].index) == oos_source_hyd
+    assert (set(net.res_source.loc[net.res_source.mdot_kg_per_s.notnull()].index) ==
+            set(net.source.index) - oos_source_hyd)
 
     assert np.allclose(net.res_ext_grid.mdot_kg_per_s.sum(),
                        -net.res_sink.mdot_kg_per_s.sum() + net.res_source.mdot_kg_per_s.sum(),
@@ -352,6 +365,42 @@ def test_connectivity_heat3(complex_heat_connectivity_grid, use_numba):
                        rtol=1e-10, atol=0)
 
 
+@pytest.mark.parametrize("use_numba", [True, False])
+def test_connectivity_heat4(complex_heat_connectivity_grid, use_numba):
+    net = copy.deepcopy(complex_heat_connectivity_grid)
+
+    net.pipe.in_service.loc[[7, 8]] = True
+    j_new = pandapipes.create_junction(net, 1, 320.15)
+    pandapipes.create_pipe_from_parameters(net, 8, j_new, 0.1, 0.1, alpha_w_per_m2k=5)
+
+    net2 = copy.deepcopy(net)
+
+    pandapipes.pipeflow(net, mode="all", check_connectivity=True, use_numba=use_numba)
+    pandapipes.pipeflow(net2, mode="all", check_connectivity=False, use_numba=use_numba)
+
+    assert pandapipes.nets_equal(net, net2, check_only_results=True)
+
+
+@pytest.mark.parametrize("use_numba", [True, False])
+def test_connectivity_heat5(complex_heat_connectivity_grid, use_numba):
+    net = copy.deepcopy(complex_heat_connectivity_grid)
+    net.pipe.in_service.loc[[7, 8]] = True
+
+    j_from, j_to = pandapipes.create_junctions(net, 2, 1, 320.15)
+
+    pandapipes.create_pipe_from_parameters(net, j_from, j_to, 0.1, 0.1, alpha_w_per_m2k=5)
+    pandapipes.create_sink(net, j_to, 0.1)
+    pandapipes.create_ext_grid(net, j_from, 1, 320.15)
+
+    net.ext_grid.loc[2, 'in_service'] = False
+    net.ext_grid.loc[1, 'type'] = 'p'
+
+    pandapipes.pipeflow(net, check_connectivity=True, mode='all', use_numba=use_numba)
+
+    with pytest.raises(PipeflowNotConverged):
+        pandapipes.pipeflow(net, check_connectivity=False, mode='all', use_numba=use_numba)
+
+
 @pytest.mark.parametrize("use_numba", [True, False])
 def test_exclude_unconnected_junction(use_numba):
     """
@@ -390,6 +439,12 @@ def get_oos_branch(net, tbl, oosj=()):
     return get_oos(net, tbl) | net[tbl].from_junction.isin(oosj) | net[tbl].to_junction.isin(oosj)
 
 
+def get_col_slice_null(tbl):
+    if tbl == "junction":
+        return "p_bar"
+    return slice(None)
+
+
 all_tbls_funcs = {"junction": get_oos, "pipe": get_oos_branch, "sink": get_oos_node_elem,
                   "source": get_oos_node_elem, "ext_grid": get_oos_node_elem,
                   "press_control": get_oos_branch}
@@ -440,7 +495,8 @@ def test_mixed_indexing_oos2(create_mixed_indexing_grid, use_numba):
     pandapipes.pipeflow(net, mode="hydraulics", use_numba=use_numba, check_connectivity=True)
     assert all(np.all(net["res_" + tbl].loc[~oos_func(net, tbl, oos_juncs)].notnull())
                for tbl, oos_func in all_tbls_funcs.items())
-    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs)].isnull())
+    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs),
+                                            get_col_slice_null(tbl)].isnull())
                for tbl, oos_func in all_tbls_funcs.items())
     assert check_mass_flows(net)
 
@@ -457,7 +513,8 @@ def test_mixed_indexing_oos3(create_mixed_indexing_grid, use_numba):
     pandapipes.pipeflow(net, mode="hydraulics", use_numba=use_numba, check_connectivity=True)
     assert all(np.all(net["res_" + tbl].loc[~oos_func(net, tbl, oos_juncs)].notnull())
                for tbl, oos_func in all_tbls_funcs.items())
-    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs)].isnull())
+    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs),
+                                            get_col_slice_null(tbl)].isnull())
                for tbl, oos_func in all_tbls_funcs.items())
     assert check_mass_flows(net)
 
@@ -474,7 +531,8 @@ def test_mixed_indexing_oos4(create_mixed_indexing_grid, use_numba):
     pandapipes.pipeflow(net, mode="hydraulics", use_numba=use_numba, check_connectivity=True)
     assert all(np.all(net["res_" + tbl].loc[~oos_func(net, tbl, oos_juncs)].notnull())
                for tbl, oos_func in all_tbls_funcs.items())
-    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs)].isnull())
+    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs),
+                                            get_col_slice_null(tbl)].isnull())
                for tbl, oos_func in all_tbls_funcs.items())
     assert check_mass_flows(net)
 
@@ -491,7 +549,8 @@ def test_mixed_indexing_oos5(create_mixed_indexing_grid, use_numba):
     pandapipes.pipeflow(net, mode="hydraulics", use_numba=use_numba, check_connectivity=True)
     assert all(np.all(net["res_" + tbl].loc[~oos_func(net, tbl, oos_juncs)].notnull())
                for tbl, oos_func in all_tbls_funcs.items())
-    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs)].isnull())
+    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs),
+                                            get_col_slice_null(tbl)].isnull())
                for tbl, oos_func in all_tbls_funcs.items())
     assert check_mass_flows(net)
 
@@ -508,7 +567,8 @@ def test_mixed_indexing_oos6(create_mixed_indexing_grid, use_numba):
     pandapipes.pipeflow(net, mode="hydraulics", use_numba=use_numba, check_connectivity=True)
     assert all(np.all(net["res_" + tbl].loc[~oos_func(net, tbl, oos_juncs)].notnull())
                for tbl, oos_func in all_tbls_funcs.items())
-    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs)].isnull())
+    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs),
+                                            get_col_slice_null(tbl)].isnull())
                for tbl, oos_func in all_tbls_funcs.items())
     assert check_mass_flows(net)
 
@@ -516,7 +576,8 @@ def test_mixed_indexing_oos6(create_mixed_indexing_grid, use_numba):
     pandapipes.pipeflow(net, mode="hydraulics", use_numba=use_numba, check_connectivity=True)
     assert all(np.all(net["res_" + tbl].loc[~oos_func(net, tbl, oos_juncs)].notnull())
                for tbl, oos_func in all_tbls_funcs.items())
-    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs)].isnull())
+    assert all(np.all(net["res_" + tbl].loc[oos_func(net, tbl, oos_juncs),
+                                            get_col_slice_null(tbl)].isnull())
                for tbl, oos_func in all_tbls_funcs.items())
     assert check_mass_flows(net)