docs: def block_cm(cm):
      """Return whether ``cm`` can be arranged as a block connectivity matrix.
  
      If so, the corresponding mechanism/purview is trivially reducible.
      Technically, only square matrices are "block diagonal", but the notion of
      connectivity carries over.
  
      We test for block connectivity by trying to grow a block of nodes such
      that:
  
      - 'source' nodes only input to nodes in the block
      - 'sink' nodes only receive inputs from source nodes in the block
  
      For example, the following connectivity matrix represents connections from
      ``nodes1 = A, B, C`` to ``nodes2 = D, E, F, G`` (without loss of
      generality, note that ``nodes1`` and ``nodes2`` may share elements)::
  
           D  E  F  G
        A [1, 1, 0, 0]
        B [1, 1, 0, 0]
        C [0, 0, 1, 1]
  
      Since nodes |AB| only connect to nodes |DE|, and node |C| only connects to
      nodes |FG|, the subgraph is reducible, because the cut ::
  
        A,B    C
        ─── ✕ ───
        D,E   F,G
  
      does not change the structure of the graph.
      """
      if np.any(cm.sum(1) == 0):
          return True
      if np.all(cm.sum(1) == 1):
          return True
  
      outputs = list(range(cm.shape[1]))
  
      # CM helpers:
      def outputs_of(nodes):
          """Return all nodes that `nodes` connect to (output to)."""
          return np.where(cm[nodes, :].sum(0))[0]
  
      def inputs_to(nodes):
          """Return all nodes which connect to (input to) `nodes`."""
          return np.where(cm[:, nodes].sum(1))[0]
  
      # Start: source node with most outputs
      sources = [np.argmax(cm.sum(1))]
      sinks = outputs_of(sources)
      sink_inputs = inputs_to(sinks)
  
      while True:
          if np.array_equal(sink_inputs, sources):
              # sources exclusively connect to sinks.
              # There are no other nodes which connect sink nodes,
              # hence set(sources) + set(sinks) form a component
              # which is not connected to the rest of the graph
              return True
  
          # Recompute sources, sinks, and sink_inputs
          sources = sink_inputs
          sinks = outputs_of(sources)
          sink_inputs = inputs_to(sinks)
  
          # Considering all output nodes?
          if np.array_equal(sinks, outputs):
              return False