Repository: networkx/networkx

Issue:
weakly_connected_components: Editing yielded sets during iteration has buggy side effects
When modifying the yielded component sets during iteration from `weakly_connected_components`, strange things happen:

```python
import networkx as nx
g = nx.MultiDiGraph()

g.add_edge(0, 1)                                                                      
g.add_edge(1, 0)                                                                      
g.add_edge(0, 3)                                                                      
g.add_edge(0, 3)                                                                      

print('Okay:')
print(list(nx.weakly_connected_components(g)))
# Everything alright, gives: [{0, 1, 3}]

print('Still okay:')
for n in nx.weakly_connected_components(g):
    print(n)
    n.remove(0)
# Will print once: {0, 1, 3}

print('NOT okay:')
# Deleting the wrong node during the iteration:
for n in nx.weakly_connected_components(g):
    print(n)
    n.remove(1)
# Will print the set twice:
# {0, 1, 3}
# {0, 1, 3}
```

In case this is actually desired (which I doubt), that should be at least documented that the connected component set is not to be edited.

Discussion:

Patch:
diff --git a/.appveyor.yml b/.appveyor.yml
index f4fd737f0..a0ef5629b 100644
--- a/.appveyor.yml
+++ b/.appveyor.yml
@@ -22,7 +22,10 @@ environment:
     - PYTHON: "C:\\Python38-x64"
       PYTHON_VERSION: "3.8"
       PYTHON_ARCH: "64"
-      PIP_FLAGS: "--pre"
+
+    - PYTHON: "C:\\Python39-x64"
+      PYTHON_VERSION: "3.9.0"
+      PYTHON_ARCH: "64"
 
 install:
   # Install Python (from the official .msi of http://python.org) and pip when
diff --git a/.travis.yml b/.travis.yml
index 83c18cbd7..09ef85474 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -44,7 +44,7 @@ matrix:
           - libgdal-dev
           - graphviz
     - python: 3.8
-    - python: 3.9-dev
+    - python: 3.9
     - python: pypy3
 
 before_install:
diff --git a/INSTALL.rst b/INSTALL.rst
index d51286758..3fbee695b 100644
--- a/INSTALL.rst
+++ b/INSTALL.rst
@@ -1,7 +1,7 @@
 Install
 =======
 
-NetworkX requires Python 3.6, 3.7, or 3.8.  If you do not already
+NetworkX requires Python 3.6, 3.7, 3.8, or 3.9.  If you do not already
 have a Python environment configured on your computer, please see the
 instructions for installing the full `scientific Python stack
 <https://scipy.org/install.html>`_.
diff --git a/doc/release/index.rst b/doc/release/index.rst
index 4744dfc85..48f76df2a 100644
--- a/doc/release/index.rst
+++ b/doc/release/index.rst
@@ -18,6 +18,7 @@ period.
    :maxdepth: 2
 
    release_dev
+   release_2.5
    release_2.4
    release_2.3
    release_2.2
diff --git a/examples/algorithms/plot_dedensification.py b/examples/algorithms/plot_dedensification.py
index 97aa04a0a..309e35273 100644
--- a/examples/algorithms/plot_dedensification.py
+++ b/examples/algorithms/plot_dedensification.py
@@ -1,7 +1,7 @@
 """
-=============
+===============
 Dedensification
-=============
+===============
 
 Examples of dedensification of a graph.  Dedensification retains the structural
 pattern of the original graph and will only add compressor nodes when doing so
diff --git a/examples/drawing/plot_degree.py b/examples/drawing/plot_degree.py
index 7670927f4..6fc5994bf 100644
--- a/examples/drawing/plot_degree.py
+++ b/examples/drawing/plot_degree.py
@@ -9,9 +9,9 @@ nodes with two common techniques: a *degree-rank plot* and a
 
 In this example, a random Graph is generated with 100 nodes. The degree of
 each node is determined, and a figure is generated showing three things:
-  1. The subgraph of connected components
-  2. The degree-rank plot for the Graph, and
-  3. The degree histogram
+1. The subgraph of connected components
+2. The degree-rank plot for the Graph, and
+3. The degree histogram
 """
 import networkx as nx
 import numpy as np
diff --git a/examples/drawing/plot_rainbow_coloring.py b/examples/drawing/plot_rainbow_coloring.py
new file mode 100644
index 000000000..1d9e6075c
--- /dev/null
+++ b/examples/drawing/plot_rainbow_coloring.py
@@ -0,0 +1,68 @@
+"""
+================
+Rainbow Coloring
+================
+
+Generate a complete graph with 13 nodes in a circular layout with the
+edges colored by node distance. The node distance is given by the minimum
+number of nodes traversed along an arc between any two nodes on the circle.
+
+Such graphs are the subject of Ringel's conjecture, which states: any complete
+graph with ``2n + 1`` nodes can be tiled by any tree with ``n + 1`` nodes
+(i.e. copies of the tree can be placed over the complete graph such that each
+edge in the complete graph is covered exactly once). The edge coloring is
+helpful in determining how to place the tree copies.
+
+References
+----------
+https://www.quantamagazine.org/mathematicians-prove-ringels-graph-theory-conjecture-20200219/
+"""
+import matplotlib.pyplot as plt
+import networkx as nx
+
+# A rainbow color mapping using matplotlib's tableau colors
+node_dist_to_color = {
+    1: "tab:red",
+    2: "tab:orange",
+    3: "tab:olive",
+    4: "tab:green",
+    5: "tab:blue",
+    6: "tab:purple",
+}
+
+# Create a complete graph with an odd number of nodes
+nnodes = 13
+G = nx.complete_graph(nnodes)
+
+# A graph with (2n + 1) nodes requires n colors for the edges
+n = (nnodes - 1) // 2
+ndist_iter = list(range(1, n + 1))
+
+# Take advantage of circular symmetry in determining node distances
+ndist_iter += ndist_iter[::-1]
+
+
+def cycle(nlist, n):
+    return nlist[-n:] + nlist[:-n]
+
+
+# Rotate nodes around the circle and assign colors for each edge based on
+# node distance
+nodes = list(G.nodes())
+for i, nd in enumerate(ndist_iter):
+    for u, v in zip(nodes, cycle(nodes, i + 1)):
+        G[u][v]["color"] = node_dist_to_color[nd]
+
+pos = nx.circular_layout(G)
+# Create a figure with 1:1 aspect ratio to preserve the circle.
+fig, ax = plt.subplots(figsize=(8, 8))
+node_opts = {"node_size": 500, "node_color": "w", "edgecolors": "k", "linewidths": 2.0}
+nx.draw_networkx_nodes(G, pos, **node_opts)
+nx.draw_networkx_labels(G, pos, font_size=14)
+# Extract color from edge data
+edge_colors = [edgedata["color"] for _, _, edgedata in G.edges(data=True)]
+nx.draw_networkx_edges(G, pos, width=2.0, edge_color=edge_colors)
+
+ax.set_axis_off()
+fig.tight_layout()
+plt.show()
diff --git a/networkx/algorithms/components/strongly_connected.py b/networkx/algorithms/components/strongly_connected.py
index 8df56c5fb..5167eb833 100644
--- a/networkx/algorithms/components/strongly_connected.py
+++ b/networkx/algorithms/components/strongly_connected.py
@@ -166,8 +166,8 @@ def kosaraju_strongly_connected_components(G, source=None):
             continue
         c = nx.dfs_preorder_nodes(G, r)
         new = {v for v in c if v not in seen}
-        yield new
         seen.update(new)
+        yield new
 
 
 @not_implemented_for("undirected")
diff --git a/networkx/algorithms/components/weakly_connected.py b/networkx/algorithms/components/weakly_connected.py
index ae38d107b..aae275f4d 100644
--- a/networkx/algorithms/components/weakly_connected.py
+++ b/networkx/algorithms/components/weakly_connected.py
@@ -60,8 +60,8 @@ def weakly_connected_components(G):
     for v in G:
         if v not in seen:
             c = set(_plain_bfs(G, v))
-            yield c
             seen.update(c)
+            yield c
 
 
 @not_implemented_for("undirected")
@@ -162,7 +162,7 @@ def _plain_bfs(G, source):
         nextlevel = set()
         for v in thislevel:
             if v not in seen:
-                yield v
                 seen.add(v)
                 nextlevel.update(Gsucc[v])
                 nextlevel.update(Gpred[v])
+                yield v
diff --git a/networkx/algorithms/minors.py b/networkx/algorithms/minors.py
index d8f269f6e..a6223f0b2 100644
--- a/networkx/algorithms/minors.py
+++ b/networkx/algorithms/minors.py
@@ -62,12 +62,12 @@ def equivalence_classes(iterable, relation):
     `{2, 5, 8}`: `0`, `3`, `6`, `9` are all divisible by `3` and leave zero
     remainder; `1`, `4`, `7` leave remainder `1`; while `2`, `5` and `8` leave
     remainder `2`. We can see this by calling `equivalence_classes` with
-    with `X` and a function implementation of `R`.::
+    with `X` and a function implementation of `R`.
 
-        >>> X = set(range(10))
-        >>> def mod3(x, y): return (x - y) % 3 == 0
-        >>> equivalence_classes(X, mod3)    # doctest: +SKIP
-        {frozenset({1, 4, 7}), frozenset({8, 2, 5}), frozenset({0, 9, 3, 6})} # doctest: +SKIP
+    >>> X = set(range(10))
+    >>> def mod3(x, y): return (x - y) % 3 == 0
+    >>> equivalence_classes(X, mod3)    # doctest: +SKIP
+    {frozenset({1, 4, 7}), frozenset({8, 2, 5}), frozenset({0, 9, 3, 6})} # doctest: +SKIP
 
     """
     # For simplicity of implementation, we initialize the return value as a
@@ -192,83 +192,82 @@ def quotient_graph(
     --------
     The quotient graph of the complete bipartite graph under the "same
     neighbors" equivalence relation is `K_2`. Under this relation, two nodes
-    are equivalent if they are not adjacent but have the same neighbor set::
-
-        >>> G = nx.complete_bipartite_graph(2, 3)
-        >>> same_neighbors = lambda u, v: (
-        ...     u not in G[v] and v not in G[u] and G[u] == G[v]
-        ... )
-        >>> Q = nx.quotient_graph(G, same_neighbors)
-        >>> K2 = nx.complete_graph(2)
-        >>> nx.is_isomorphic(Q, K2)
-        True
+    are equivalent if they are not adjacent but have the same neighbor set.
+
+    >>> G = nx.complete_bipartite_graph(2, 3)
+    >>> same_neighbors = lambda u, v: (
+    ...     u not in G[v] and v not in G[u] and G[u] == G[v]
+    ... )
+    >>> Q = nx.quotient_graph(G, same_neighbors)
+    >>> K2 = nx.complete_graph(2)
+    >>> nx.is_isomorphic(Q, K2)
+    True
 
     The quotient graph of a directed graph under the "same strongly connected
     component" equivalence relation is the condensation of the graph (see
     :func:`condensation`). This example comes from the Wikipedia article
-    *`Strongly connected component`_*::
-
-        >>> G = nx.DiGraph()
-        >>> edges = [
-        ...     "ab",
-        ...     "be",
-        ...     "bf",
-        ...     "bc",
-        ...     "cg",
-        ...     "cd",
-        ...     "dc",
-        ...     "dh",
-        ...     "ea",
-        ...     "ef",
-        ...     "fg",
-        ...     "gf",
-        ...     "hd",
-        ...     "hf",
-        ... ]
-        >>> G.add_edges_from(tuple(x) for x in edges)
-        >>> components = list(nx.strongly_connected_components(G))
-        >>> sorted(sorted(component) for component in components)
-        [['a', 'b', 'e'], ['c', 'd', 'h'], ['f', 'g']]
-        >>>
-        >>> C = nx.condensation(G, components)
-        >>> component_of = C.graph["mapping"]
-        >>> same_component = lambda u, v: component_of[u] == component_of[v]
-        >>> Q = nx.quotient_graph(G, same_component)
-        >>> nx.is_isomorphic(C, Q)
-        True
+    *`Strongly connected component`_*.
+
+    >>> G = nx.DiGraph()
+    >>> edges = [
+    ...     "ab",
+    ...     "be",
+    ...     "bf",
+    ...     "bc",
+    ...     "cg",
+    ...     "cd",
+    ...     "dc",
+    ...     "dh",
+    ...     "ea",
+    ...     "ef",
+    ...     "fg",
+    ...     "gf",
+    ...     "hd",
+    ...     "hf",
+    ... ]
+    >>> G.add_edges_from(tuple(x) for x in edges)
+    >>> components = list(nx.strongly_connected_components(G))
+    >>> sorted(sorted(component) for component in components)
+    [['a', 'b', 'e'], ['c', 'd', 'h'], ['f', 'g']]
+    >>>
+    >>> C = nx.condensation(G, components)
+    >>> component_of = C.graph["mapping"]
+    >>> same_component = lambda u, v: component_of[u] == component_of[v]
+    >>> Q = nx.quotient_graph(G, same_component)
+    >>> nx.is_isomorphic(C, Q)
+    True
 
     Node identification can be represented as the quotient of a graph under the
     equivalence relation that places the two nodes in one block and each other
-    node in its own singleton block::
-
-        >>> K24 = nx.complete_bipartite_graph(2, 4)
-        >>> K34 = nx.complete_bipartite_graph(3, 4)
-        >>> C = nx.contracted_nodes(K34, 1, 2)
-        >>> nodes = {1, 2}
-        >>> is_contracted = lambda u, v: u in nodes and v in nodes
-        >>> Q = nx.quotient_graph(K34, is_contracted)
-        >>> nx.is_isomorphic(Q, C)
-        True
-        >>> nx.is_isomorphic(Q, K24)
-        True
+    node in its own singleton block.
+
+    >>> K24 = nx.complete_bipartite_graph(2, 4)
+    >>> K34 = nx.complete_bipartite_graph(3, 4)
+    >>> C = nx.contracted_nodes(K34, 1, 2)
+    >>> nodes = {1, 2}
+    >>> is_contracted = lambda u, v: u in nodes and v in nodes
+    >>> Q = nx.quotient_graph(K34, is_contracted)
+    >>> nx.is_isomorphic(Q, C)
+    True
+    >>> nx.is_isomorphic(Q, K24)
+    True
 
     The blockmodeling technique described in [1]_ can be implemented as a
-    quotient graph::
+    quotient graph.
 
-        >>> G = nx.path_graph(6)
-        >>> partition = [{0, 1}, {2, 3}, {4, 5}]
-        >>> M = nx.quotient_graph(G, partition, relabel=True)
-        >>> list(M.edges())
-        [(0, 1), (1, 2)]
+    >>> G = nx.path_graph(6)
+    >>> partition = [{0, 1}, {2, 3}, {4, 5}]
+    >>> M = nx.quotient_graph(G, partition, relabel=True)
+    >>> list(M.edges())
+    [(0, 1), (1, 2)]
 
     Here is the sample example but using partition as a dict of block sets.
-    ::
 
-        >>> G = nx.path_graph(6)
-        >>> partition = {0: {0, 1}, 2: {2, 3}, 4: {4, 5}}
-        >>> M = nx.quotient_graph(G, partition, relabel=True)
-        >>> list(M.edges())
-        [(0, 1), (1, 2)]
+    >>> G = nx.path_graph(6)
+    >>> partition = {0: {0, 1}, 2: {2, 3}, 4: {4, 5}}
+    >>> M = nx.quotient_graph(G, partition, relabel=True)
+    >>> list(M.edges())
+    [(0, 1), (1, 2)]
 
     Partitions can be represented in various ways:
     ::
@@ -413,31 +412,32 @@ def contracted_nodes(G, u, v, self_loops=True, copy=True):
     Parameters
     ----------
     G : NetworkX graph
-       The graph whose nodes will be contracted.
+        The graph whose nodes will be contracted.
 
     u, v : nodes
-       Must be nodes in `G`.
+        Must be nodes in `G`.
 
     self_loops : Boolean
-       If this is True, any edges joining `u` and `v` in `G` become
-       self-loops on the new node in the returned graph.
+        If this is True, any edges joining `u` and `v` in `G` become
+        self-loops on the new node in the returned graph.
 
     copy : Boolean
         If this is True (default True), make a copy of
         `G` and return that instead of directly changing `G`.
 
+
     Returns
     -------
     Networkx graph
-       If Copy is True:
-       A new graph object of the same type as `G` (leaving `G` unmodified)
-       with `u` and `v` identified in a single node. The right node `v`
-       will be merged into the node `u`, so only `u` will appear in the
-       returned graph.
-       if Copy is False:
-       Modifies `G` with `u` and `v` identified in a single node.
-       The right node `v` will be merged into the node `u`, so
-       only `u` will appear in the returned graph.
+        If Copy is True,
+        A new graph object of the same type as `G` (leaving `G` unmodified)
+        with `u` and `v` identified in a single node. The right node `v`
+        will be merged into the node `u`, so only `u` will appear in the
+        returned graph.
+        If copy is False,
+        Modifies `G` with `u` and `v` identified in a single node.
+        The right node `v` will be merged into the node `u`, so
+        only `u` will appear in the returned graph.
 
     Notes
     -----
@@ -450,7 +450,7 @@ def contracted_nodes(G, u, v, self_loops=True, copy=True):
     Examples
     --------
     Contracting two nonadjacent nodes of the cycle graph on four nodes `C_4`
-    yields the path graph (ignoring parallel edges)::
+    yields the path graph (ignoring parallel edges):
 
     >>> G = nx.cycle_graph(4)
     >>> M = nx.contracted_nodes(G, 1, 3)
@@ -557,22 +557,22 @@ def contracted_edge(G, edge, self_loops=True):
 
     Examples
     --------
-    Attempting to contract two nonadjacent nodes yields an error::
+    Attempting to contract two nonadjacent nodes yields an error:
 
-        >>> G = nx.cycle_graph(4)
-        >>> nx.contracted_edge(G, (1, 3))
-        Traceback (most recent call last):
-          ...
-        ValueError: Edge (1, 3) does not exist in graph G; cannot contract it
+    >>> G = nx.cycle_graph(4)
+    >>> nx.contracted_edge(G, (1, 3))
+    Traceback (most recent call last):
+      ...
+    ValueError: Edge (1, 3) does not exist in graph G; cannot contract it
 
     Contracting two adjacent nodes in the cycle graph on *n* nodes yields the
-    cycle graph on *n - 1* nodes::
+    cycle graph on *n - 1* nodes:
 
-        >>> C5 = nx.cycle_graph(5)
-        >>> C4 = nx.cycle_graph(4)
-        >>> M = nx.contracted_edge(C5, (0, 1), self_loops=False)
-        >>> nx.is_isomorphic(M, C4)
-        True
+    >>> C5 = nx.cycle_graph(5)
+    >>> C4 = nx.cycle_graph(4)
+    >>> M = nx.contracted_edge(C5, (0, 1), self_loops=False)
+    >>> nx.is_isomorphic(M, C4)
+    True
 
     See also
     --------
diff --git a/networkx/release.py b/networkx/release.py
index 17dc97714..30f20ea2d 100644
--- a/networkx/release.py
+++ b/networkx/release.py
@@ -219,6 +219,7 @@ classifiers = [
     "Programming Language :: Python :: 3.6",
     "Programming Language :: Python :: 3.7",
     "Programming Language :: Python :: 3.8",
+    "Programming Language :: Python :: 3.9",
     "Programming Language :: Python :: 3 :: Only",
     "Topic :: Software Development :: Libraries :: Python Modules",
     "Topic :: Scientific/Engineering :: Bio-Informatics",
diff --git a/requirements/default.txt b/requirements/default.txt
index 0bb350521..1a8aed1ad 100644
--- a/requirements/default.txt
+++ b/requirements/default.txt
@@ -1,5 +1,5 @@
 decorator>=4.4
-numpy>=1.19; platform_python_implementation!='PyPy' and python_version<'3.9'
-scipy>=1.5; platform_python_implementation!='PyPy' and python_version<'3.9'
-matplotlib>=3.3; platform_python_implementation!='PyPy' and python_version<'3.9'
-pandas>=1.1; platform_python_implementation!='PyPy' and python_version<'3.9'
+numpy>=1.19; platform_python_implementation!='PyPy' and python_version<'3.10'
+scipy>=1.5; platform_python_implementation!='PyPy' and python_version<'3.10'
+matplotlib>=3.3; platform_python_implementation!='PyPy' and python_version<'3.10'
+pandas>=1.1; platform_python_implementation!='PyPy' and python_version<'3.10'
diff --git a/requirements/doc.txt b/requirements/doc.txt
index 54a77b1f2..3fc3093f0 100644
--- a/requirements/doc.txt
+++ b/requirements/doc.txt
@@ -1,6 +1,6 @@
-sphinx==3.2.1
+sphinx==3.3.1
 sphinx_rtd_theme==0.5.0
-sphinx-gallery==0.7.0
-pillow>=7.2.0
+sphinx-gallery==0.8.1
+pillow>=8.0
 nb2plots>=0.6
 texext>=0.6.6


Test patch:
diff --git a/.github/workflows/testing.yml b/.github/workflows/testing.yml
index 8828c4839..2cc00567f 100644
--- a/.github/workflows/testing.yml
+++ b/.github/workflows/testing.yml
@@ -12,7 +12,7 @@ jobs:
     strategy:
       max-parallel: 3
       matrix:
-        python-version: [3.6, 3.7, 3.8]
+        python-version: [3.6, 3.7, 3.8, 3.9]
     steps:
     - uses: actions/checkout@v2
     - name: Set up Python ${{ matrix.python-version }}
diff --git a/networkx/algorithms/components/tests/test_strongly_connected.py b/networkx/algorithms/components/tests/test_strongly_connected.py
index 959b333c8..54d74cdc2 100644
--- a/networkx/algorithms/components/tests/test_strongly_connected.py
+++ b/networkx/algorithms/components/tests/test_strongly_connected.py
@@ -188,6 +188,22 @@ class TestStronglyConnected:
         )
         pytest.raises(NetworkXNotImplemented, nx.condensation, G)
 
+    strong_cc_methods = (
+        nx.strongly_connected_components,
+        nx.kosaraju_strongly_connected_components,
+        nx.strongly_connected_components_recursive,
+    )
+
+    @pytest.mark.parametrize("get_components", strong_cc_methods)
+    def test_connected_mutability(self, get_components):
+        DG = nx.path_graph(5, create_using=nx.DiGraph)
+        G = nx.disjoint_union(DG, DG)
+        seen = set()
+        for component in get_components(G):
+            assert len(seen & component) == 0
+            seen.update(component)
+            component.clear()
+
 
 #    Commented out due to variability on Travis-CI hardware/operating systems
 #    def test_linear_time(self):
diff --git a/networkx/algorithms/components/tests/test_weakly_connected.py b/networkx/algorithms/components/tests/test_weakly_connected.py
index 393f688d5..9964236e5 100644
--- a/networkx/algorithms/components/tests/test_weakly_connected.py
+++ b/networkx/algorithms/components/tests/test_weakly_connected.py
@@ -76,3 +76,12 @@ class TestWeaklyConnected:
         pytest.raises(NetworkXNotImplemented, nx.weakly_connected_components, G)
         pytest.raises(NetworkXNotImplemented, nx.number_weakly_connected_components, G)
         pytest.raises(NetworkXNotImplemented, nx.is_weakly_connected, G)
+
+    def test_connected_mutability(self):
+        DG = nx.path_graph(5, create_using=nx.DiGraph)
+        G = nx.disjoint_union(DG, DG)
+        seen = set()
+        for component in nx.weakly_connected_components(G):
+            assert len(seen & component) == 0
+            seen.update(component)
+            component.clear()
diff --git a/networkx/classes/tests/test_multigraph.py b/networkx/classes/tests/test_multigraph.py
index 5fee1897c..ec07c229b 100644
--- a/networkx/classes/tests/test_multigraph.py
+++ b/networkx/classes/tests/test_multigraph.py
@@ -373,3 +373,14 @@ class TestEdgeSubgraph:
 
         """
         assert self.G.graph is self.H.graph
+
+
+def test_multigraph_add_edges_from_four_tuple_misordered():
+    """add_edges_from expects 4-tuples of the format (u, v, key, data_dict).
+
+    Ensure 4-tuples of form (u, v, data_dict, key) raises exception.
+    """
+    G = nx.MultiGraph()
+    with pytest.raises(TypeError):
+        # key/data values flipped in 4-tuple
+        G.add_edges_from([(0, 1, {"color": "red"}, 0)])