Fast compiled graphs

This is a Cython implementation of the base class for sparse and dense graphs in Sage. It is not intended for use on its own. Specific graph types should extend this base class and implement missing functionalities. Whenever possible, specific methods should also be overridden with implementations that suit the graph type under consideration.

Data structure

The class CGraph maintains the following variables:

  • cdef int num_verts
  • cdef int num_arcs
  • cdef int *in_degrees
  • cdef int *out_degrees
  • cdef bitset_t active_vertices

The bitset active_vertices is a list of all available vertices for use, but only the ones which are set are considered to actually be in the graph. The variables num_verts and num_arcs are self-explanatory. Note that num_verts is the number of bits set in active_vertices, not the full length of the bitset. The arrays in_degrees and out_degrees are of the same length as the bitset.

For more information about active vertices, see the documentation for the method realloc.

class sage.graphs.base.c_graph.CGraph

Bases: object

Compiled sparse and dense graphs.

add_arc(u, v)

Add the given arc to this graph.

INPUT:

  • u – integer; the tail of an arc.
  • v – integer; the head of an arc.

OUTPUT:

  • Raise NotImplementedError. This method is not implemented at the CGraph level. A child class should provide a suitable implementation.

See also

  • add_arcadd_arc method for sparse graphs.
  • add_arcadd_arc method for dense graphs.

EXAMPLE:

sage: from sage.graphs.base.c_graph import CGraph
sage: G = CGraph()
sage: G.add_arc(0, 1)
Traceback (most recent call last):
...
NotImplementedError
add_vertex(k=-1)

Adds vertex k to the graph.

INPUT:

  • k – nonnegative integer or -1 (default: -1). If k = -1, a new vertex is added and the integer used is returned. That is, for k = -1, this function will find the first available vertex that is not in self and add that vertex to this graph.

OUTPUT:

  • -1 – indicates that no vertex was added because the current allocation is already full or the vertex is out of range.
  • nonnegative integer – this vertex is now guaranteed to be in the graph.

See also

  • add_vertex_unsafe() – add a vertex to a graph. This method is potentially unsafe. You should instead use add_vertex().
  • add_vertices() – add a bunch of vertices to a graph.

EXAMPLES:

Adding vertices to a sparse graph:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: G = SparseGraph(3, extra_vertices=3)
sage: G.add_vertex(3)
3
sage: G.add_arc(2, 5)
Traceback (most recent call last):
...
RuntimeError: Vertex (5) is not a vertex of the graph.
sage: G.add_arc(1, 3)
sage: G.has_arc(1, 3)
True
sage: G.has_arc(2, 3)
False

Adding vertices to a dense graph:

sage: from sage.graphs.base.dense_graph import DenseGraph
sage: G = DenseGraph(3, extra_vertices=3)
sage: G.add_vertex(3)
3
sage: G.add_arc(2,5)
Traceback (most recent call last):
...
RuntimeError: Vertex (5) is not a vertex of the graph.
sage: G.add_arc(1, 3)
sage: G.has_arc(1, 3)
True
sage: G.has_arc(2, 3)
False

Repeatedly adding a vertex using k = -1 will allocate more memory as required:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: G = SparseGraph(3, extra_vertices=0)
sage: G.verts()
[0, 1, 2]
sage: for i in range(10):
...       _ = G.add_vertex(-1);
...
sage: G.verts()
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12]
sage: from sage.graphs.base.dense_graph import DenseGraph
sage: G = DenseGraph(3, extra_vertices=0)
sage: G.verts()
[0, 1, 2]
sage: for i in range(12):
...       _ = G.add_vertex(-1);
...
sage: G.verts()
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14]

TESTS:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: G = SparseGraph(3, extra_vertices=0)
sage: G.add_vertex(6)
Traceback (most recent call last):
...
RuntimeError: Requested vertex is past twice the allocated range: use realloc.
sage: from sage.graphs.base.dense_graph import DenseGraph
sage: G = DenseGraph(3, extra_vertices=0)
sage: G.add_vertex(6)
Traceback (most recent call last):
...
RuntimeError: Requested vertex is past twice the allocated range: use realloc.
add_vertices(verts)

Adds vertices from the iterable verts.

INPUT:

  • verts – an iterable of vertices.

OUTPUT:

See also

EXAMPLE:

Adding vertices for sparse graphs:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: S = SparseGraph(nverts=4, extra_vertices=4)
sage: S.verts()
[0, 1, 2, 3]
sage: S.add_vertices([3,5,7,9])
sage: S.verts()
[0, 1, 2, 3, 5, 7, 9]
sage: S.realloc(20)
sage: S.verts()
[0, 1, 2, 3, 5, 7, 9]

Adding vertices for dense graphs:

sage: from sage.graphs.base.dense_graph import DenseGraph
sage: D = DenseGraph(nverts=4, extra_vertices=4)
sage: D.verts()
[0, 1, 2, 3]
sage: D.add_vertices([3,5,7,9])
sage: D.verts()
[0, 1, 2, 3, 5, 7, 9]
sage: D.realloc(20)
sage: D.verts()
[0, 1, 2, 3, 5, 7, 9]
all_arcs(u, v)

Return the labels of all arcs from u to v.

INPUT:

  • u – integer; the tail of an arc.
  • v – integer; the head of an arc.

OUTPUT:

  • Raise NotImplementedError. This method is not implemented at the CGraph level. A child class should provide a suitable implementation.

See also

  • all_arcsall_arcs method for sparse graphs.

EXAMPLE:

sage: from sage.graphs.base.c_graph import CGraph
sage: G = CGraph()
sage: G.all_arcs(0, 1)
Traceback (most recent call last):
...
NotImplementedError
check_vertex(n)

Checks that n is a vertex of self.

This method is different from has_vertex(). The current method raises an error if n is not a vertex of this graph. On the other hand, has_vertex() returns a boolean to signify whether or not n is a vertex of this graph.

INPUT:

  • n – a nonnegative integer representing a vertex.

OUTPUT:

  • Raise an error if n is not a vertex of this graph.

See also

  • has_vertex() – determine whether this graph has a specific vertex.

EXAMPLES:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: S = SparseGraph(nverts=10, expected_degree=3, extra_vertices=10)
sage: S.check_vertex(4)
sage: S.check_vertex(12)
Traceback (most recent call last):
...
RuntimeError: Vertex (12) is not a vertex of the graph.
sage: S.check_vertex(24)
Traceback (most recent call last):
...
RuntimeError: Vertex (24) is not a vertex of the graph.
sage: S.check_vertex(-19)
Traceback (most recent call last):
...
RuntimeError: Vertex (-19) is not a vertex of the graph.
sage: from sage.graphs.base.dense_graph import DenseGraph
sage: D = DenseGraph(nverts=10, extra_vertices=10)
sage: D.check_vertex(4)
sage: D.check_vertex(12)
Traceback (most recent call last):
...
RuntimeError: Vertex (12) is not a vertex of the graph.
sage: D.check_vertex(24)
Traceback (most recent call last):
...
RuntimeError: Vertex (24) is not a vertex of the graph.
sage: D.check_vertex(-19)
Traceback (most recent call last):
...
RuntimeError: Vertex (-19) is not a vertex of the graph.
current_allocation()

Report the number of vertices allocated.

INPUT:

  • None.

OUTPUT:

  • The number of vertices allocated. This number is usually different from the order of a graph. We may have allocated enough memory for a graph to hold n > 0 vertices, but the order (actual number of vertices) of the graph could be less than n.

EXAMPLES:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: S = SparseGraph(nverts=4, extra_vertices=4)
sage: S.current_allocation()
8
sage: S.add_vertex(6)
6
sage: S.current_allocation()
8
sage: S.add_vertex(10)
10
sage: S.current_allocation()
16
sage: S.add_vertex(40)
Traceback (most recent call last):
...
RuntimeError: Requested vertex is past twice the allocated range: use realloc.
sage: S.realloc(50)
sage: S.add_vertex(40)
40
sage: S.current_allocation()
50
sage: S.realloc(30)
-1
sage: S.current_allocation()
50
sage: S.del_vertex(40)
sage: S.realloc(30)
sage: S.current_allocation()
30

The actual number of vertices in a graph might be less than the number of vertices allocated for the graph:

sage: from sage.graphs.base.dense_graph import DenseGraph
sage: G = DenseGraph(nverts=3, extra_vertices=2)
sage: order = len(G.verts())
sage: order
3
sage: G.current_allocation()
5
sage: order < G.current_allocation()
True
del_all_arcs(u, v)

Delete all arcs from u to v.

INPUT:

  • u – integer; the tail of an arc.
  • v – integer; the head of an arc.

OUTPUT:

  • Raise NotImplementedError. This method is not implemented at the CGraph level. A child class should provide a suitable implementation.

See also

EXAMPLE:

sage: from sage.graphs.base.c_graph import CGraph
sage: G = CGraph()
sage: G.del_all_arcs(0,1)
Traceback (most recent call last):
...
NotImplementedError
del_vertex(v)

Deletes the vertex v, along with all edges incident to it. If v is not in self, fails silently.

INPUT:

  • v – a nonnegative integer representing a vertex.

OUTPUT:

  • None.

See also

  • del_vertex_unsafe() – delete a vertex from a graph. This method is potentially unsafe. Use del_vertex() instead.

EXAMPLES:

Deleting vertices of sparse graphs:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: G = SparseGraph(3)
sage: G.add_arc(0, 1)
sage: G.add_arc(0, 2)
sage: G.add_arc(1, 2)
sage: G.add_arc(2, 0)
sage: G.del_vertex(2)
sage: for i in range(2):
...       for j in range(2):
...           if G.has_arc(i, j):
...               print i, j
0 1
sage: G = SparseGraph(3)
sage: G.add_arc(0, 1)
sage: G.add_arc(0, 2)
sage: G.add_arc(1, 2)
sage: G.add_arc(2, 0)
sage: G.del_vertex(1)
sage: for i in xrange(3):
...       for j in xrange(3):
...           if G.has_arc(i, j):
...               print i, j
0 2
2 0

Deleting vertices of dense graphs:

sage: from sage.graphs.base.dense_graph import DenseGraph
sage: G = DenseGraph(4)
sage: G.add_arc(0, 1); G.add_arc(0, 2)
sage: G.add_arc(3, 1); G.add_arc(3, 2)
sage: G.add_arc(1, 2)
sage: G.verts()
[0, 1, 2, 3]
sage: G.del_vertex(3); G.verts()
[0, 1, 2]
sage: for i in range(3):
...       for j in range(3):
...           if G.has_arc(i, j):
...               print i, j
...
0 1
0 2
1 2

If the vertex to be deleted is not in this graph, then fail silently:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: G = SparseGraph(3)
sage: G.verts()
[0, 1, 2]
sage: G.has_vertex(3)
False
sage: G.del_vertex(3)
sage: G.verts()
[0, 1, 2]
sage: from sage.graphs.base.dense_graph import DenseGraph
sage: G = DenseGraph(5)
sage: G.verts()
[0, 1, 2, 3, 4]
sage: G.has_vertex(6)
False
sage: G.del_vertex(6)
sage: G.verts()
[0, 1, 2, 3, 4]
has_arc(u, v)

Determine whether or not the given arc is in this graph.

INPUT:

  • u – integer; the tail of an arc.
  • v – integer; the head of an arc.

OUTPUT:

  • Print a Not Implemented! message. This method is not implemented at the CGraph level. A child class should provide a suitable implementation.

See also

  • has_archas_arc method for sparse graphs.
  • has_archas_arc method for dense graphs.

EXAMPLE:

sage: from sage.graphs.base.c_graph import CGraph
sage: G = CGraph()
sage: G.has_arc(0, 1)
Not Implemented!
False
has_vertex(n)

Determine whether the vertex n is in self.

This method is different from check_vertex(). The current method returns a boolean to signify whether or not n is a vertex of this graph. On the other hand, check_vertex() raises an error if n is not a vertex of this graph.

INPUT:

  • n – a nonnegative integer representing a vertex.

OUTPUT:

  • True if n is a vertex of this graph; False otherwise.

See also

  • check_vertex() – raise an error if this graph does not contain a specific vertex.

EXAMPLES:

Upon initialization, a SparseGraph or DenseGraph has the first nverts vertices:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: S = SparseGraph(nverts=10, expected_degree=3, extra_vertices=10)
sage: S.has_vertex(6)
True
sage: S.has_vertex(12)
False
sage: S.has_vertex(24)
False
sage: S.has_vertex(-19)
False
sage: from sage.graphs.base.dense_graph import DenseGraph
sage: D = DenseGraph(nverts=10, extra_vertices=10)
sage: D.has_vertex(6)
True
sage: D.has_vertex(12)
False
sage: D.has_vertex(24)
False
sage: D.has_vertex(-19)
False
in_neighbors(v)

Gives the in-neighbors of the vertex v.

INPUT:

  • v – integer representing a vertex of this graph.

OUTPUT:

  • Raise NotImplementedError. This method is not implemented at the CGraph level. A child class should provide a suitable implementation.

See also

EXAMPLE:

sage: from sage.graphs.base.c_graph import CGraph
sage: G = CGraph()
sage: G.in_neighbors(0)
Traceback (most recent call last):
...
NotImplementedError
out_neighbors(u)

Gives the out-neighbors of the vertex u.

INPUT:

  • u – integer representing a vertex of this graph.

OUTPUT:

  • Raise NotImplementedError. This method is not implemented at the CGraph level. A child class should provide a suitable implementation.

See also

  • out_neighborsout_neighbors implementation for sparse graphs.
  • out_neighborsout_neighbors implementation for dense graphs.

EXAMPLE:

sage: from sage.graphs.base.c_graph import CGraph
sage: G = CGraph()
sage: G.out_neighbors(0)
Traceback (most recent call last):
...
NotImplementedError
realloc(total)

Reallocate the number of vertices to use, without actually adding any.

INPUT:

  • total – integer; the total size to make the array of vertices.

OUTPUT:

  • Raise a NotImplementedError. This method is not implemented in this base class. A child class should provide a suitable implementation.

See also

  • realloc – a realloc implementation for sparse graphs.
  • realloc – a realloc implementation for dense graphs.

EXAMPLES:

First, note that realloc() is implemented for SparseGraph and DenseGraph differently, and is not implemented at the CGraph level:

sage: from sage.graphs.base.c_graph import CGraph
sage: G = CGraph()
sage: G.realloc(20)
Traceback (most recent call last):
...
NotImplementedError

The realloc implementation for sparse graphs:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: S = SparseGraph(nverts=4, extra_vertices=4)
sage: S.current_allocation()
8
sage: S.add_vertex(6)
6
sage: S.current_allocation()
8
sage: S.add_vertex(10)
10
sage: S.current_allocation()
16
sage: S.add_vertex(40)
Traceback (most recent call last):
...
RuntimeError: Requested vertex is past twice the allocated range: use realloc.
sage: S.realloc(50)
sage: S.add_vertex(40)
40
sage: S.current_allocation()
50
sage: S.realloc(30)
-1
sage: S.current_allocation()
50
sage: S.del_vertex(40)
sage: S.realloc(30)
sage: S.current_allocation()
30

The realloc implementation for dense graphs:

sage: from sage.graphs.base.dense_graph import DenseGraph
sage: D = DenseGraph(nverts=4, extra_vertices=4)
sage: D.current_allocation()
8
sage: D.add_vertex(6)
6
sage: D.current_allocation()
8
sage: D.add_vertex(10)
10
sage: D.current_allocation()
16
sage: D.add_vertex(40)
Traceback (most recent call last):
...
RuntimeError: Requested vertex is past twice the allocated range: use realloc.
sage: D.realloc(50)
sage: D.add_vertex(40)
40
sage: D.current_allocation()
50
sage: D.realloc(30)
-1
sage: D.current_allocation()
50
sage: D.del_vertex(40)
sage: D.realloc(30)
sage: D.current_allocation()
30
verts()

Returns a list of the vertices in self.

INPUT:

  • None.

OUTPUT:

  • A list of all vertices in this graph.

EXAMPLE:

sage: from sage.graphs.base.sparse_graph import SparseGraph
sage: S = SparseGraph(nverts=4, extra_vertices=4)
sage: S.verts()
[0, 1, 2, 3]
sage: S.add_vertices([3,5,7,9])
sage: S.verts()
[0, 1, 2, 3, 5, 7, 9]
sage: S.realloc(20)
sage: S.verts()
[0, 1, 2, 3, 5, 7, 9]
sage: from sage.graphs.base.dense_graph import DenseGraph
sage: G = DenseGraph(3, extra_vertices=2)
sage: G.verts()
[0, 1, 2]
sage: G.del_vertex(0)
sage: G.verts()
[1, 2]
class sage.graphs.base.c_graph.CGraphBackend

Bases: sage.graphs.base.graph_backends.GenericGraphBackend

Base class for sparse and dense graph backends.

sage: from sage.graphs.base.c_graph import CGraphBackend

This class is extended by SparseGraphBackend and DenseGraphBackend, which are fully functional backends. This class is mainly just for vertex functions, which are the same for both. A CGraphBackend will not work on its own:

sage: from sage.graphs.base.c_graph import CGraphBackend
sage: CGB = CGraphBackend()
sage: CGB.degree(0, True)
Traceback (most recent call last):
...
AttributeError: 'CGraphBackend' object has no attribute 'vertex_ints'

The appropriate way to use these backends is via Sage graphs:

sage: G = Graph(30, implementation="c_graph")
sage: G.add_edges([(0,1), (0,3), (4,5), (9, 23)])
sage: G.edges(labels=False)
[(0, 1), (0, 3), (4, 5), (9, 23)]

See also

add_vertex(name)

Add a vertex to self.

INPUT:

  • name – the vertex to be added (must be hashable).

OUTPUT:

  • None.

See also

  • add_vertices() – add a bunch of vertices of this graph.
  • has_vertex() – returns whether or not this graph has a specific vertex.

EXAMPLE:

sage: D = sage.graphs.base.dense_graph.DenseGraphBackend(9)
sage: D.add_vertex(10)
sage: D.add_vertex([])
Traceback (most recent call last):
...
TypeError: unhashable type: 'list'
sage: S = sage.graphs.base.sparse_graph.SparseGraphBackend(9)
sage: S.add_vertex(10)
sage: S.add_vertex([])
Traceback (most recent call last):
...
TypeError: unhashable type: 'list'
add_vertices(vertices)

Add vertices to self.

INPUT:

  • vertices – iterator of vertex labels.

OUTPUT:

  • None.

See also

EXAMPLE:

sage: D = sage.graphs.base.sparse_graph.SparseGraphBackend(1)
sage: D.add_vertices([1,2,3])
sage: G = sage.graphs.base.sparse_graph.SparseGraphBackend(0)
sage: G.add_vertices([0,1])
sage: list(G.iterator_verts(None))
[0, 1]
sage: list(G.iterator_edges([0,1], True))
[]
sage: import sage.graphs.base.dense_graph
sage: D = sage.graphs.base.dense_graph.DenseGraphBackend(9)
sage: D.add_vertices([10,11,12])
bidirectional_dijkstra(x, y)

Returns the shortest path between x and y using a bidirectional version of Dijkstra’s algorithm.

INPUT:

  • x – the starting vertex in the shortest path from x to y.
  • y – the end vertex in the shortest path from x to y.

OUTPUT:

  • A list of vertices in the shortest path from x to y.

EXAMPLE:

sage: G = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: for (u,v) in G.edges(labels=None):
...      G.set_edge_label(u,v,1)
sage: G.shortest_path(0, 1, by_weight=True)
[0, 1]

TEST:

Bugfix from #7673

sage: G = Graph(implementation="networkx")
sage: G.add_edges([(0,1,9),(0,2,8),(1,2,7)])
sage: Gc = G.copy(implementation='c_graph')
sage: sp = G.shortest_path_length(0,1,by_weight=True)
sage: spc = Gc.shortest_path_length(0,1,by_weight=True)
sage: sp == spc
True

Returns a breadth-first search from vertex v.

INPUT:

  • v – a vertex from which to start the breadth-first search.
  • reverse – boolean (default: False). This is only relevant to digraphs. If this is a digraph, consider the reversed graph in which the out-neighbors become the in-neighbors and vice versa.
  • ignore_direction – boolean (default: False). This is only relevant to digraphs. If this is a digraph, ignore all orientations and consider the graph as undirected.

ALGORITHM:

Below is a general template for breadth-first search.

  • Input: A directed or undirected graph G = (V, E) of order n > 0. A vertex s from which to start the search. The vertices are numbered from 1 to n = |V|, i.e. V = \{1, 2, \dots, n\}.
  • Output: A list D of distances of all vertices from s. A tree T rooted at s.
  1. Q \leftarrow [s] # a queue of nodes to visit
  2. D \leftarrow [\infty, \infty, \dots, \infty] # n copies of \infty
  3. D[s] \leftarrow 0
  4. T \leftarrow [\,]
  5. while \text{length}(Q) > 0 do
    1. v \leftarrow \text{dequeue}(Q)
    2. for each w \in \text{adj}(v) do # for digraphs, use out-neighbor set \text{oadj}(v)
      1. if D[w] = \infty then
        1. D[w] \leftarrow D[v] + 1
        2. \text{enqueue}(Q, w)
        3. \text{append}(T, vw)
  6. return (D, T)

See also

EXAMPLES:

Breadth-first search of the Petersen graph starting at vertex 0:

sage: G = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: list(G.breadth_first_search(0))
[0, 1, 4, 5, 2, 6, 3, 9, 7, 8]

Visiting German cities using breadth-first search:

sage: G = Graph({"Mannheim": ["Frankfurt","Karlsruhe"],
...   "Frankfurt": ["Mannheim","Wurzburg","Kassel"],
...   "Kassel": ["Frankfurt","Munchen"],
...   "Munchen": ["Kassel","Nurnberg","Augsburg"],
...   "Augsburg": ["Munchen","Karlsruhe"],
...   "Karlsruhe": ["Mannheim","Augsburg"],
...   "Wurzburg": ["Frankfurt","Erfurt","Nurnberg"],
...   "Nurnberg": ["Wurzburg","Stuttgart","Munchen"],
...   "Stuttgart": ["Nurnberg"],
...   "Erfurt": ["Wurzburg"]}, implementation="c_graph")
sage: list(G.breadth_first_search("Frankfurt"))
['Frankfurt', 'Mannheim', 'Kassel', 'Wurzburg', 'Karlsruhe', 'Munchen', 'Erfurt', 'Nurnberg', 'Augsburg', 'Stuttgart']
degree(v, directed)

Return the degree of the vertex v.

INPUT:

  • v – a vertex of the graph.
  • directed – boolean; whether to take into account the orientation of this graph in counting the degree of v.

OUTPUT:

  • The degree of vertex v.

EXAMPLE:

sage: from sage.graphs.base.sparse_graph import SparseGraphBackend
sage: B = SparseGraphBackend(7)
sage: B.degree(3, False)
0

TESTS:

Ensure that ticket #8395 is fixed.

sage: def my_add_edges(G, m, n):
...       for i in range(m):
...           u = randint(0, n)
...           v = randint(0, n)
...           G.add_edge(u, v)
sage: G = Graph({1:[1]}); G
Looped graph on 1 vertex
sage: G.edges(labels=False)
[(1, 1)]
sage: G.degree(); G.size()
[2]
1
sage: sum(G.degree()) == 2 * G.size()
True
sage: G = Graph({1:[1,2], 2:[2,3]}, loops=True); G
Looped graph on 3 vertices
sage: my_add_edges(G, 100, 50)
sage: sum(G.degree()) == 2 * G.size()
True
sage: G = Graph({1:[2,2], 2:[3]}); G
Multi-graph on 3 vertices
sage: G.edges(labels=False)
[(1, 2), (1, 2), (2, 3)]
sage: G.degree(); G.size()
[2, 3, 1]
3
sage: sum(G.degree()) == 2 * G.size()
True
sage: G.allow_loops(True); G
Looped multi-graph on 3 vertices
sage: my_add_edges(G, 100, 50)
sage: sum(G.degree()) == 2 * G.size()
True
sage: D = DiGraph({1:[2], 2:[1,3]}); D
Digraph on 3 vertices
sage: D.edges(labels=False)
[(1, 2), (2, 1), (2, 3)]
sage: D.degree(); D.size()
[2, 3, 1]
3
sage: sum(D.degree()) == 2 * D.size()
True
sage: D.allow_loops(True); D
Looped digraph on 3 vertices
sage: my_add_edges(D, 100, 50)
sage: sum(D.degree()) == 2 * D.size()
True
sage: D.allow_multiple_edges(True)
sage: my_add_edges(D, 200, 50)
sage: sum(D.degree()) == 2 * D.size()
True
sage: G = Graph({1:[2,2,2]})
sage: G.allow_loops(True)
sage: G.add_edge(1,1)
sage: G.add_edge(1,1)
sage: G.edges(labels=False)
[(1, 1), (1, 1), (1, 2), (1, 2), (1, 2)]
sage: G.degree(1)
7
sage: G.allow_loops(False)
sage: G.edges(labels=False)
[(1, 2), (1, 2), (1, 2)]
sage: G.degree(1)
3
sage: G = Graph({1:{2:['a','a','a']}})
sage: G.allow_loops(True)
sage: G.add_edge(1,1,'b')
sage: G.add_edge(1,1,'b')
sage: G.add_edge(1,1)
sage: G.add_edge(1,1)
sage: G.edges()
[(1, 1, None), (1, 1, None), (1, 1, 'b'), (1, 1, 'b'), (1, 2, 'a'), (1, 2, 'a'), (1, 2, 'a')]
sage: G.degree(1)
11
sage: G.allow_loops(False)
sage: G.edges()
[(1, 2, 'a'), (1, 2, 'a'), (1, 2, 'a')]
sage: G.degree(1)
3
sage: G = Graph({1:{2:['a','a','a']}})
sage: G.allow_loops(True)
sage: G.add_edge(1,1,'b')
sage: G.add_edge(1,1,'b')
sage: G.edges()
[(1, 1, 'b'), (1, 1, 'b'), (1, 2, 'a'), (1, 2, 'a'), (1, 2, 'a')]
sage: G.degree(1)
7
sage: G.allow_loops(False)
sage: G.edges()
[(1, 2, 'a'), (1, 2, 'a'), (1, 2, 'a')]
sage: G.degree(1)
3
del_vertex(v)

Delete a vertex in self, failing silently if the vertex is not in the graph.

INPUT:

  • v – vertex to be deleted.

OUTPUT:

  • None.

See also

EXAMPLE:

sage: D = sage.graphs.base.dense_graph.DenseGraphBackend(9)
sage: D.del_vertex(0)
sage: D.has_vertex(0)
False
sage: S = sage.graphs.base.sparse_graph.SparseGraphBackend(9)
sage: S.del_vertex(0)
sage: S.has_vertex(0)
False
del_vertices(vertices)

Delete vertices from an iterable container.

INPUT:

  • vertices – iterator of vertex labels.

OUTPUT:

See also

EXAMPLE:

sage: import sage.graphs.base.dense_graph
sage: D = sage.graphs.base.dense_graph.DenseGraphBackend(9)
sage: D.del_vertices([7,8])
sage: D.has_vertex(7)
False
sage: D.has_vertex(6)
True
sage: D = sage.graphs.base.sparse_graph.SparseGraphBackend(9)
sage: D.del_vertices([1,2,3])
sage: D.has_vertex(1)
False
sage: D.has_vertex(0)
True

Returns a depth-first search from vertex v.

INPUT:

  • v – a vertex from which to start the depth-first search.
  • reverse – boolean (default: False). This is only relevant to digraphs. If this is a digraph, consider the reversed graph in which the out-neighbors become the in-neighbors and vice versa.
  • ignore_direction – boolean (default: False). This is only relevant to digraphs. If this is a digraph, ignore all orientations and consider the graph as undirected.

ALGORITHM:

Below is a general template for depth-first search.

  • Input: A directed or undirected graph G = (V, E) of order n > 0. A vertex s from which to start the search. The vertices are numbered from 1 to n = |V|, i.e. V = \{1, 2, \dots, n\}.
  • Output: A list D of distances of all vertices from s. A tree T rooted at s.
  1. S \leftarrow [s] # a stack of nodes to visit
  2. D \leftarrow [\infty, \infty, \dots, \infty] # n copies of \infty
  3. D[s] \leftarrow 0
  4. T \leftarrow [\,]
  5. while \text{length}(S) > 0 do
    1. v \leftarrow \text{pop}(S)
    2. for each w \in \text{adj}(v) do # for digraphs, use out-neighbor set \text{oadj}(v)
      1. if D[w] = \infty then
        1. D[w] \leftarrow D[v] + 1
        2. \text{push}(S, w)
        3. \text{append}(T, vw)
  6. return (D, T)

See also

EXAMPLES:

Traversing the Petersen graph using depth-first search:

sage: G = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: list(G.depth_first_search(0))
[0, 5, 8, 6, 9, 7, 2, 3, 4, 1]

Visiting German cities using depth-first search:

sage: G = Graph({"Mannheim": ["Frankfurt","Karlsruhe"],
...   "Frankfurt": ["Mannheim","Wurzburg","Kassel"],
...   "Kassel": ["Frankfurt","Munchen"],
...   "Munchen": ["Kassel","Nurnberg","Augsburg"],
...   "Augsburg": ["Munchen","Karlsruhe"],
...   "Karlsruhe": ["Mannheim","Augsburg"],
...   "Wurzburg": ["Frankfurt","Erfurt","Nurnberg"],
...   "Nurnberg": ["Wurzburg","Stuttgart","Munchen"],
...   "Stuttgart": ["Nurnberg"],
...   "Erfurt": ["Wurzburg"]}, implementation="c_graph")
sage: list(G.depth_first_search("Frankfurt"))
['Frankfurt', 'Wurzburg', 'Nurnberg', 'Munchen', 'Kassel', 'Augsburg', 'Karlsruhe', 'Mannheim', 'Stuttgart', 'Erfurt']
has_vertex(v)

Returns whether v is a vertex of self.

INPUT:

  • v – any object.

OUTPUT:

  • True if v is a vertex of this graph; False otherwise.

EXAMPLE:

sage: from sage.graphs.base.sparse_graph import SparseGraphBackend
sage: B = SparseGraphBackend(7)
sage: B.has_vertex(6)
True
sage: B.has_vertex(7)
False
is_connected()

Returns whether the graph is connected.

INPUT:

  • None.

OUTPUT:

  • True if this graph is connected; False otherwise.

EXAMPLES:

Petersen’s graph is connected:

sage: DiGraph(graphs.PetersenGraph(),implementation="c_graph").is_connected()
True

While the disjoint union of two of them is not:

sage: DiGraph(2*graphs.PetersenGraph(),implementation="c_graph").is_connected()
False
A graph with non-integer vertex labels::
sage: Graph(graphs.CubeGraph(3), implementation=’c_graph’).is_connected() True
is_directed_acyclic(certificate=False)

Returns whether the graph is both directed and acylic (possibly with a certificate)

INPUT:

  • certificate – whether to return a certificate (False by default).

OUTPUT:

When certificate=False, returns a boolean value. When certificate=True :

  • If the graph is acyclic, returns a pair (True, ordering) where ordering is a list of the vertices such that u appears before v in ordering if u, v is an edge.
  • Else, returns a pair (False, cycle) where cycle is a list of vertices representing a circuit in the graph.

Note

The graph is assumed to be directed. An exception is raised if it is not.

EXAMPLES:

At first, the following graph is acyclic:

sage: D = DiGraph({ 0:[1,2,3], 4:[2,5], 1:[8], 2:[7], 3:[7], 5:[6,7], 7:[8], 6:[9], 8:[10], 9:[10] })
sage: D.plot(layout='circular').show()
sage: D.is_directed_acyclic()
True

Adding an edge from 9 to 7 does not change it:

sage: D.add_edge(9,7)
sage: D.is_directed_acyclic()
True

We can obtain as a proof an ordering of the vertices such that u appears before v if uv is an edge of the graph:

sage: D.is_directed_acyclic(certificate = True)
(True, [4, 5, 6, 9, 0, 1, 2, 3, 7, 8, 10])

Adding an edge from 7 to 4, though, makes a difference:

sage: D.add_edge(7,4)
sage: D.is_directed_acyclic()
False

Indeed, it creates a circuit 7, 4, 5:

sage: D.is_directed_acyclic(certificate = True)
(False, [7, 4, 5])

Checking acyclic graphs are indeed acyclic

sage: def random_acyclic(n, p):
...    g = graphs.RandomGNP(n, p)
...    h = DiGraph()
...    h.add_edges([ ((u,v) if u<v else (v,u)) for u,v,_ in g.edges() ])
...    return h
...
sage: all( random_acyclic(100, .2).is_directed_acyclic()    # long time
...        for i in range(50))                              # long time
True
is_strongly_connected()

Returns whether the graph is strongly connected.

INPUT:

  • None.

OUTPUT:

  • True if this graph is strongly connected; False otherwise.

EXAMPLES:

The circuit on 3 vertices is obviously strongly connected:

sage: g = DiGraph({0: [1], 1: [2], 2: [0]}, implementation="c_graph")
sage: g.is_strongly_connected()
True

But a transitive triangle is not:

sage: g = DiGraph({0: [1,2], 1: [2]}, implementation="c_graph")
sage: g.is_strongly_connected()
False
iterator_in_nbrs(v)

Returns an iterator over the incoming neighbors of v.

INPUT:

  • v – a vertex of this graph.

OUTPUT:

  • An iterator over the in-neighbors of the vertex v.

See also

EXAMPLE:

sage: P = DiGraph(graphs.PetersenGraph().to_directed(), implementation="c_graph")
sage: list(P._backend.iterator_in_nbrs(0))
[1, 4, 5]
iterator_nbrs(v)

Returns an iterator over the neighbors of v.

INPUT:

  • v – a vertex of this graph.

OUTPUT:

  • An iterator over the neighbors the vertex v.

See also

EXAMPLE:

sage: P = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: list(P._backend.iterator_nbrs(0))
[1, 4, 5]
iterator_out_nbrs(v)

Returns an iterator over the outgoing neighbors of v.

INPUT:

  • v – a vertex of this graph.

OUTPUT:

  • An iterator over the out-neighbors of the vertex v.

See also

EXAMPLE:

sage: P = DiGraph(graphs.PetersenGraph().to_directed(), implementation="c_graph")
sage: list(P._backend.iterator_out_nbrs(0))
[1, 4, 5]
iterator_verts(verts=None)

Returns an iterator over the vertices of self intersected with verts.

INPUT:

  • verts – an iterable container of objects (default: None).

OUTPUT:

  • If verts=None, return an iterator over all vertices of this graph.
  • If verts is an iterable container of vertices, find the intersection of verts with the vertex set of this graph and return an iterator over the resulting intersection.

See also

EXAMPLE:

sage: P = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: list(P._backend.iterator_verts(P))
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9]
sage: list(P._backend.iterator_verts())
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9]
sage: list(P._backend.iterator_verts([1, 2, 3]))
[1, 2, 3]
sage: list(P._backend.iterator_verts([1, 2, 10]))
[1, 2]
loops(new=None)

Returns whether loops are allowed in this graph.

INPUT:

  • new – (default: None); boolean (to set) or None (to get).

OUTPUT:

  • If new=None, return True if this graph allows self-loops or False if self-loops are not allowed.
  • If new is a boolean, set the self-loop permission of this graph according to the boolean value of new.

EXAMPLE:

sage: G = Graph(implementation='c_graph')
sage: G._backend.loops()
False
sage: G._backend.loops(True)
sage: G._backend.loops()
True
name(new=None)

Returns the name of this graph.

INPUT:

  • new – (default: None); boolean (to set) or None (to get).

OUTPUT:

  • If new=None, return the name of this graph. Otherwise, set the name of this graph to the value of new.

EXAMPLE:

sage: G = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: G._backend.name()
'Petersen graph'
sage: G._backend.name("Peter Pan's graph")
sage: G._backend.name()
"Peter Pan's graph"
num_edges(directed)

Returns the number of edges in self.

INPUT:

  • directed – boolean; whether to count (u,v) and (v,u) as one or two edges.

OUTPUT:

  • If directed=True, counts the number of directed edges in this graph. Otherwise, return the size of this graph.

See also

EXAMPLE:

sage: G = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: G._backend.num_edges(False)
15

TESTS:

Ensure that ticket #8395 is fixed.

sage: G = Graph({1:[1]}); G
Looped graph on 1 vertex
sage: G.edges(labels=False)
[(1, 1)]
sage: G.size()
1
sage: G = Graph({1:[2,2]}); G
Multi-graph on 2 vertices
sage: G.edges(labels=False)
[(1, 2), (1, 2)]
sage: G.size()
2
sage: G = Graph({1:[1,1]}); G
Looped multi-graph on 1 vertex
sage: G.edges(labels=False)
[(1, 1), (1, 1)]
sage: G.size()
2
sage: D = DiGraph({1:[1]}); D
Looped digraph on 1 vertex
sage: D.edges(labels=False)
[(1, 1)]
sage: D.size()
1
sage: D = DiGraph({1:[2,2], 2:[1,1]}); D
Multi-digraph on 2 vertices
sage: D.edges(labels=False)
[(1, 2), (1, 2), (2, 1), (2, 1)]
sage: D.size()
4
sage: D = DiGraph({1:[1,1]}); D
Looped multi-digraph on 1 vertex
sage: D.edges(labels=False)
[(1, 1), (1, 1)]
sage: D.size()
2
sage: from sage.graphs.base.sparse_graph import SparseGraphBackend
sage: S = SparseGraphBackend(7)
sage: S.num_edges(directed=False)
0
sage: S.loops(True)
sage: S.add_edge(1, 1, None, directed=False)
sage: S.num_edges(directed=False)
1
sage: S.multiple_edges(True)
sage: S.add_edge(1, 1, None, directed=False)
sage: S.num_edges(directed=False)
2
sage: from sage.graphs.base.dense_graph import DenseGraphBackend
sage: D = DenseGraphBackend(7)
sage: D.num_edges(directed=False)
0
sage: D.loops(True)
sage: D.add_edge(1, 1, None, directed=False)
sage: D.num_edges(directed=False)
1
num_verts()

Returns the number of vertices in self.

INPUT:

  • None.

OUTPUT:

  • The order of this graph.

See also

  • num_edges() – return the number of (directed) edges in this graph.

EXAMPLE:

sage: G = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: G._backend.num_verts()
10
relabel(perm, directed)

Relabels the graph according to perm.

INPUT:

  • perm – anything which represents a permutation as v --> perm[v], for example a dict or a list.
  • directed – ignored (this is here for compatibility with other backends).

EXAMPLES:

sage: G = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: G._backend.relabel(range(9,-1,-1), False)
sage: G.edges()
[(0, 2, None),
 (0, 3, None),
 (0, 5, None),
 (1, 3, None),
 (1, 4, None),
 (1, 6, None),
 (2, 4, None),
 (2, 7, None),
 (3, 8, None),
 (4, 9, None),
 (5, 6, None),
 (5, 9, None),
 (6, 7, None),
 (7, 8, None),
 (8, 9, None)]
shortest_path(x, y)

Returns the shortest path between x and y.

INPUT:

  • x – the starting vertex in the shortest path from x to y.
  • y – the end vertex in the shortest path from x to y.

OUTPUT:

  • A list of vertices in the shortest path from x to y.

EXAMPLE:

sage: G = Graph(graphs.PetersenGraph(), implementation="c_graph")
sage: G.shortest_path(0, 1)
[0, 1]
shortest_path_all_vertices(v, cutoff=None)

Returns for each vertex u a shortest v-u path.

INPUT:

  • v – a starting vertex in the shortest path.
  • cutoff – maximal distance. Longer paths will not be returned.

OUTPUT:

  • A list which associates to each vertex u the shortest path between u and v if there is one.

Note

The weight of edges is not taken into account.

ALGORITHM:

This is just a breadth-first search.

EXAMPLES:

On the Petersen Graph:

sage: g = graphs.PetersenGraph()
sage: paths = g._backend.shortest_path_all_vertices(0)
sage: all([ len(paths[v]) == 0 or len(paths[v])-1 == g.distance(0,v) for v in g])
True

On a disconnected graph

sage: g = 2*graphs.RandomGNP(20,.3)
sage: paths = g._backend.shortest_path_all_vertices(0)
sage: all([ (not paths.has_key(v) and g.distance(0,v) == +Infinity) or len(paths[v])-1 == g.distance(0,v) for v in g])
True
strongly_connected_component_containing_vertex(v)

Returns the strongly connected component containing the given vertex.

INPUT:

  • v – a vertex

EXAMPLES:

The digraph obtained from the PetersenGraph has an unique strongly connected component:

sage: g = DiGraph(graphs.PetersenGraph())
sage: g.strongly_connected_component_containing_vertex(0)
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9]

In the Butterfly DiGraph, each vertex is a strongly connected component:

sage: g = digraphs.ButterflyGraph(3)
sage: all([[v] == g.strongly_connected_component_containing_vertex(v) for v in g])
True
class sage.graphs.base.c_graph.Search_iterator

Bases: object

An iterator for traversing a (di)graph.

This class is commonly used to perform a depth-first or breadth-first search. The class does not build all at once in memory the whole list of visited vertices. The class maintains the following variables:

  • graph – a graph whose vertices are to be iterated over.
  • direction – integer; this determines the position at which vertices to be visited are removed from the list stack. For breadth-first search (BFS), element removal occurs at the start of the list, as signified by the value direction=0. This is because in implementations of BFS, the list of vertices to visit are usually maintained by a queue, so element insertion and removal follow a first-in first-out (FIFO) protocol. For depth-first search (DFS), element removal occurs at the end of the list, as signified by the value direction=-1. The reason is that DFS is usually implemented using a stack to maintain the list of vertices to visit. Hence, element insertion and removal follow a last-in first-out (LIFO) protocol.
  • stack – a list of vertices to visit.
  • seen – a list of vertices that are already visited.
  • test_out – boolean; whether we want to consider the out-neighbors of the graph to be traversed. For undirected graphs, we consider both the in- and out-neighbors. However, for digraphs we only traverse along out-neighbors.
  • test_in – boolean; whether we want to consider the in-neighbors of the graph to be traversed. For undirected graphs, we consider both the in- and out-neighbors.

EXAMPLE:

sage: g = graphs.PetersenGraph()
sage: list(g.breadth_first_search(0))
[0, 1, 4, 5, 2, 6, 3, 9, 7, 8]
next()

x.next() -> the next value, or raise StopIteration

sage.graphs.base.c_graph.all_pairs_shortest_path_BFS(gg)

Computes the shortest path and the distances between each pair of vertices through successive breadth-first-searches

OUTPUT:

This function return a pair of dictionaries (distances,paths) where distance[u][v] denotes the distance of a shortest path from u to v and paths[u][v] denotes an inneighbor w of v such that dist(u,v)= 1 + dist(u,w).

Warning

Because this function works on matrices whose size is quadratic compared to the number of vertices, it uses short variables instead of long ones to divide by 2 the size in memory. This means that the current implementation does not run on a graph of more than 65536 nodes (this can be easily changed if necessary, but would require much more memory. It may be worth writing two versions). For information, the current version of the algorithm on a graph with 65536=2^{16} nodes creates in memory 2 tables on 2^{32} short elements (2bytes each), for a total of 2^{34} bytes or 16 gigabytes.

EXAMPLE:

On a grid:

sage: g = graphs.Grid2dGraph(10,10)
sage: from sage.graphs.base.c_graph import all_pairs_shortest_path_BFS
sage: dist, path = all_pairs_shortest_path_BFS(g)
sage: all( dist[u][v] == g.distance(u,v) for u in g for v in g )
True

TESTS:

Too large graphs:

sage: all_pairs_shortest_path_BFS(Graph(65536))
Traceback (most recent call last):
...
ValueError: The graph backend contains more than 65535 nodes
sage.graphs.base.c_graph.floyd_warshall(gg, paths=True, distances=False)

Computes the shortest path/distances between all pairs of vertices.

For more information on the Floyd-Warshall algorithm, see the Wikipedia article on Floyd-Warshall.

INPUT:

  • gg – the graph on which to work.
  • paths (boolean) – whether to return the dictionary of shortest paths. Set to True by default.
  • distances (boolean) – whether to return the dictionary of distances. Set to False by default.

OUTPUT:

Depending on the input, this function return the dictionary of paths, the dictionary of distances, or a pair of dictionaries (distances, paths) where distance[u][v] denotes the distance of a shortest path from u to v and paths[u][v] denotes an inneighbor w of v such that dist(u,v)= 1 + dist(u,w).

Warning

Because this function works on matrices whose size is quadratic compared to the number of vertices, it uses short variables instead of long ones to divide by 2 the size in memory. This means that the current implementation does not run on a graph of more than 65536 nodes (this can be easily changed if necessary, but would require much more memory. It may be worth writing two versions). For information, the current version of the algorithm on a graph with 65536=2^{16} nodes creates in memory 2 tables on 2^{32} short elements (2bytes each), for a total of 2^{34} bytes or 16 gigabytes. Let us also remember that if the memory size is quadratic, the algorithm runs in cubic time.

Note

When paths = False the algorithm saves roughly half of the memory as it does not have to maintain the matrix of predecessors. However, setting distances=False produces no such effect as the algorithm can not run without computing them. They will not be returned, but they will be stored while the method is running.

EXAMPLES:

Shortest paths in a small grid

sage: g = graphs.Grid2dGraph(2,2)
sage: from sage.graphs.base.c_graph import floyd_warshall
sage: print floyd_warshall(g)
{(0, 1): {(0, 1): None, (1, 0): (0, 0), (0, 0): (0, 1), (1, 1): (0, 1)},
(1, 0): {(0, 1): (0, 0), (1, 0): None, (0, 0): (1, 0), (1, 1): (1, 0)},
(0, 0): {(0, 1): (0, 0), (1, 0): (0, 0), (0, 0): None, (1, 1): (0, 1)},
(1, 1): {(0, 1): (1, 1), (1, 0): (1, 1), (0, 0): (0, 1), (1, 1): None}}

Checking the distances are correct

sage: g = graphs.Grid2dGraph(5,5)
sage: dist,path = floyd_warshall(g, distances = True)
sage: all( dist[u][v] == g.distance(u,v) for u in g for v in g )
True

Checking a random path is valid

sage: u,v = g.random_vertex(), g.random_vertex()
sage: p = [v]
sage: while p[0] != None:
...     p.insert(0,path[u][p[0]])
sage: len(p) == dist[u][v] + 2
True

Distances for all pairs of vertices in a diamond:

sage: g = graphs.DiamondGraph()
sage: floyd_warshall(g, paths = False, distances = True)
{0: {0: 0, 1: 1, 2: 1, 3: 2}, 
 1: {0: 1, 1: 0, 2: 1, 3: 1}, 
 2: {0: 1, 1: 1, 2: 0, 3: 1}, 
 3: {0: 2, 1: 1, 2: 1, 3: 0}}

TESTS:

Too large graphs:

sage: from sage.graphs.base.c_graph import floyd_warshall
sage: floyd_warshall(Graph(65536))
Traceback (most recent call last):
...
ValueError: The graph backend contains more than 65535 nodes

Table Of Contents

Previous topic

Graph database

Next topic

Fast sparse graphs

This Page