code: def tutte_polynomial(G):
      r"""Returns the Tutte polynomial of `G`
  
      This function computes the Tutte polynomial via an iterative version of
      the deletion-contraction algorithm.
  
      The Tutte polynomial `T_G(x, y)` is a fundamental graph polynomial invariant in
      two variables. It encodes a wide array of information related to the
      edge-connectivity of a graph; "Many problems about graphs can be reduced to
      problems of finding and evaluating the Tutte polynomial at certain values" [1]_.
      In fact, every deletion-contraction-expressible feature of a graph is a
      specialization of the Tutte polynomial [2]_ (see Notes for examples).
  
      There are several equivalent definitions; here are three:
  
      Def 1 (rank-nullity expansion): For `G` an undirected graph, `n(G)` the
      number of vertices of `G`, `E` the edge set of `G`, `V` the vertex set of
      `G`, and `c(A)` the number of connected components of the graph with vertex
      set `V` and edge set `A` [3]_:
  
      .. math::
  
          T_G(x, y) = \sum_{A \in E} (x-1)^{c(A) - c(E)} (y-1)^{c(A) + |A| - n(G)}
  
      Def 2 (spanning tree expansion): Let `G` be an undirected graph, `T` a spanning
      tree of `G`, and `E` the edge set of `G`. Let `E` have an arbitrary strict
      linear order `L`. Let `B_e` be the unique minimal nonempty edge cut of
      $E \setminus T \cup {e}$. An edge `e` is internally active with respect to
      `T` and `L` if `e` is the least edge in `B_e` according to the linear order
      `L`. The internal activity of `T` (denoted `i(T)`) is the number of edges
      in $E \setminus T$ that are internally active with respect to `T` and `L`.
      Let `P_e` be the unique path in $T \cup {e}$ whose source and target vertex
      are the same. An edge `e` is externally active with respect to `T` and `L`
      if `e` is the least edge in `P_e` according to the linear order `L`. The
      external activity of `T` (denoted `e(T)`) is the number of edges in
      $E \setminus T$ that are externally active with respect to `T` and `L`.
      Then [4]_ [5]_:
  
      .. math::
  
          T_G(x, y) = \sum_{T \text{ a spanning tree of } G} x^{i(T)} y^{e(T)}
  
      Def 3 (deletion-contraction recurrence): For `G` an undirected graph, `G-e`
      the graph obtained from `G` by deleting edge `e`, `G/e` the graph obtained
      from `G` by contracting edge `e`, `k(G)` the number of cut-edges of `G`,
      and `l(G)` the number of self-loops of `G`:
  
      .. math::
          T_G(x, y) = \begin{cases}
      	   x^{k(G)} y^{l(G)}, & \text{if all edges are cut-edges or self-loops} \\
             T_{G-e}(x, y) + T_{G/e}(x, y), & \text{otherwise, for an arbitrary edge $e$ not a cut-edge or loop}
          \end{cases}
  
      Parameters
      ----------
      G : NetworkX graph
  
      Returns
      -------
      instance of `sympy.core.add.Add`
          A Sympy expression representing the Tutte polynomial for `G`.
  
      Examples
      --------
      >>> C = nx.cycle_graph(5)
      >>> nx.tutte_polynomial(C)
      x**4 + x**3 + x**2 + x + y
  
      >>> D = nx.diamond_graph()
      >>> nx.tutte_polynomial(D)
      x**3 + 2*x**2 + 2*x*y + x + y**2 + y
  
      Notes
      -----
      Some specializations of the Tutte polynomial:
  
      - `T_G(1, 1)` counts the number of spanning trees of `G`
      - `T_G(1, 2)` counts the number of connected spanning subgraphs of `G`
      - `T_G(2, 1)` counts the number of spanning forests in `G`
      - `T_G(0, 2)` counts the number of strong orientations of `G`
      - `T_G(2, 0)` counts the number of acyclic orientations of `G`
  
      Edge contraction is defined and deletion-contraction is introduced in [6]_.
      Combinatorial meaning of the coefficients is introduced in [7]_.
      Universality, properties, and applications are discussed in [8]_.
  
      Practically, up-front computation of the Tutte polynomial may be useful when
      users wish to repeatedly calculate edge-connectivity-related information
      about one or more graphs.
  
      References
      ----------
      .. [1] M. Brandt,
         "The Tutte Polynomial."
         Talking About Combinatorial Objects Seminar, 2015
         https://math.berkeley.edu/~brandtm/talks/tutte.pdf
      .. [2] A. Björklund, T. Husfeldt, P. Kaski, M. Koivisto,
         "Computing the Tutte polynomial in vertex-exponential time"
         49th Annual IEEE Symposium on Foundations of Computer Science, 2008
         https://ieeexplore.ieee.org/abstract/document/4691000
      .. [3] Y. Shi, M. Dehmer, X. Li, I. Gutman,
         "Graph Polynomials," p. 14
      .. [4] Y. Shi, M. Dehmer, X. Li, I. Gutman,
         "Graph Polynomials," p. 46
      .. [5] A. Nešetril, J. Goodall,
         "Graph invariants, homomorphisms, and the Tutte polynomial"
         https://iuuk.mff.cuni.cz/~andrew/Tutte.pdf
      .. [6] D. B. West,
         "Introduction to Graph Theory," p. 84
      .. [7] G. Coutinho,
         "A brief introduction to the Tutte polynomial"
         Structural Analysis of Complex Networks, 2011
         https://homepages.dcc.ufmg.br/~gabriel/seminars/coutinho_tuttepolynomial_seminar.pdf
      .. [8] J. A. Ellis-Monaghan, C. Merino,
         "Graph polynomials and their applications I: The Tutte polynomial"
         Structural Analysis of Complex Networks, 2011
         https://arxiv.org/pdf/0803.3079.pdf
      """
      import sympy
  
      x = sympy.Symbol("x")
      y = sympy.Symbol("y")
      stack = deque()
      stack.append(nx.MultiGraph(G))
  
      polynomial = 0
      while stack:
          G = stack.pop()
          bridges = set(nx.bridges(G))
  
          e = None
          for i in G.edges:
              if (i[0], i[1]) not in bridges and i[0] != i[1]:
                  e = i
                  break
          if not e:
              loops = list(nx.selfloop_edges(G, keys=True))
              polynomial += x ** len(bridges) * y ** len(loops)
          else:
              # deletion-contraction
              C = nx.contracted_edge(G, e, self_loops=True)
              C.remove_edge(e[0], e[0])
              G.remove_edge(*e)
              stack.append(G)
              stack.append(C)
      return sympy.simplify(polynomial)
  
imports: from collections import deque
  import networkx
  from networkx.utils import not_implemented_for
label_desc: Use image node_id 1 for calling a global function with example usage: tutte_polynomial(G) and returns: sympy