#include "misc.h"
#include "sets.h"
#include "smallgraph.h"
#include "queue.h"
#include <assert.h>

/*************************************************************************
**
**                            The Basics
**
*************************************************************************/

SMALL_GRAPH *SmallGraphAlloc(unsigned int n)
{
    static Boolean startup = 1;
    SMALL_GRAPH *G = Calloc(1, sizeof(SMALL_GRAPH));
    assert(n <= MAX_SSET);
    if(startup)
    {
	startup = 0;
	SetStartup();
    }
    G->n = n;
    return G;
}

SMALL_GRAPH *SmallGraphConnect(SMALL_GRAPH *G, int i, int j)
{   
    if(SmallGraphAreConnected(G, i, j))
	return G;
    SSetAdd(G->A[i], j);
    SSetAdd(G->A[j], i);
    ++G->degree[i];
    ++G->degree[j];
    return G;
}

SMALL_GRAPH *SmallGraphEdgesAllDelete(SMALL_GRAPH *G)
{
    int i;
    for(i=0; i < MAX_SSET; i++)
    {
	G->degree[i] = 0;
	SSetEmpty(G->A[i]);
    }
    return G;
}

SMALL_GRAPH *SmallGraphDisconnect(SMALL_GRAPH *G, int i, int j)
{   
    if(!SmallGraphAreConnected(G, i, j))
	return G;
    SSetDelete(G->A[i], j);
    SSetDelete(G->A[j], i);
    --G->degree[i];
    --G->degree[j];
    return G;
}

#ifndef SmallGraphAreConnected
Boolean SmallGraphAreConnected(SMALL_GRAPH *G, int i, int j)
{
    if(SSetIn(G->A[i],j))
    {
	assert(SSetIn(G->A[j],i));
	return true;
    }
    else
	return false;
}
#endif

void SmallGraphPrintAdjMatrix(FILE *fp, SMALL_GRAPH *G)
{
    int i, j;
    for(i=0; i<G->n; i++)
    {
	for(j=0; j<G->n - 1; j++)
	    fprintf(fp, "%d ", !!SmallGraphAreConnected(G,i,j));
	fprintf(fp, "%d\n", !!SmallGraphAreConnected(G,i,j));
    }
}


SMALL_GRAPH *SmallGraphReadAdjMatrix(FILE *fp)
{
    int i,j,n;
    SMALL_GRAPH *G;
    fscanf(fp, "%d", &n);
    G = SmallGraphAlloc(n);
    for(i=0; i<n; i++) for(j=0; j<n; j++)
    {
	int connected;
	if(fscanf(fp, "%d", &connected) != 1)
	    Fatal("too few vertices listed in file");
	if(connected)
	    SmallGraphConnect(G,i,j);
    }
    return G;
}


SMALL_GRAPH *SmallGraphComplement(SMALL_GRAPH *Gbar, SMALL_GRAPH *G)
{
    int i;
    SSET full = 0;
    if(!Gbar)
	Gbar = SmallGraphAlloc(G->n);

    for(i=0; i < G->n; i++)
	SSetAdd(full, i);

    Gbar->n = G->n;

    /* A tad tricky: we don't want to include self-loops, so we add them,
    ** and then XOR with the full row.
    */
    for(i=0; i < G->n; i++)
    {
	Gbar->A[i] = (G->A[i] | (SSET1 << i)) ^ full;
	Gbar->degree[i] = G->n - 1 - G->degree[i];
    }
    return Gbar;
}


SMALL_GRAPH *SmallGraphUnion(SMALL_GRAPH *dest, SMALL_GRAPH *G1, SMALL_GRAPH *G2)
{
    int i;

    if(G1->n != G2->n)
	return NULL;
    
    if(dest)
	dest->n = G1->n;
    else
	dest = SmallGraphAlloc(G1->n);

    for(i=0; i < G1->n; i++)
    {
	dest->A[i] = G1->A[i] | G2->A[i];
	dest->degree[i] = SSetCardinality(dest->A[i]);
    }

    return dest;
}


int SmallGraphBFS(SMALL_GRAPH *G, int root, int distance, int *nodeArray, int *distArray)
{
    QUEUE *BFSQ;
    int i, count = 0;

    assert(0 <= root && root < G->n);
    assert(distance >= 0);
    assert(nodeArray != NULL);
    assert(distArray != NULL);

    if(distance == 0) /* We could let the rest of the routine run, but why bother? */
    {
	nodeArray[0] = root;
	distArray[root] = 0;
	return 1;
    }

    for(i=0; i<G->n; i++)
	nodeArray[i] = distArray[i] = -1;

    distArray[root] = 0;
    BFSQ = QueueAlloc(G->n);
    QueuePut(BFSQ, (foint)root);
    while(QueueSize(BFSQ) > 0)
    {
	int v = QueueGet(BFSQ).i;

	/* At this point, distArray[v] should be assigned (when v was appended
	 * onto the queue), but v hasn't been "visited" or "counted" yet.
	 */

	assert(0 <= v && v < G->n);
	assert(0 <= distArray[v] && distArray[v] < G->n);

	assert(nodeArray[count] == -1);
	nodeArray[count] = v;
	count++;

	if(distArray[v] < distance) /* v's neighbors will be within BFS distance */
	{
	    int neighbor[MAX_SSET];
	    int j, numNeighbors = SSetToArray(neighbor, G->A[v]); /* This is the slow part, O(n) */
	    for(j=0; j < numNeighbors; j++)
		if(distArray[neighbor[j]] == -1) /* some of the neighbors might have already been visited */
		{
		    distArray[neighbor[j]] = distArray[v] + 1;
		    QueuePut(BFSQ, (foint)neighbor[j]);
		}
	}
    }
    QueueFree(BFSQ);
    return count;
}

SMALL_GRAPH *SmallGraphInduced(SMALL_GRAPH *Gv, SMALL_GRAPH *G, SSET V)
{
    unsigned array[MAX_SSET], nV = SSetToArray(array, V), i, j;
    SMALL_GRAPH GGv;
    if(Gv)
    {
	if(Gv == G) /* be careful, destination is same as source */
	    Gv = &GGv;
	Gv->n = nV;
	SmallGraphEdgesAllDelete(Gv);
    }
    else
	Gv = SmallGraphAlloc(nV);

    for(i=0; i < nV; i++) for(j=i+1; j < nV; j++)
	if(SmallGraphAreConnected(G, array[i], array[j]))
	    SmallGraphConnect(Gv, i, j);

    if(Gv == &GGv)
	*(Gv = G) = GGv;
    return Gv;
}


SMALL_GRAPH *SmallGraphInduced_NoVertexDelete(SMALL_GRAPH *Gv, SMALL_GRAPH *G, SSET V)
{
    unsigned array[MAX_SSET], nV = SSetToArray(array, V), i, j;
    SMALL_GRAPH GGv;
    if(Gv)
    {
	if(Gv == G) /* be careful, destination is same as source */
	    Gv = &GGv;
	Gv->n = G->n;
	SmallGraphEdgesAllDelete(Gv);
    }
    else
	Gv = SmallGraphAlloc(G->n);

    for(i=0; i < nV; i++) for(j=i+1; j < nV; j++)
	if(SmallGraphAreConnected(G, array[i], array[j]))
	    SmallGraphConnect(Gv, array[i], array[j]);
    if(Gv == &GGv)
	*(Gv = G) = GGv;
    return Gv;
}

/*
** A reasonably fast search for triangles.  If there exists a triangle,
** it'll return the SSET representation of it, otherwise return 0.
*/
#if 0
static int SmallGraphContainsK3(SMALL_GRAPH *G)
{
    int i,j;
    for(i=0; i < G->n; i++)
	for(j=i+1; j < G->n; j++)
	    if(SmallGraphAreConnected(G,i,j) && SmallGraphConnectingCausesK3(G,i,j))
	    {
		int k;
		for(k=0; k < G->n; k++)
		    if(SSetIn(SSetIntersect(G->A[i], G->A[j]), k))
			break;
		return (SSET1<<i)|(SSET1<<j)|(SSET1<<k);
	    }
    return 0;
}
#endif



/**************************************************************************
**
** These are the Clique and Indep set (exponential time) algorithms.
** They work for general graphs.
**
**************************************************************************/

SSET SmallGraph_IsCombClique(CLIQUE *c);

CLIQUE *SmallGraphKnFirst(SMALL_GRAPH *g, int k)
{
    int i, nDegk = g->n, prevNDegk;
    SSET setDegk;
    CLIQUE *c;

    assert(k <= g->n);
    if(k == 0)
	return NULL;
    
    c = (CLIQUE*)Calloc(1,sizeof(CLIQUE));
    c->G = SmallGraphAlloc(g->n);
    *(c->G) = *g;
    c->cliqueSize = k;

#if 1
    /*
    ** First reduce the potential number of vertices needed to check. A
    ** necessary condition for a vertex to be a member of a Kk is to have
    ** degree >= k-1, and all edges go to other vertices of degree >= k-1.
    ** Iterate inducing subgraphs 'til we cain't induce no more...
    */
    while(nDegk >= k)
    {
	prevNDegk = nDegk;  /* just to check for changes */
	nDegk = 0;
	SSetEmpty(setDegk);
	for(i=0; i < c->G->n; i++)
	    if(c->G->degree[i] >= k-1)
	    {
		++nDegk;
		SSetAdd(setDegk, i);
	    }
	if(nDegk == prevNDegk)  /* nobody eliminated */
	    break;
	SmallGraphInduced_NoVertexDelete(c->G, c->G, setDegk);
    }

    if(nDegk < k)
    {
	SmallGraphFree(c->G);
	Free(c);
	return NULL;
    }
    SSetToArray(c->inducedArray, setDegk);
#else
    nDegk = g->n;
    for(i=0; i < nDegk; i++)
	c->inducedArray[i] = i;
#endif
    c->combin = CombinZeroth(nDegk, k, c->combArray);
    if(SmallGraph_IsCombClique(c) || SmallGraphKnNext(c))
	return c;
    SmallGraphCliqueFree(c);
    return NULL;
}

CLIQUE *SmallGraphInFirst(SMALL_GRAPH *G, int n)
{
    static SMALL_GRAPH Gbar;  /* non re-entrant */
    SmallGraphComplement(&Gbar, G);
    return SmallGraphKnFirst(&Gbar, n);
}


/*
** We check to see if these particular k vertices are a clique by
** bitwise-anding together their adjacency lists (with self-loop added
** to each).  Return 0 of it's not a clique, otherwise return it's
** SSET representation.  Note that if the graph has changed since the
** previous call, we may miss some cliques.
*/
SSET SmallGraph_IsCombClique(CLIQUE *c)
{
    SSET result;
    int i;

    SSetEmpty(c->sset);
    for(i=0; i < c->cliqueSize; i++)
	SSetAdd(c->sset, c->inducedArray[c->combArray[i]]);

    result = c->sset;
    for(i=0; i < c->cliqueSize; i++)
    {
	result &= c->G->A[c->inducedArray[c->combArray[i]]] |
	    (SSET1 << c->inducedArray[c->combArray[i]]);    /* add self-loop */
	if(result != c->sset)
	    return 0;
    }
    return c->sset;
}

SSET SmallGraphKnNext(CLIQUE *c)
{
    while(CombinNext(c->combin))
	if((c->sset = SmallGraph_IsCombClique(c)))
	    return c->sset;
    return 0;
}

void SmallGraphCliqueFree(CLIQUE *c)
{
    SmallGraphFree(c->G);
    CombinFree(c->combin);
    Free(c);
}


Boolean SmallGraphKnContains(SMALL_GRAPH *G, int n)
{
    CLIQUE *c;
    assert(n <= G->n);
    if(n < 2)
	return true;
    c = SmallGraphKnFirst(G, n);
    if(c)
	SmallGraphCliqueFree(c);
    return c != NULL;
}

Boolean SmallGraphInContains(SMALL_GRAPH *G, int n)
{
    CLIQUE *c;
    assert(n <= G->n);
    if(n < 2)
	return true;
    c = SmallGraphInFirst(G, n);
    if(c)
	SmallGraphCliqueFree(c);
    return c != NULL;
}


/**************************************************************************
**
**  Graph Isomorphism
**
**************************************************************************/

static SMALL_GRAPH *isoG1, *isoG2;

static int _permutationIdentical(int n, int perm[n])
{
    int i, j;
    for(i=0; i<n; i++)
	if(isoG1->degree[i] != isoG2->degree[perm[i]])
	    return 0;

    for(i=0; i<n; i++) for(j=i+1; j<n; j++)
	/* The !GraphAreConnected is just to turn a bitstring into a boolean */
	if(!SmallGraphAreConnected(isoG1, i,j) !=
	    !SmallGraphAreConnected(isoG2, perm[i], perm[j]))
	    return 0;   /* non-isomorphic */
    return 1;   /* isomorphic! */
}

Boolean SmallGraphsIsomorphic(SMALL_GRAPH *G1, SMALL_GRAPH *G2)
{
    int i, n = G1->n, workArray1[n], workArray2[n];
    SSET degreeOnce;

    /*
    ** First some simple tests.
    */
    if(G1->n != G2->n)
	return false;
    
    if(n < 2)
	return true;

    /*
    ** Ensure each degree occurs the same number of times in each.
    */
    for(i=0; i<n; i++)
	workArray1[i] = workArray2[i] = 0;
    for(i=0; i<n; i++)
    {
	++workArray1[G1->degree[i]];
	++workArray2[G2->degree[i]];
    }
    for(i=0; i<n; i++)
	if(workArray1[i] != workArray2[i])
	    return false;

    /*
    ** Let degree d appear only once.  Then there is exactly one vertex
    ** v1 in G1 with degree d, and exactly one vertex v2 in G2 with degree d.
    ** G1 and G2 are isomorphic only if the neighborhood of v1 is isomorphic
    ** to the neighborhood of v2.
    */
    SSetEmpty(degreeOnce);
    for(i=0; i<n; i++)
	if(workArray1[i] == 1)
	    SSetAdd(degreeOnce, i);

    for(i=0; i<n; i++)
    {
	/* Find out if the degree of vertex i in G1 appears only once */
	if(SSetIn(degreeOnce, G1->degree[i]))
	{
	    int j, degree = G1->degree[i];
	    SMALL_GRAPH neighG1i, neighG2j;

	    /* find the (unique) vertex in G2 that has the same degree */
	    for(j=0; j < n; j++)
		if(G2->degree[j] == degree)
		    break;
	    assert(j < n);

	    SmallGraphInduced(&neighG1i, G1, G1->A[i]);
	    SmallGraphInduced(&neighG2j, G2, G2->A[j]);

	    /*
	    ** Note: this recursion works only as long as
	    ** _permutationIdentical doesn't call GraphsIsomorphic.
	    ** (if it does, isoG1 and isoG2 get messed up).
	    */

	    if(!SmallGraphsIsomorphic(&neighG1i, &neighG2j))
		return false;
	    /* Otherwise they *might* be isomorphic, so keep going */
	}
    }

    /*
    ** Oh well, fire up the exponential search.
    ** CombinAllPermutations will return 0 iff all permutations were
    ** tried; the function _permutationIdentical should return non-zero
    ** when it finds an identical permutation, and that non-zero value
    ** will be returned here, indicating an identical permutation was
    ** found, ie, that the graphs are isomorphic.
    */
    isoG1 = G1; isoG2 = G2;
    return !!CombinAllPermutations(n, _permutationIdentical);
}
