Speed up from_* and add a quantities filter

.gitignore

@@ -155,4 +155,6 @@ docs/api/*.qmd
 
 uv.lock
 
-tests/testdata/confidential/*
+tests/testdata/confidential/*
+
+tests/profiling/output/

docs/user-guide/network.qmd

@@ -250,12 +250,13 @@ A MIKE 1D network contains multiple levels that are unified into a generic netwo
 
 Large result files can contain thousands of nodes and gridpoints. Loading all of that data into memory is slow and may cause memory issues — especially when you only need the timeseries at a handful of nodes where observations exist.
 
-Both constructors accept the same two optional arguments to restrict what gets loaded:
+Both constructors accept the same three optional arguments to restrict what gets loaded:
 
 | Argument | Type | Effect |
 |---|---|---|
 | `nodes` | `None` \| `str` \| `list[str]` | Control which nodes have timeseries data loaded. `None` (default) loads all nodes; `[]` skips all node data; a name or list loads only those nodes. |
 | `reaches` | `None` \| `str` \| `list[str]` | Control which reaches have intermediate gridpoint data populated. `None` (default) loads everything; `[]` skips all gridpoints; a name or list of names loads only those reaches. |
+| `quantities` | `None` \| `str` \| `list[str]` | Control which quantities are read at each selected location. `None` (default) reads everything; `[]` reads nothing; a name or list reads only those quantities. |
 
 ::: {.callout-note}
 Selective loading only controls **which timeseries are held in memory**. The full network topology (nodes, reaches, lengths) is always constructed so that `find()`, `recall()`, and graph algorithms still work on the complete network.
@@ -283,6 +284,13 @@ network_subset = Network.from_mike(
 network_subset
 ```
 
+`quantities` cuts the read down a third way, which helps when a calibration loop rebuilds the network for every trial but only scores one quantity. A location that does not carry a requested quantity simply stays topology-only:
+
+```{python}
+network_discharge = Network.from_mike(path_to_res1d, quantities="Discharge")
+network_discharge.quantities
+```
+
 When only some nodes are loaded, `to_dataframe()` and `to_dataset()` only contain columns for those nodes — the rest are graph-connected but data-free:
 
 ```{python}

pyproject.toml

@@ -46,7 +46,14 @@ classifiers = [
 networks = ["mikeio1d", "networkx"]
 
 [dependency-groups]
-dev = ["pytest", "plotly >= 4.5", "ruff==0.6.2", "netCDF4", "dask"]
+dev = [
+    "pytest",
+    "plotly >= 4.5",
+    "ruff==0.6.2",
+    "netCDF4",
+    "dask",
+    "snakeviz>=2.2.2",
+]
 
 docs = ["quartodoc==0.11.1", "nbformat", "nbconvert", "ipykernel", "griffe<2"]
 

src/modelskill/model/adapters/_res1d.py

@@ -10,7 +10,15 @@
 from modelskill.network import NetworkNode, ReachBreakPoint, NetworkReach
 
 
-def _simplify_colnames(node: ResultNode | ResultGridPoint) -> pd.DataFrame:
+# Topology-only nodes and gridpoints all share this frame instead of each
+# allocating its own. A large network has two per reach, which profiling showed
+# to be the biggest single cost of a filtered load. Never mutate it in place.
+_EMPTY_DATA = pd.DataFrame()
+
+
+def _simplify_colnames(
+    node: ResultNode | ResultGridPoint, quantities: set[str] | None = None
+) -> pd.DataFrame:
     # We remove suffixes and indexes so the columns contain only the quantity names
 
     # Some formats keep no timeseries at all on some locations - MIKE 11, for instance,
@@ -22,9 +30,24 @@ def _simplify_colnames(node: ResultNode | ResultGridPoint) -> pd.DataFrame:
     # The columns in a Res1D dataframe follow the convention "Quantity:Location:Sublocation"
     # where Location refers to the node id or the reach id followed by the chainage.
     RES1D_NAME_SEP = ":"
-    df = node.to_dataframe()
+
+    available = list(node.quantities)
+    wanted = (
+        available if quantities is None else [q for q in available if q in quantities]
+    )
+
+    if not wanted:
+        # A location need not carry every requested quantity; it stays topology-only.
+        return _EMPTY_DATA
+
+    if len(wanted) == len(available):
+        # Reading the whole location is one interop call rather than one per quantity.
+        df = node.to_dataframe()
+    else:
+        df = pd.concat([getattr(node, q).to_dataframe() for q in wanted], axis=1)
+
     renamer_dict = {}
-    for quantity in node.quantities:
+    for quantity in wanted:
         column_pairs = [
             (col, quantity)
             for col in df.columns
@@ -86,7 +109,7 @@ def __init__(
         boundary: dict[str, pd.DataFrame] | None = None,
     ):
         self._id = id
-        self._data = pd.DataFrame() if data is None else data
+        self._data = _EMPTY_DATA if data is None else data
         self._boundary = {} if boundary is None else boundary
 
     @property
@@ -107,7 +130,7 @@ def __init__(
         self, reach_id: str, chainage: float, data: pd.DataFrame | None = None
     ):
         self._id = (reach_id, chainage)
-        self._data = pd.DataFrame() if data is None else data
+        self._data = _EMPTY_DATA if data is None else data
 
     @property
     def id(self) -> tuple[str, float]:
@@ -129,6 +152,7 @@ def __init__(
         *,
         populate_gridpoints: bool = True,
         length: float | None = None,
+        quantities: set[str] | None = None,
     ):
         self._id = reach.name
 
@@ -163,7 +187,9 @@ def __init__(
             GridPoint(
                 gridpoint.reach_name,
                 gridpoint.chainage,
-                _simplify_colnames(gridpoint) if populate_gridpoints else None,
+                _simplify_colnames(gridpoint, quantities)
+                if populate_gridpoints
+                else None,
             )
             for gridpoint in intermediate_gridpoints
         ]

src/modelskill/network.py

@@ -422,6 +422,7 @@ def from_mike(
         *,
         nodes: str | list[str] | None = None,
         reaches: str | list[str] | None = None,
+        quantities: str | list[str] | None = None,
     ) -> Network:
         """Create a Network from a MIKE 1D or MIKE 11 result file.
 
@@ -448,6 +449,17 @@ def from_mike(
             * A single reach name or a list of reach names — only those reaches
               get gridpoint data; others are topology-only.
             * ``[]`` (empty list) — no gridpoint data is loaded at all.
+        quantities : str, list of str, or None, optional
+            Controls which quantities are read at each selected location.
+
+            * ``None`` *(default)* — every quantity is read.
+            * A single quantity name or a list of names — only those are read.
+            * ``[]`` (empty list) — no data is read at all.
+
+            A location that does not carry a requested quantity becomes
+            topology-only rather than an error, so this composes with ``nodes``
+            and ``reaches`` on files where nodes and reaches hold different
+            quantities.
 
         Returns
         -------
@@ -484,6 +496,11 @@ def from_mike(
         ...     reaches=["reach_1"],
         ... )
 
+        Read a single quantity, for a calibration loop that only scores
+        discharge:
+
+        >>> network = Network.from_mike("model.res1d", quantities="Discharge")
+
         Notes
         -----
         MIKE 11 keeps its timeseries on reach gridpoints rather than on nodes,
@@ -498,6 +515,7 @@ def from_mike(
             res,
             nodes=nodes,
             reaches=reaches,
+            quantities=quantities,
             allowed=_MIKE_EXTENSIONS,
             caller="from_mike",
         )
@@ -511,6 +529,7 @@ def from_epanet(
         inp: str | Path | None = None,
         nodes: str | list[str] | None = None,
         reaches: str | list[str] | None = None,
+        quantities: str | list[str] | None = None,
     ) -> Network:
         """Create a Network from an EPANET result file and its companions.
 
@@ -539,6 +558,9 @@ def from_epanet(
             Which reaches get their gridpoint data loaded. See
             :meth:`from_mike`. EPANET results have no intermediate gridpoints,
             so this argument has no effect.
+        quantities : str, list of str, or None, optional
+            Which quantities are read at each selected location. See
+            :meth:`from_mike`.
 
         Returns
         -------
@@ -596,6 +618,7 @@ def from_epanet(
             res,
             nodes=nodes,
             reaches=reaches,
+            quantities=quantities,
             allowed=_EPANET_EXTENSIONS,
             caller="from_epanet",
             resx=resx,
@@ -611,6 +634,7 @@ def _from_mikeio1d(
         reaches: str | list[str] | None,
         allowed: frozenset[str],
         caller: str,
+        quantities: str | list[str] | None = None,
         resx: str | Path | Res1D | None = None,
         inp: str | Path | None = None,
     ) -> Network:
@@ -664,8 +688,22 @@ def _from_mikeio1d(
         extra = None if resx is None else cls._open_companion_result(res, resx)
         lengths = None if inp is None else cls._read_companion_lengths(inp)
 
+        # None is threaded through as "read everything" rather than expanded to
+        # res.quantities, which would only cost a lookup for the same result.
+        if quantities is None:
+            quantities_set: set[str] | None = None
+        elif isinstance(quantities, str):
+            quantities_set = {quantities}
+        else:
+            quantities_set = set(quantities)
+
         list_of_reaches = cls._load_res1d_network(
-            res, nodes_list, reaches_list, extra=extra, lengths=lengths
+            res,
+            nodes_list,
+            reaches_list,
+            extra=extra,
+            lengths=lengths,
+            quantities=quantities_set,
         )
         return cls(list_of_reaches)
 
@@ -796,6 +834,7 @@ def _load_res1d_network(
         *,
         extra: Res1D | None = None,
         lengths: dict[str, float] | None = None,
+        quantities: set[str] | None = None,
     ) -> list[Res1DReach]:
         from modelskill.model.adapters._res1d import (
             Res1DReach,
@@ -812,21 +851,31 @@ def _load_res1d_network(
         # potential memory issues. For this reason, we create this intermediate step that populates
         # only the data in the passed nodes
 
+        # A node shared by several reaches is visited once per reach endpoint, but
+        # _generate_graph keeps only the first copy of its data, so read it once.
+        # The boundary is per-reach and must stay outside the cache.
+        node_data: dict[str, pd.DataFrame] = {}
+
         def _init_node(reach: ResultReach, is_end: bool) -> Res1DNode:
             id = reach.end_node if is_end else reach.start_node
             gpt_idx = -1 if is_end else 0
             if id in nodes_set:
-                node = res.nodes[id]
-                df = _simplify_colnames(node)
-                # Merged here rather than up front so selective loading still
-                # decides what is held in memory.
-                if extra is not None and id in extra.nodes:
-                    df = _merge_extra_quantities(
-                        df, _simplify_colnames(extra.nodes[id]), node_id=id
-                    )
+                if id not in node_data:
+                    df = _simplify_colnames(res.nodes[id], quantities)
+                    # Merged here rather than up front so selective loading still
+                    # decides what is held in memory.
+                    if extra is not None and id in extra.nodes:
+                        df = _merge_extra_quantities(
+                            df,
+                            _simplify_colnames(extra.nodes[id], quantities),
+                            node_id=id,
+                        )
+                    node_data[id] = df
                 overlapping_gridpoint = reach.gridpoints[gpt_idx]
-                boundary = _simplify_colnames(overlapping_gridpoint)
-                return Res1DNode(id, data=df, boundary={reach.name: boundary})
+                boundary = _simplify_colnames(overlapping_gridpoint, quantities)
+                return Res1DNode(
+                    id, data=node_data[id], boundary={reach.name: boundary}
+                )
             else:
                 return Res1DNode(id)
 
@@ -837,6 +886,7 @@ def _init_node(reach: ResultReach, is_end: bool) -> Res1DNode:
                 _init_node(reach, True),
                 populate_gridpoints=reach.name in reaches_set,
                 length=lengths.get(reach.name),
+                quantities=quantities,
             )
             for reach in res.reaches.values()
         ]