1// Copyright 2014 Sonia Keys
2// License MIT: http://opensource.org/licenses/MIT
3
4package graph
5
6import (
7	"errors"
8	"fmt"
9
10	"github.com/soniakeys/bits"
11)
12
13// undir_RO.go is code generated from undir_cg.go by directives in graph.go.
14// Editing undir_cg.go is okay.  It is the code generation source.
15// DO NOT EDIT undir_RO.go.
16// The RO means read only and it is upper case RO to slow you down a bit
17// in case you start to edit the file.
18//-------------------
19
20// Bipartite constructs an object indexing the bipartite structure of a graph.
21//
22// In a bipartite component, nodes can be partitioned into two sets, or
23// "colors," such that every edge in the component goes from one set to the
24// other.
25//
26// If the graph is bipartite, the method constructs and returns a new
27// Bipartite object as b and returns ok = true.
28//
29// If the component is not bipartite, a representative odd cycle as oc and
30// returns ok = false.
31//
32// In the case of a graph with mulitiple connected components, this method
33// provides no control over the color orientation by component.  See
34// Undirected.BipartiteComponent if this control is needed.
35//
36// There are equivalent labeled and unlabeled versions of this method.
37func (g LabeledUndirected) Bipartite() (b *LabeledBipartite, oc []NI, ok bool) {
38	c1 := bits.New(g.Order())
39	c2 := bits.New(g.Order())
40	r, _, _ := g.ConnectedComponentReps()
41	// accumulate n2 number of zero bits in c2 as number of one bits in n1
42	var n, n2 int
43	for _, r := range r {
44		ok, n, _, oc = g.BipartiteComponent(r, c1, c2)
45		if !ok {
46			return
47		}
48		n2 += n
49	}
50	return &LabeledBipartite{g, c2, n2}, nil, true
51}
52
53// BipartiteComponent analyzes the bipartite structure of a connected component
54// of an undirected graph.
55//
56// In a bipartite component, nodes can be partitioned into two sets, or
57// "colors," such that every edge in the component goes from one set to the
58// other.
59//
60// Argument n can be any representative node of the component to be analyzed.
61// Arguments c1 and c2 must be separate bits.Bits objects constructed to be
62// of length of the number of nodes of g.  These bitmaps are used in the
63// component traversal and the bits of the component must be zero when the
64// method is called.
65//
66// If the component is bipartite, BipartiteComponent populates bitmaps
67// c1 and c2 with the two-coloring of the component, always assigning the set
68// with representative node n to bitmap c1.  It returns b = true,
69// and also returns the number of bits set in c1 and c2 as n1 and n2
70// respectively.
71//
72// If the component is not bipartite, BipartiteComponent returns b = false
73// and a representative odd cycle as oc.
74//
75// See also method Bipartite.
76//
77// There are equivalent labeled and unlabeled versions of this method.
78func (g LabeledUndirected) BipartiteComponent(n NI, c1, c2 bits.Bits) (b bool, n1, n2 int, oc []NI) {
79	a := g.LabeledAdjacencyList
80	b = true
81	var open bool
82	var df func(n NI, c1, c2 *bits.Bits, n1, n2 *int)
83	df = func(n NI, c1, c2 *bits.Bits, n1, n2 *int) {
84		c1.SetBit(int(n), 1)
85		*n1++
86		for _, nb := range a[n] {
87			if c1.Bit(int(nb.To)) == 1 {
88				b = false
89				oc = []NI{nb.To, n}
90				open = true
91				return
92			}
93			if c2.Bit(int(nb.To)) == 1 {
94				continue
95			}
96			df(nb.To, c2, c1, n2, n1)
97			if b {
98				continue
99			}
100			switch {
101			case !open:
102			case n == oc[0]:
103				open = false
104			default:
105				oc = append(oc, n)
106			}
107			return
108		}
109	}
110	df(n, &c1, &c2, &n1, &n2)
111	if b {
112		return b, n1, n2, nil
113	}
114	return b, 0, 0, oc
115}
116
117// BronKerbosch1 finds maximal cliques in an undirected graph.
118//
119// The graph must not contain parallel edges or loops.
120//
121// See https://en.wikipedia.org/wiki/Clique_(graph_theory) and
122// https://en.wikipedia.org/wiki/Bron%E2%80%93Kerbosch_algorithm for background.
123//
124// This method implements the BronKerbosch1 algorithm of WP; that is,
125// the original algorithm without improvements.
126//
127// The method calls the emit argument for each maximal clique in g, as long
128// as emit returns true.  If emit returns false, BronKerbosch1 returns
129// immediately.
130//
131// There are equivalent labeled and unlabeled versions of this method.
132//
133// See also more sophisticated variants BronKerbosch2 and BronKerbosch3.
134func (g LabeledUndirected) BronKerbosch1(emit func(bits.Bits) bool) {
135	a := g.LabeledAdjacencyList
136	var f func(R, P, X bits.Bits) bool
137	f = func(R, P, X bits.Bits) bool {
138		switch {
139		case !P.AllZeros():
140			r2 := bits.New(len(a))
141			p2 := bits.New(len(a))
142			x2 := bits.New(len(a))
143			pf := func(n int) bool {
144				r2.Set(R)
145				r2.SetBit(n, 1)
146				p2.ClearAll()
147				x2.ClearAll()
148				for _, to := range a[n] {
149					if P.Bit(int(to.To)) == 1 {
150						p2.SetBit(int(to.To), 1)
151					}
152					if X.Bit(int(to.To)) == 1 {
153						x2.SetBit(int(to.To), 1)
154					}
155				}
156				if !f(r2, p2, x2) {
157					return false
158				}
159				P.SetBit(n, 0)
160				X.SetBit(n, 1)
161				return true
162			}
163			if !P.IterateOnes(pf) {
164				return false
165			}
166		case X.AllZeros():
167			return emit(R)
168		}
169		return true
170	}
171	var R, P, X bits.Bits
172	R = bits.New(len(a))
173	P = bits.New(len(a))
174	X = bits.New(len(a))
175	P.SetAll()
176	f(R, P, X)
177}
178
179// BKPivotMaxDegree is a strategy for BronKerbosch methods.
180//
181// To use it, take the method value (see golang.org/ref/spec#Method_values)
182// and pass it as the argument to BronKerbosch2 or 3.
183//
184// The strategy is to pick the node from P or X with the maximum degree
185// (number of edges) in g.  Note this is a shortcut from evaluating degrees
186// in P.
187//
188// There are equivalent labeled and unlabeled versions of this method.
189func (g LabeledUndirected) BKPivotMaxDegree(P, X bits.Bits) (p NI) {
190	// choose pivot u as highest degree node from P or X
191	a := g.LabeledAdjacencyList
192	maxDeg := -1
193	P.IterateOnes(func(n int) bool { // scan P
194		if d := len(a[n]); d > maxDeg {
195			p = NI(n)
196			maxDeg = d
197		}
198		return true
199	})
200	X.IterateOnes(func(n int) bool { // scan X
201		if d := len(a[n]); d > maxDeg {
202			p = NI(n)
203			maxDeg = d
204		}
205		return true
206	})
207	return
208}
209
210// BKPivotMinP is a strategy for BronKerbosch methods.
211//
212// To use it, take the method value (see golang.org/ref/spec#Method_values)
213// and pass it as the argument to BronKerbosch2 or 3.
214//
215// The strategy is to simply pick the first node in P.
216//
217// There are equivalent labeled and unlabeled versions of this method.
218func (g LabeledUndirected) BKPivotMinP(P, X bits.Bits) NI {
219	return NI(P.OneFrom(0))
220}
221
222// BronKerbosch2 finds maximal cliques in an undirected graph.
223//
224// The graph must not contain parallel edges or loops.
225//
226// See https://en.wikipedia.org/wiki/Clique_(graph_theory) and
227// https://en.wikipedia.org/wiki/Bron%E2%80%93Kerbosch_algorithm for background.
228//
229// This method implements the BronKerbosch2 algorithm of WP; that is,
230// the original algorithm plus pivoting.
231//
232// The argument is a pivot function that must return a node of P or X.
233// P is guaranteed to contain at least one node.  X is not.
234// For example see BKPivotMaxDegree.
235//
236// The method calls the emit argument for each maximal clique in g, as long
237// as emit returns true.  If emit returns false, BronKerbosch1 returns
238// immediately.
239//
240// There are equivalent labeled and unlabeled versions of this method.
241//
242// See also simpler variant BronKerbosch1 and more sophisticated variant
243// BronKerbosch3.
244func (g LabeledUndirected) BronKerbosch2(pivot func(P, X bits.Bits) NI, emit func(bits.Bits) bool) {
245	a := g.LabeledAdjacencyList
246	var f func(R, P, X bits.Bits) bool
247	f = func(R, P, X bits.Bits) bool {
248		switch {
249		case !P.AllZeros():
250			r2 := bits.New(len(a))
251			p2 := bits.New(len(a))
252			x2 := bits.New(len(a))
253			pnu := bits.New(len(a))
254			// compute P \ N(u).  next 5 lines are only difference from BK1
255			pnu.Set(P)
256			for _, to := range a[pivot(P, X)] {
257				pnu.SetBit(int(to.To), 0)
258			}
259			// remaining code like BK1
260			pf := func(n int) bool {
261				r2.Set(R)
262				r2.SetBit(n, 1)
263				p2.ClearAll()
264				x2.ClearAll()
265				for _, to := range a[n] {
266					if P.Bit(int(to.To)) == 1 {
267						p2.SetBit(int(to.To), 1)
268					}
269					if X.Bit(int(to.To)) == 1 {
270						x2.SetBit(int(to.To), 1)
271					}
272				}
273				if !f(r2, p2, x2) {
274					return false
275				}
276				P.SetBit(n, 0)
277				X.SetBit(n, 1)
278				return true
279			}
280			if !pnu.IterateOnes(pf) {
281				return false
282			}
283		case X.AllZeros():
284			return emit(R)
285		}
286		return true
287	}
288	R := bits.New(len(a))
289	P := bits.New(len(a))
290	X := bits.New(len(a))
291	P.SetAll()
292	f(R, P, X)
293}
294
295// BronKerbosch3 finds maximal cliques in an undirected graph.
296//
297// The graph must not contain parallel edges or loops.
298//
299// See https://en.wikipedia.org/wiki/Clique_(graph_theory) and
300// https://en.wikipedia.org/wiki/Bron%E2%80%93Kerbosch_algorithm for background.
301//
302// This method implements the BronKerbosch3 algorithm of WP; that is,
303// the original algorithm with pivoting and degeneracy ordering.
304//
305// The argument is a pivot function that must return a node of P or X.
306// P is guaranteed to contain at least one node.  X is not.
307// For example see BKPivotMaxDegree.
308//
309// The method calls the emit argument for each maximal clique in g, as long
310// as emit returns true.  If emit returns false, BronKerbosch1 returns
311// immediately.
312//
313// There are equivalent labeled and unlabeled versions of this method.
314//
315// See also simpler variants BronKerbosch1 and BronKerbosch2.
316func (g LabeledUndirected) BronKerbosch3(pivot func(P, X bits.Bits) NI, emit func(bits.Bits) bool) {
317	a := g.LabeledAdjacencyList
318	var f func(R, P, X bits.Bits) bool
319	f = func(R, P, X bits.Bits) bool {
320		switch {
321		case !P.AllZeros():
322			r2 := bits.New(len(a))
323			p2 := bits.New(len(a))
324			x2 := bits.New(len(a))
325			pnu := bits.New(len(a))
326			// compute P \ N(u).  next lines are only difference from BK1
327			pnu.Set(P)
328			for _, to := range a[pivot(P, X)] {
329				pnu.SetBit(int(to.To), 0)
330			}
331			// remaining code like BK2
332			pf := func(n int) bool {
333				r2.Set(R)
334				r2.SetBit(n, 1)
335				p2.ClearAll()
336				x2.ClearAll()
337				for _, to := range a[n] {
338					if P.Bit(int(to.To)) == 1 {
339						p2.SetBit(int(to.To), 1)
340					}
341					if X.Bit(int(to.To)) == 1 {
342						x2.SetBit(int(to.To), 1)
343					}
344				}
345				if !f(r2, p2, x2) {
346					return false
347				}
348				P.SetBit(n, 0)
349				X.SetBit(n, 1)
350				return true
351			}
352			if !pnu.IterateOnes(pf) {
353				return false
354			}
355		case X.AllZeros():
356			return emit(R)
357		}
358		return true
359	}
360	R := bits.New(len(a))
361	P := bits.New(len(a))
362	X := bits.New(len(a))
363	P.SetAll()
364	// code above same as BK2
365	// code below new to BK3
366	ord, _ := g.DegeneracyOrdering()
367	p2 := bits.New(len(a))
368	x2 := bits.New(len(a))
369	for _, n := range ord {
370		R.SetBit(int(n), 1)
371		p2.ClearAll()
372		x2.ClearAll()
373		for _, to := range a[n] {
374			if P.Bit(int(to.To)) == 1 {
375				p2.SetBit(int(to.To), 1)
376			}
377			if X.Bit(int(to.To)) == 1 {
378				x2.SetBit(int(to.To), 1)
379			}
380		}
381		if !f(R, p2, x2) {
382			return
383		}
384		R.SetBit(int(n), 0)
385		P.SetBit(int(n), 0)
386		X.SetBit(int(n), 1)
387	}
388}
389
390// ConnectedComponentBits returns a function that iterates over connected
391// components of g, returning a member bitmap for each.
392//
393// Each call of the returned function returns the order, arc size,
394// and bits of a connected component.  The underlying bits allocation is
395// the same for each call and is overwritten on subsequent calls.  Use or
396// save the bits before calling the function again.  The function returns
397// zeros after returning all connected components.
398//
399// There are equivalent labeled and unlabeled versions of this method.
400//
401// See also ConnectedComponentInts, ConnectedComponentReps, and
402// ConnectedComponentReps.
403func (g LabeledUndirected) ConnectedComponentBits() func() (order, arcSize int, bits bits.Bits) {
404	a := g.LabeledAdjacencyList
405	vg := bits.New(len(a)) // nodes visited in graph
406	vc := bits.New(len(a)) // nodes visited in current component
407	var order, arcSize int
408	var df func(NI)
409	df = func(n NI) {
410		vg.SetBit(int(n), 1)
411		vc.SetBit(int(n), 1)
412		order++
413		arcSize += len(a[n])
414		for _, nb := range a[n] {
415			if vg.Bit(int(nb.To)) == 0 {
416				df(nb.To)
417			}
418		}
419		return
420	}
421	var n int
422	return func() (o, ma int, b bits.Bits) {
423		for ; n < len(a); n++ {
424			if vg.Bit(n) == 0 {
425				vc.ClearAll()
426				order, arcSize = 0, 0
427				df(NI(n))
428				return order, arcSize, vc
429			}
430		}
431		return // return zeros signalling no more components
432	}
433}
434
435// ConnectedComponenInts returns a list of component numbers (ints) for each
436// node of graph g.
437//
438// The method assigns numbers to components 1-based, 1 through the number of
439// components.  Return value ci contains the component number for each node.
440// Return value nc is the number of components.
441//
442// There are equivalent labeled and unlabeled versions of this method.
443//
444// See also ConnectedComponentBits, ConnectedComponentLists, and
445// ConnectedComponentReps.
446func (g LabeledUndirected) ConnectedComponentInts() (ci []int, nc int) {
447	a := g.LabeledAdjacencyList
448	ci = make([]int, len(a))
449	var df func(NI)
450	df = func(nd NI) {
451		ci[nd] = nc
452		for _, to := range a[nd] {
453			if ci[to.To] == 0 {
454				df(to.To)
455			}
456		}
457		return
458	}
459	for nd := range a {
460		if ci[nd] == 0 {
461			nc++
462			df(NI(nd))
463		}
464	}
465	return
466}
467
468// ConnectedComponentLists returns a function that iterates over connected
469// components of g, returning the member list of each.
470//
471// Each call of the returned function returns a node list of a connected
472// component and the arc size of the component.  The returned function returns
473// nil, 0 after returning all connected components.
474//
475// There are equivalent labeled and unlabeled versions of this method.
476//
477// See also ConnectedComponentBits, ConnectedComponentInts, and
478// ConnectedComponentReps.
479func (g LabeledUndirected) ConnectedComponentLists() func() (nodes []NI, arcSize int) {
480	a := g.LabeledAdjacencyList
481	vg := bits.New(len(a)) // nodes visited in graph
482	var l []NI             // accumulated node list of current component
483	var ma int             // accumulated arc size of current component
484	var df func(NI)
485	df = func(n NI) {
486		vg.SetBit(int(n), 1)
487		l = append(l, n)
488		ma += len(a[n])
489		for _, nb := range a[n] {
490			if vg.Bit(int(nb.To)) == 0 {
491				df(nb.To)
492			}
493		}
494		return
495	}
496	var n int
497	return func() ([]NI, int) {
498		for ; n < len(a); n++ {
499			if vg.Bit(n) == 0 {
500				l, ma = nil, 0
501				df(NI(n))
502				return l, ma
503			}
504		}
505		return nil, 0
506	}
507}
508
509// ConnectedComponentReps returns a representative node from each connected
510// component of g.
511//
512// Returned is a slice with a single representative node from each connected
513// component and also parallel slices with the orders and arc sizes
514// in the corresponding components.
515//
516// This is fairly minimal information describing connected components.
517// From a representative node, other nodes in the component can be reached
518// by depth first traversal for example.
519//
520// There are equivalent labeled and unlabeled versions of this method.
521//
522// See also ConnectedComponentBits and ConnectedComponentLists which can
523// collect component members in a single traversal, and IsConnected which
524// is an even simpler boolean test.
525func (g LabeledUndirected) ConnectedComponentReps() (reps []NI, orders, arcSizes []int) {
526	a := g.LabeledAdjacencyList
527	c := bits.New(len(a))
528	var o, ma int
529	var df func(NI)
530	df = func(n NI) {
531		c.SetBit(int(n), 1)
532		o++
533		ma += len(a[n])
534		for _, nb := range a[n] {
535			if c.Bit(int(nb.To)) == 0 {
536				df(nb.To)
537			}
538		}
539		return
540	}
541	for n := range a {
542		if c.Bit(n) == 0 {
543			o, ma = 0, 0
544			df(NI(n))
545			reps = append(reps, NI(n))
546			orders = append(orders, o)
547			arcSizes = append(arcSizes, ma)
548		}
549	}
550	return
551}
552
553// Copy makes a deep copy of g.
554// Copy also computes the arc size ma, the number of arcs.
555//
556// There are equivalent labeled and unlabeled versions of this method.
557func (g LabeledUndirected) Copy() (c LabeledUndirected, ma int) {
558	l, s := g.LabeledAdjacencyList.Copy()
559	return LabeledUndirected{l}, s
560}
561
562// Degeneracy is a measure of dense subgraphs within a graph.
563//
564// See Wikipedia https://en.wikipedia.org/wiki/Degeneracy_(graph_theory)
565//
566// See also method DegeneracyOrdering which returns a degeneracy node
567// ordering and k-core breaks.
568//
569// There are equivalent labeled and unlabeled versions of this method.
570func (g LabeledUndirected) Degeneracy() (k int) {
571	a := g.LabeledAdjacencyList
572	// WP algorithm, attributed to Matula and Beck.
573	L := bits.New(len(a))
574	d := make([]int, len(a))
575	var D [][]NI
576	for v, nb := range a {
577		dv := len(nb)
578		d[v] = dv
579		for len(D) <= dv {
580			D = append(D, nil)
581		}
582		D[dv] = append(D[dv], NI(v))
583	}
584	for range a {
585		// find a non-empty D
586		i := 0
587		for len(D[i]) == 0 {
588			i++
589		}
590		// k is max(i, k)
591		if i > k {
592			k = i
593		}
594		// select from D[i]
595		Di := D[i]
596		last := len(Di) - 1
597		v := Di[last]
598		// Add v to ordering, remove from Di
599		L.SetBit(int(v), 1)
600		D[i] = Di[:last]
601		// move neighbors
602		for _, nb := range a[v] {
603			if L.Bit(int(nb.To)) == 1 {
604				continue
605			}
606			dn := d[nb.To] // old number of neighbors of nb
607			Ddn := D[dn]   // nb is in this list
608			// remove it from the list
609			for wx, w := range Ddn {
610				if w == nb.To {
611					last := len(Ddn) - 1
612					Ddn[wx], Ddn[last] = Ddn[last], Ddn[wx]
613					D[dn] = Ddn[:last]
614				}
615			}
616			dn-- // new number of neighbors
617			d[nb.To] = dn
618			// re--add it to it's new list
619			D[dn] = append(D[dn], nb.To)
620		}
621	}
622	return
623}
624
625// DegeneracyOrdering computes degeneracy node ordering and k-core breaks.
626//
627// See Wikipedia https://en.wikipedia.org/wiki/Degeneracy_(graph_theory)
628//
629// In return value ordering, nodes are ordered by their "coreness" as
630// defined at https://en.wikipedia.org/wiki/Degeneracy_(graph_theory)#k-Cores.
631//
632// Return value kbreaks indexes ordering by coreness number.  len(kbreaks)
633// will be one more than the graph degeneracy as returned by the Degeneracy
634// method.  If degeneracy is d, d = len(kbreaks) - 1, kbreaks[d] is the last
635// value in kbreaks and ordering[:kbreaks[d]] contains nodes of the d-cores
636// of the graph.  kbreaks[0] is always the number of nodes in g as all nodes
637// are in in a 0-core.
638//
639// Note that definitions of "k-core" differ on whether a k-core must be a
640// single connected component.  This method does not resolve individual
641// connected components.
642//
643// See also method Degeneracy which returns just the degeneracy number.
644//
645// There are equivalent labeled and unlabeled versions of this method.
646func (g LabeledUndirected) DegeneracyOrdering() (ordering []NI, kbreaks []int) {
647	a := g.LabeledAdjacencyList
648	// WP algorithm
649	k := 0
650	ordering = make([]NI, len(a))
651	kbreaks = []int{len(a)}
652	L := bits.New(len(a))
653	d := make([]int, len(a))
654	var D [][]NI
655	for v, nb := range a {
656		dv := len(nb)
657		d[v] = dv
658		for len(D) <= dv {
659			D = append(D, nil)
660		}
661		D[dv] = append(D[dv], NI(v))
662	}
663	for ox := len(a) - 1; ox >= 0; ox-- {
664		// find a non-empty D
665		i := 0
666		for len(D[i]) == 0 {
667			i++
668		}
669		// k is max(i, k)
670		if i > k {
671			for len(kbreaks) <= i {
672				kbreaks = append(kbreaks, ox+1)
673			}
674			k = i
675		}
676		// select from D[i]
677		Di := D[i]
678		last := len(Di) - 1
679		v := Di[last]
680		// Add v to ordering, remove from Di
681		ordering[ox] = v
682		L.SetBit(int(v), 1)
683		D[i] = Di[:last]
684		// move neighbors
685		for _, nb := range a[v] {
686			if L.Bit(int(nb.To)) == 1 {
687				continue
688			}
689			dn := d[nb.To] // old number of neighbors of nb
690			Ddn := D[dn]   // nb is in this list
691			// remove it from the list
692			for wx, w := range Ddn {
693				if w == nb.To {
694					last := len(Ddn) - 1
695					Ddn[wx], Ddn[last] = Ddn[last], Ddn[wx]
696					D[dn] = Ddn[:last]
697				}
698			}
699			dn-- // new number of neighbors
700			d[nb.To] = dn
701			// re--add it to it's new list
702			D[dn] = append(D[dn], nb.To)
703		}
704	}
705	//for i, j := 0, k; i < j; i, j = i+1, j-1 {
706	//	kbreaks[i], kbreaks[j] = kbreaks[j], kbreaks[i]
707	//}
708	return
709}
710
711// Degree for undirected graphs, returns the degree of a node.
712//
713// The degree of a node in an undirected graph is the number of incident
714// edges, where loops count twice.
715//
716// If g is known to be loop-free, the result is simply equivalent to len(g[n]).
717// See handshaking lemma example at AdjacencyList.ArcSize.
718//
719// There are equivalent labeled and unlabeled versions of this method.
720func (g LabeledUndirected) Degree(n NI) int {
721	to := g.LabeledAdjacencyList[n]
722	d := len(to) // just "out" degree,
723	for _, to := range to {
724		if to.To == n {
725			d++ // except loops count twice
726		}
727	}
728	return d
729}
730
731// DegreeCentralization returns the degree centralization metric of a graph.
732//
733// Degree of a node is one measure of node centrality and is directly
734// available from the adjacency list representation.  This allows degree
735// centralization for the graph to be very efficiently computed.
736//
737// The value returned is from 0 to 1 inclusive for simple graphs of three or
738// more nodes.  As a special case, 0 is returned for graphs of two or fewer
739// nodes.  The value returned can be > 1 for graphs with loops or parallel
740// edges.
741//
742// There are equivalent labeled and unlabeled versions of this method.
743func (g LabeledUndirected) DegreeCentralization() float64 {
744	a := g.LabeledAdjacencyList
745	if len(a) <= 2 {
746		return 0
747	}
748	var max, sum int
749	for _, to := range a {
750		if len(to) > max {
751			max = len(to)
752		}
753		sum += len(to)
754	}
755	return float64(len(a)*max-sum) / float64((len(a)-1)*(len(a)-2))
756}
757
758// Density returns density for a simple graph.
759//
760// See also Density function.
761//
762// There are equivalent labeled and unlabeled versions of this method.
763func (g LabeledUndirected) Density() float64 {
764	return Density(g.Order(), g.Size())
765}
766
767// Eulerian scans an undirected graph to determine if it is Eulerian.
768//
769// If the graph represents an Eulerian cycle, it returns -1, -1, nil.
770//
771// If the graph does not represent an Eulerian cycle but does represent an
772// Eulerian path, it returns the two end nodes of the path, and nil.
773//
774// Otherwise it returns an error.
775//
776// See also method EulerianStart, which short-circuits as soon as it finds
777// a node that must be a start or end node of an Eulerian path.
778//
779// There are equivalent labeled and unlabeled versions of this method.
780func (g LabeledUndirected) Eulerian() (end1, end2 NI, err error) {
781	end1 = -1
782	end2 = -1
783	for n := range g.LabeledAdjacencyList {
784		switch {
785		case g.Degree(NI(n))%2 == 0:
786		case end1 < 0:
787			end1 = NI(n)
788		case end2 < 0:
789			end2 = NI(n)
790		default:
791			err = errors.New("non-Eulerian")
792			return
793		}
794	}
795	return
796}
797
798// EulerianCycle finds an Eulerian cycle in an undirected multigraph.
799//
800// * If g has no nodes, result is nil, nil.
801//
802// * If g is Eulerian, result is an Eulerian cycle with err = nil.
803// The first element of the result represents only a start node.
804// The remaining elements represent the half arcs of the cycle.
805//
806// * Otherwise, result is nil, with a non-nil error giving a reason the graph
807// is not Eulerian.
808//
809// Internally, EulerianCycle copies the entire graph g.
810// See EulerianCycleD for a more space efficient version.
811//
812// There are equivalent labeled and unlabeled versions of this method.
813func (g LabeledUndirected) EulerianCycle() ([]Half, error) {
814	c, _ := g.Copy()
815	return c.EulerianCycleD(c.Size())
816}
817
818// EulerianCycleD finds an Eulerian cycle in an undirected multigraph.
819//
820// EulerianCycleD is destructive on its receiver g.  See EulerianCycle for
821// a non-destructive version.
822//
823// Parameter m must be the size of the undirected graph -- the
824// number of edges.  Use Undirected.Size if the size is unknown.
825//
826// * If g has no nodes, result is nil, nil.
827//
828// * If g is Eulerian, result is an Eulerian cycle with err = nil.
829// The first element of the result represents only a start node.
830// The remaining elements represent the half arcs of the cycle.
831//
832// * Otherwise, result is nil, with a non-nil error giving a reason the graph
833// is not Eulerian.
834//
835// There are equivalent labeled and unlabeled versions of this method.
836func (g LabeledUndirected) EulerianCycleD(m int) ([]Half, error) {
837	if g.Order() == 0 {
838		return nil, nil
839	}
840	e := newLabEulerian(g.LabeledAdjacencyList, m)
841	e.p[0] = Half{0, -1}
842	for e.s >= 0 {
843		v := e.top()
844		if err := e.pushUndir(); err != nil {
845			return nil, err
846		}
847		if e.top().To != v.To {
848			return nil, errors.New("not Eulerian")
849		}
850		e.keep()
851	}
852	if !e.uv.AllZeros() {
853		return nil, errors.New("not strongly connected")
854	}
855	return e.p, nil
856}
857
858// EulerianPath finds an Eulerian path in an undirected multigraph.
859//
860// * If g has no nodes, result is nil, nil.
861//
862// * If g has an Eulerian path, result is an Eulerian path with err = nil.
863// The first element of the result represents only a start node.
864// The remaining elements represent the half arcs of the path.
865//
866// * Otherwise, result is nil, with a non-nil error giving a reason the graph
867// is not Eulerian.
868//
869// Internally, EulerianPath copies the entire graph g.
870// See EulerianPathD for a more space efficient version.
871//
872// There are equivalent labeled and unlabeled versions of this method.
873func (g LabeledUndirected) EulerianPath() ([]Half, error) {
874	c, _ := g.Copy()
875	start := c.EulerianStart()
876	if start < 0 {
877		start = 0
878	}
879	return c.EulerianPathD(c.Size(), start)
880}
881
882// EulerianPathD finds an Eulerian path in a undirected multigraph.
883//
884// EulerianPathD is destructive on its receiver g.  See EulerianPath for
885// a non-destructive version.
886//
887// Argument m must be the correct size, or number of edges in g.
888// Argument start must be a valid start node for the path.
889//
890// * If g has no nodes, result is nil, nil.
891//
892// * If g has an Eulerian path starting at start, result is an Eulerian path
893// with err = nil.
894// The first element of the result represents only a start node.
895// The remaining elements represent the half arcs of the path.
896//
897// * Otherwise, result is nil, with a non-nil error giving a reason the graph
898// is not Eulerian.
899//
900// There are equivalent labeled and unlabeled versions of this method.
901func (g LabeledUndirected) EulerianPathD(m int, start NI) ([]Half, error) {
902	if g.Order() == 0 {
903		return nil, nil
904	}
905	e := newLabEulerian(g.LabeledAdjacencyList, m)
906	e.p[0] = Half{start, -1}
907	// unlike EulerianCycle, the first path doesn't have to be a cycle.
908	if err := e.pushUndir(); err != nil {
909		return nil, err
910	}
911	e.keep()
912	for e.s >= 0 {
913		start = e.top().To
914		e.push()
915		// paths after the first must be cycles though
916		// (as long as there are nodes on the stack)
917		if e.top().To != start {
918			return nil, errors.New("no Eulerian path")
919		}
920		e.keep()
921	}
922	if !e.uv.AllZeros() {
923		return nil, errors.New("no Eulerian path")
924	}
925	return e.p, nil
926}
927
928// EulerianStart finds a candidate start node for an Eulerian path.
929//
930// A graph representing an Eulerian path can have two nodes with odd degree.
931// If it does, these must be the end nodes of the path.  EulerianEnd scans
932// for a node with an odd degree, returning immediately with the first one
933// it finds.
934//
935// If the scan completes without finding a node with odd degree the method
936// returns -1.
937//
938// See also method Eulerian, which completely validates a graph as representing
939// an Eulerian path.
940//
941// There are equivalent labeled and unlabeled versions of this method.
942func (g LabeledUndirected) EulerianStart() NI {
943	for n := range g.LabeledAdjacencyList {
944		if g.Degree(NI(n))%2 != 0 {
945			return NI(n)
946		}
947	}
948	return -1
949}
950
951// AddNode maps a node in a supergraph to a subgraph node.
952//
953// Argument p must be an NI in supergraph s.Super.  AddNode panics if
954// p is not a valid node index of s.Super.
955//
956// AddNode is idempotent in that it does not add a new node to the subgraph if
957// a subgraph node already exists mapped to supergraph node p.
958//
959// The mapped subgraph NI is returned.
960func (s *LabeledUndirectedSubgraph) AddNode(p NI) (b NI) {
961	if int(p) < 0 || int(p) >= s.Super.Order() {
962		panic(fmt.Sprint("AddNode: NI ", p, " not in supergraph"))
963	}
964	if b, ok := s.SubNI[p]; ok {
965		return b
966	}
967	a := s.LabeledUndirected.LabeledAdjacencyList
968	b = NI(len(a))
969	s.LabeledUndirected.LabeledAdjacencyList = append(a, nil)
970	s.SuperNI = append(s.SuperNI, p)
971	s.SubNI[p] = b
972	return
973}
974
975// InduceList constructs a node-induced subgraph.
976//
977// The subgraph is induced on receiver graph g.  Argument l must be a list of
978// NIs in receiver graph g.  Receiver g becomes the supergraph of the induced
979// subgraph.
980//
981// Duplicate NIs are allowed in list l.  The duplicates are effectively removed
982// and only a single corresponding node is created in the subgraph.  Subgraph
983// NIs are mapped in the order of list l, execpt for ignoring duplicates.
984// NIs in l that are not in g will panic.
985//
986// Returned is the constructed Subgraph object containing the induced subgraph
987// and the mappings to the supergraph.
988func (g *LabeledUndirected) InduceList(l []NI) *LabeledUndirectedSubgraph {
989	sub, sup := mapList(l)
990	return &LabeledUndirectedSubgraph{
991		Super:   g,
992		SubNI:   sub,
993		SuperNI: sup,
994		LabeledUndirected: LabeledUndirected{
995			g.LabeledAdjacencyList.induceArcs(sub, sup),
996		}}
997}
998
999// InduceBits constructs a node-induced subgraph.
1000//
1001// The subgraph is induced on receiver graph g.  Argument t must be a bitmap
1002// representing NIs in receiver graph g.  Receiver g becomes the supergraph
1003// of the induced subgraph.  NIs in t that are not in g will panic.
1004//
1005// Returned is the constructed Subgraph object containing the induced subgraph
1006// and the mappings to the supergraph.
1007func (g *LabeledUndirected) InduceBits(t bits.Bits) *LabeledUndirectedSubgraph {
1008	sub, sup := mapBits(t)
1009	return &LabeledUndirectedSubgraph{
1010		Super:   g,
1011		SubNI:   sub,
1012		SuperNI: sup,
1013		LabeledUndirected: LabeledUndirected{
1014			g.LabeledAdjacencyList.induceArcs(sub, sup),
1015		}}
1016}
1017
1018// IsConnected tests if an undirected graph is a single connected component.
1019//
1020// There are equivalent labeled and unlabeled versions of this method.
1021//
1022// See also ConnectedComponentReps for a method returning more information.
1023func (g LabeledUndirected) IsConnected() bool {
1024	a := g.LabeledAdjacencyList
1025	if len(a) == 0 {
1026		return true
1027	}
1028	b := bits.New(len(a))
1029	var df func(NI)
1030	df = func(n NI) {
1031		b.SetBit(int(n), 1)
1032		for _, to := range a[n] {
1033			if b.Bit(int(to.To)) == 0 {
1034				df(to.To)
1035			}
1036		}
1037	}
1038	df(0)
1039	return b.AllOnes()
1040}
1041
1042// IsTree identifies trees in undirected graphs.
1043//
1044// Return value isTree is true if the connected component reachable from root
1045// is a tree.  Further, return value allTree is true if the entire graph g is
1046// connected.
1047//
1048// There are equivalent labeled and unlabeled versions of this method.
1049func (g LabeledUndirected) IsTree(root NI) (isTree, allTree bool) {
1050	a := g.LabeledAdjacencyList
1051	v := bits.New(len(a))
1052	v.SetAll()
1053	var df func(NI, NI) bool
1054	df = func(fr, n NI) bool {
1055		if v.Bit(int(n)) == 0 {
1056			return false
1057		}
1058		v.SetBit(int(n), 0)
1059		for _, to := range a[n] {
1060			if to.To != fr && !df(n, to.To) {
1061				return false
1062			}
1063		}
1064		return true
1065	}
1066	v.SetBit(int(root), 0)
1067	for _, to := range a[root] {
1068		if !df(root, to.To) {
1069			return false, false
1070		}
1071	}
1072	return true, v.AllZeros()
1073}
1074
1075// Size returns the number of edges in g.
1076//
1077// See also ArcSize and AnyLoop.
1078func (g LabeledUndirected) Size() int {
1079	m2 := 0
1080	for fr, to := range g.LabeledAdjacencyList {
1081		m2 += len(to)
1082		for _, to := range to {
1083			if to.To == NI(fr) {
1084				m2++
1085			}
1086		}
1087	}
1088	return m2 / 2
1089}
1090
1091// Density returns edge density of a bipartite graph.
1092//
1093// Edge density is number of edges over maximum possible number of edges.
1094// Maximum possible number of edges in a bipartite graph is number of
1095// nodes of one color times number of nodes of the other color.
1096func (g LabeledBipartite) Density() float64 {
1097	a := g.LabeledUndirected.LabeledAdjacencyList
1098	s := 0
1099	g.Color.IterateOnes(func(n int) bool {
1100		s += len(a[n])
1101		return true
1102	})
1103	return float64(s) / float64(g.N0*(len(a)-g.N0))
1104}
1105
1106// PermuteBiadjacency permutes a bipartite graph in place so that a prefix
1107// of the adjacency list encodes a biadjacency matrix.
1108//
1109// The permutation applied is returned.  This would be helpful in referencing
1110// any externally stored node information.
1111//
1112// The biadjacency matrix is encoded as the prefix AdjacencyList[:g.N0].
1113// Note though that this slice does not represent a valid complete
1114// AdjacencyList.  BoundsOk would return false, for example.
1115//
1116// In adjacency list terms, the result of the permutation is that nodes of
1117// the prefix only have arcs to the suffix and nodes of the suffix only have
1118// arcs to the prefix.
1119func (g LabeledBipartite) PermuteBiadjacency() []int {
1120	p := make([]int, g.Order())
1121	i := 0
1122	g.Color.IterateZeros(func(n int) bool {
1123		p[n] = i
1124		i++
1125		return true
1126	})
1127	g.Color.IterateOnes(func(n int) bool {
1128		p[n] = i
1129		i++
1130		return true
1131	})
1132	g.Permute(p)
1133	g.Color.ClearAll()
1134	for i := g.N0; i < g.Order(); i++ {
1135		g.Color.SetBit(i, 1)
1136	}
1137	return p
1138}
1139