def remove(G, loops=None):
   if loops is None:
     loops = nx.cycle_basis(G)
   for l in loops:
     s1, s2 = l[0], l[1]
     if (s1, s2) in G.edges():
       G.remove_edge(s1, s2)
     elif (s2, s1) in G.edges():
       G.remove_edge(s2, s1)
   loops2 = nx.cycle_basis(G)
   print('New graph has %i loops' %len(loops2))
   return G

def correctLoops(code, loops_graph, charge_type):
    while nx.cycle_basis(loops_graph) != []:
        cycle = nx.cycle_basis(loops_graph)[0]
        loop = path.Path(cycle)
        for data in code.Primal.nodes():
            if loop.contains_points([data]) == [True]:
                charge = code.Primal.node[data]['charge'][charge_type]
                code.Primal.node[data]['charge'][charge_type] = (charge + 1)%2

        l = len(cycle)
        for i in range(l):
            n1, n2 = cycle[i], cycle[(i+1)%l]
            loops_graph.remove_edge(*(n1,n2))

    return code, loops_graph

def cycle_space (T,C):
  """Return a list of cycle graphs representing the fundamental cycles of the
     spanning tree T extended by the chords edges of the graph C

  Parameters
  ----------
  T: a tree graph
     A NetworkX graph.
  C: graph representing chords for the tree T.
    A NetworkX (multidi)graph.

  Returns
  -------
  Z : a list of cycles

  """
  Z     = list ();
  edges = C.edges_iter ();

  for idx,e in enumerate (edges):
    if T.has_edge(*e) or T.has_edge(*e[::-1]):
      Z.append (list (e))
    else:
      T.add_edge (*e)
      Z.append (nx.cycle_basis (nx.Graph(T))[0])
      T.remove_edge (*e)

  return Z

def main(prog, argv):
	parser = argparse.ArgumentParser(prog=prog)
	parser.add_argument('rows', metavar='LENGTH', type=int)
	parser.add_argument('cols', metavar='WIDTH', type=int)
	parser.add_argument('--verbose', '-v', action='store_true')
	parser.add_argument('--output-cb', '-C', metavar='PATH', type=str, help='generate .cyclebasis file')
	parser.add_argument('--output', '-o', type=str, required=True, help='.gpickle output file')

	args = parser.parse_args(argv)

	cellrows, cellcols = args.rows, args.cols

	if args.verbose:
		print('Generating graph and attributes')
	g = make_circuit(cellrows, cellcols)
	xys = full_xy_dict(cellrows, cellcols)
	measured_edge = battery_vertices()

	if args.verbose:
		print('Saving circuit')
	save_output(args.output, g, xys, measured_edge)

	if args.output_cb is not None:
		if args.verbose:
			print('Collecting desirable cycles')
		good = collect_good_cycles(g, cellrows, cellcols)

		assert validate_paths(g, good)

		if args.verbose:
			print('Generating fallback cycles')
		fallback = nx.cycle_basis(g)

		cyclebasis = build_cyclebasis_terminal(good, fallback, thorough=False, verbose=args.verbose)
		fileio.cycles.write_cycles(cyclebasis, args.output_cb)

def order_tables(dbdir_str):
    """
    Sorts tables in constraint order, primary tables first

    Constructs a constraint graph for the tables in the supplied directory,
    sorts it in topological order, and returns the table names in that
    order in a list.
    """

    graph = nx.DiGraph()

    for file_str in os.listdir(dbdir_str):
        if file_str.endswith(".table"):
            add_table_to_graph(dbdir_str, file_str, graph)

    ordered_tables = None
    try:
        ordered_tables = nx.topological_sort(graph)
    except nx.NetworkXUnfeasible:
        print >> sys.stderr, "Tables and their constraints do not form a DAG"
        # The NetworkX package does not include a function that finds cycles
        # in a directed graph.  Lacking that, report the undirected cycles.
        # Since our tables currently contain no cycles other than self-loops,
        # no need to implement a directed cycle finder.
        cycles = nx.cycle_basis(nx.Graph(graph))
        print >> sys.stderr, "Found possible cycles:"
        print >> sys.stderr, cycles
        sys.exit(1)

    return ordered_tables

def find_sidechains(molecule_graph):
    # Identify chiral atoms
    atoms = molecule_graph.atoms
    chiral_centres = chir.get_chiral_sets(atoms)
    # Identify sidechains (Ca-Cb-X), apart from proline and glycine.
    sidechains = {}
    # Detection of sidechains requires the multiple bonds be present in the atom graph.
    chir.multi_bonds(atoms)
    for k, v in chiral_centres.items():
        carbons = [atom for atom in v if atom.element == "C"]
        amides = [
            carbon
            for carbon in carbons
            if any([type(nb) == chir.GhostAtom and nb.element == "O" for nb in nx.neighbors(atoms, carbon)])
            and any([nb.element == "N" or nb.element == "O" for nb in nx.neighbors(atoms, carbon)])
        ]
        nbs_n = [nb for nb in v if nb.element == "N"]
        if amides and nbs_n:
            amide_bond = (k, amides[0])
            n_bond = (k, nbs_n[0])
            h_bond = (k, [h for h in nx.neighbors(atoms, k) if h.element == "H"][0])
            # Now find sidechains by cutting the Ca-C, Ca-N and Ca-H bonds
            atoms.remove_edges_from([amide_bond, n_bond, h_bond])
            sidechain_atoms = [
                atom
                for atom in [comp for comp in nx.connected_components(atoms) if k in comp][0]
                if type(atom) != chir.GhostAtom
            ]
            atoms.add_edges_from([amide_bond, n_bond, h_bond])
            if not any([k in cycle for cycle in nx.cycle_basis(atoms.subgraph(sidechain_atoms))]):
                sidechains[k] = atoms.subgraph(sidechain_atoms)
    chir.remove_ghost_atoms(atoms)
    return sidechains

def hardy_cross(G, n, init_guess=None):
    cycles = nx.cycle_basis(G)
    cycles = np.array([[tuple([cycles[j][i], cycles[j][i+1]]) if (i < len(cycles[j])-1) else tuple([cycles[j][i], cycles[j][0]]) for i in range(len(cycles[j]))] for j in range(len(cycles))])

    L = [G.node[i]['demand'] for i in G.node.keys()]
    edges = np.array(G.edges())
    edge_idx = np.full((len(G), len(G)), 9999, dtype=int)
    edge_idx[edges[:,0], edges[:,1]] = np.arange(len(G.edges()))
    edge_idx[edges[:,1], edges[:,0]] = np.arange(len(G.edges()))
    
    edge_dir = np.zeros((len(G), len(G)), dtype=int)
    edge_dir[edges[:,0], edges[:,1]] = 1
    edge_dir[edges[:,1], edges[:,0]] = -1

    if init_guess == None:
        init_guess = np.linalg.lstsq(nx.incidence_matrix(G, oriented=True).toarray(), L)[0]
    A = init_guess.copy().astype(float)
    for i in range(n):
        for u in cycles:
            R = np.array([G[j[0]][j[1]]['weight'] for j in u])
            D = np.array(edge_dir[u[:,0], u[:,1]])
            C = np.array(A[edge_idx[u[:,0], u[:,1]]])
            C = C*D
            dV = (R*C).sum()
            di = dV/R.sum()
	    C = (C - di)*D
	    A[edge_idx[u[:,0], u[:,1]]] = C
    return A

def main():

    g = nx.Graph()

    g.add_cycle(range(5))
    g.add_cycle(range(5, 10))
    g.add_node(20)

    # 循環しているノードのリストを返す
    print(nx.cycle_basis(g)) #=> [[1, 2, 3, 4, 0], [9, 8, 7, 6, 5]]

    # もう一つ循環を作ってみる
    g.add_edges_from([(0, 5), (3, 8)])

    # 循環しているノードが一つ増えている
    print(nx.cycle_basis(g)) #=> [[9, 8, 7, 6, 5],

def prune_graph(G):
    """
        Return a graph describing the loopy part of G, which is
        implicitly described by the list of cycles.
        The loopy part does not contain any

        (a) tree subgraphs of G
        (b) bridges of G

        Thus pruning may disconnect the graph into several 
        connected components.
    """
    cycles = nx.cycle_basis(G)
    pruned = G.copy()
    cycle_nodes = set(chain.from_iterable(cycles))
    cycle_edges = []

    for c in cycles:
        cycle = c + [c[0]]
        a, b = tee(cycle)
        next(b, None)
        edges = izip(a, b)

        cycle_edges.append(edges)
    
    all_cycle_edges = set(tuple(sorted(e)) \
            for e in chain.from_iterable(cycle_edges))
    # remove treelike components and bridges by removing all
    # edges not belonging to loops and then all nodes not
    # belonging to loops.
    pruned.remove_edges_from(e for e in G.edges_iter() \
            if (not tuple(sorted(e)) in all_cycle_edges))
    pruned.remove_nodes_from(n for n in G if not n in cycle_nodes)

    return pruned

def is_tree(graph):
	"""
	Checks whether graph is a tree.

	A graph is a tree if it is 
		(i) undirected - always, in our setting 
		(ii) connected, and 
		(iii) contains no cycles.

	Args:
		graph: networkx.Graph instance

	Return a boolean that indicates whether this graph is a tree.
	"""

	if len(graph.nodes()) == 0:
		return True

	if not nx.is_connected(graph):
		return False

	if len(nx.cycle_basis(graph)) > 0:
		return False

	return True

    def order_tables(self):
        """ Sorts tables in constraint order, least constrained tables first

        Constructs a constraint graph for the tables in the supplied database
        or databases directory, sorts it in topological order with the least
        constained tables first, and returns the table names in that order in
        a list.
        """

        self.add_all_tables_to_graph()

        # With NetworkX v1.3, topological_sort raises NetworkXUnfeasible if
        # the graph is not acyclic.  With v1.2, it returns None.
        ordered_tables = None
        try:
            ordered_tables = nx.topological_sort(self.graph)
        except:
            ordered_tables = None
        if not ordered_tables:
            # Fails only if the graph has cycles.
            # The NetworkX package does not include a function that finds
            # cycles in a directed graph.  Lacking that, report the undirected
            # cycles.  Since our tables currently contain no cycles other than
            # self-loops, no need to implement a directed cycle finder.
            # @ToDo: Add a warning to the wiki about not putting cycles in
            # table relationships.  It's never necessary -- show a proper
            # relationship hierarchy.
            cycles = nx.cycle_basis(nx.Graph(self.graph))
            errmsg = "Tables and their constraints do not form a DAG.\n" + \
                     "Found possible cycles:\n" + str(cycles)
            raise ValueError, errmsg

        return ordered_tables

    def __init__(self, graph):

        edges = {}
        vertices = {}
        faces = {}
        half_edges = {}

        cycles = nx.cycle_basis(graph)
        fi = 0
        for cycle in cycles:
            n = len(cycle)
            for i in range(n-1):
                e = (cycle[i], cycle[i+1])
                if e not in edges:
                    edges[e] = fi
                    twin_a = e[0], e[1]
                    twin_b = e[1], e[0]
                    if twin_a not in half_edges:
                        half_edges[twin_a] = fi
                    if twin_b not in half_edges:
                        half_edges[twin_b] = None
            e = cycle[n-1], cycle[0]
            if e not in edges:
                edges[e] = fi
            faces[fi] = e


            fi += 1

        self.edges = edges
        self.faces = faces
        self.half_edges = half_edges

def VizualizeGraph(G, param):
    import os

    try:
        import matplotlib

        matplotlib.use("Agg")

        try:
            os.mkdir(param.output_directory + "/graph_regions" + str(int(param.mean_ins_size)))
        except OSError:
            # directory is already created
            pass
        # CB = nx.connected_component_subgraphs(G)
        CB = nx.cycle_basis(G)
        print "NR of SubG: ", len(CB)
        counter = 1
        for cycle in CB:
            if len(cycle) >= 6:  # at leats 6 nodes if proper cycle (haplotypic region)
                subgraph = nx.Graph()
                subgraph.add_cycle(cycle)
                nx.draw(subgraph)
                matplotlib.pyplot.savefig(
                    param.output_directory
                    + "graph_regions"
                    + str(int(param.mean_ins_size))
                    + "/"
                    + str(counter)
                    + ".png"
                )
                matplotlib.pyplot.clf()
                counter += 1
    except ImportError:
        pass
    return ()

 def _get_foot_graph(self):
     radius = self.distance/2.
     G = nx.Graph()
     query = ("""
         WITH origin AS
             (SELECT point
                 FROM rnodes_intersections
                 WHERE id = '%d')
         SELECT end1, end2, distance, run_score
             FROM routing_edges_latlon, origin
             WHERE
                 ST_Dwithin(origin.point::geography, point1::geography, %f) OR
                 ST_Dwithin(origin.point::geography, point2::geography, %f);"""
              % (self.start_node, radius, radius))
     self._cur.execute(query)
     G.add_edges_from([[int(node[0]), int(node[1]),
                        {'dist': node[2], 'run_score': node[3]}]
                       for node in self._cur.fetchall()])
     G_cycles = nx.Graph()
     for cycle in nx.cycle_basis(G):
         G_cycles.add_cycle(cycle)
     for edge in G.edges():
         if not G_cycles.has_edge(*edge):
             G.remove_edge(*edge)
     return G

def findPath(schema_graph):
	'''Determines if the schema contains a closed curve.'''
	cycles = nx.cycle_basis(schema_graph)
	if len(cycles) > 0:
		return cycles
	else:
		raise nx.exception.NetworkXNoCycle('No closed cycle found.')

 def test_cycle_basis(self):
     G=self.G
     cy=networkx.cycle_basis(G,0)
     sort_cy= sorted( sorted(c) for c in cy )
     assert_equal(sort_cy, [[0,1,2,3],[0,1,6,7,8],[0,3,4,5]])
     cy=networkx.cycle_basis(G,1)
     sort_cy= sorted( sorted(c) for c in cy )
     assert_equal(sort_cy, [[0,1,2,3],[0,1,6,7,8],[0,3,4,5]])
     cy=networkx.cycle_basis(G,9)
     sort_cy= sorted( sorted(c) for c in cy )
     assert_equal(sort_cy, [[0,1,2,3],[0,1,6,7,8],[0,3,4,5]])
     # test disconnected graphs
     G.add_cycle(list("ABC"))
     cy=networkx.cycle_basis(G,9)
     sort_cy= sorted(sorted(c) for c in cy[:-1]) + [sorted(cy[-1])]
     assert_equal(sort_cy, [[0,1,2,3],[0,1,6,7,8],[0,3,4,5],['A','B','C']])

def calculate(g, voltage):
    edges_num = nx.number_of_edges(g)
    # sort nodes in edges
    edges = [edge if edge[0] < edge[1] else (edge[1], edge[0])
             for edge in nx.edges(g)]

    a = np.zeros((edges_num, edges_num))
    b = np.zeros((edges_num, 1))
    i = 0
    # first law
    for node in [node for node in nx.nodes(g) if node != 0]:
        for neighbor in nx.all_neighbors(g, node):
            edge = tuple(sorted((node, neighbor)))
            a[i][edges.index(edge)] = 1 if neighbor < node else -1
        i += 1
    # second law
    cycles = nx.cycle_basis(g, 0)
    for cycle in cycles:
        for j in range(0, len(cycle)):
            node = cycle[j]
            next_node = cycle[(j + 1) % len(cycle)]
            edge = tuple(sorted((node, next_node)))
            resistance = g[node][next_node]['weight']
            a[i][edges.index(edge)] = resistance if node < next_node else -resistance
        if 0 in cycle:
            b[i] = voltage
        i += 1
    # solve
    x = np.linalg.solve(a, b)
    for (x1, x2), res in zip(edges, x):
        g[x1][x2]['current'] = res[0]

def check_closure(g, node):
    h = g.subgraph(g[node])
    if nx.cycle_basis(h): return True
    sorted_nodes = sorted(h.degree().items(), key=operator.itemgetter(1))
    nodes = [n[0] for n in sorted_nodes[:2] if n[1] == 1]
    if len(nodes) == 2:
        if g.node[nodes[0]]['vertices'] > g.degree(nodes[0]) and \
           g.node[nodes[1]]['vertices'] > g.degree(nodes[1]):
            g.add_edge(nodes[0], nodes[1])
            return False

        if g.node[nodes[0]]['vertices'] > g.degree(nodes[0]) and \
            g.node[nodes[1]]['vertices'] == g.degree(nodes[1]):
            h = g.subgraph(g[nodes[1]])
            sorted_nodes = sorted(h.degree().items(), key=operator.itemgetter(1))
            xnode = [n[0] for n in sorted_nodes if n[1] == 1 and n[0] not in [nodes[0], nodes[1], node]].pop()
            for n in g[xnode]:
                if n in [nodes[1], node]: continue
                g.add_edge(nodes[0], n)
            g.remove_node(xnode)
            return False

        if g.node[nodes[0]]['vertices'] == g.degree(nodes[0]) and \
            g.node[nodes[1]]['vertices'] > g.degree(nodes[1]):
            h = g.subgraph(g[nodes[0]])
            sorted_nodes = sorted(h.degree().items(), key=operator.itemgetter(1))
            xnode = [n[0] for n in sorted_nodes if n[1] == 1 and n[0] not in [nodes[0], nodes[1], node]].pop()
            for n in g[xnode]:
                if n in [nodes[0], node]: continue
                g.add_edge(nodes[1], n)
            g.remove_node(xnode)
            return False

 def f22(self):
     start = 0
     cycle_list = nx.cycle_basis(self.G)
     res = len(cycle_list)
     stop = 0
     # self.feature_time.append(stop - start)
     return res

def planar_cycle_basis_nx(g, xs, ys):
	if g.is_directed():
		raise TypeError('Not implemented for directed graphs')
	if g.is_multigraph():
		raise TypeError('Not implemented for multi-graphs')

	if not (set(g) == set(xs) == set(ys)):
		raise ValueError('g, xs, ys must all share same set of vertices')

	my_g = nx.Graph()
	my_g.add_nodes_from(iter(g))
	my_g.add_edges_from(g.edges())

	nx.set_node_attributes(my_g, 'x', xs)
	nx.set_node_attributes(my_g, 'y', ys)

	# BAND-AID (HACK)
	# This algorithm has an unresolved issue with graphs that have at least
	#   two biconnected components, one of which is "inside" another
	#   (geometrically speaking).  Luckily, it is safe to look at each
	#   biconnected component separately, as all edges in a cycle always
	#   belong to the same biconnected component.
	result = []
	for subg in nx.biconnected_component_subgraphs(my_g):
		if len(subg) >= 3: # minimum possible size for cycles to exist
			result.extend(planar_cycle_basis_impl(subg))

	# Restore some confidence in the result...
	if len(result) != len(nx.cycle_basis(g)):
		raise RuntimeError(
			'planar_cycle_basis produced a result of incorrect '
			'length on the given graph! (does it have crossing edges?)'
		)

	return result