separator.c

//
// Created by Stephane on 10/03/2020.
//

#include <assert.h>
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>


#include "components.h"
#include "graph.h"
#include "main.h"
#include "separator.h"
#include "utils.h"
#include "heap.h"



// Options separation
int chooseSeedNeighborMaxDegree = 0;

int separatorJustMemorizeReceptorSize = 1;
int separatorJustMemoSizeMin = 10000;

int verifyTheSeparator = 0;


int useMinPartitions = 2; // 0: maintain heap, 1-2: just put min elements at the beginning.
// 1: maintain mins   2: calculate mins only when necessary



//
// Datas structures
//

int *nbNeighborsInB;
int *nbNeighborsInA;
int nbEdgesInA, nbEdgesInB;

int *dateABDisconnected; // Current date for AB Diconnected vertices of C (date=number of current separe call)
int *nbNeighborsABDisconnected; // For AB Disconnected vertices : Number of neighbors in C which are AB disconnected

int *priorities;

int *verticesToThrow;

// Copies to memorize nb neighbors and not recalculate at each separation run
int * nbNInDestCopy; // To build the best separator
int * nbNInABCopy; // To memorize the number of neighbors in A an B of C vertices for best separator
int * nbNeighborsInS;
int nbEACopy, nbEdgesInS;
int minDegree, nbMinD; // minimal degree of vertices of S, initialised in initializeNbNeighbors()
int maxDegree, nbMaxD;


int nbCallsSepare = 0;

int nbCallsSearchSeparator = 0;

Separator sep = NULL;
Separator bestSep = NULL;
double bestEval;
int bestIsUninitialized;


int sizeMin;
int sizeMax;


// Connected components
int *component;


//
// Warning: the preliminary call to initializeNbNeighbors() place first in the lists of neighbors
// those that are in the current set S[]. It is then unnecessary to verify each time that the vertex
// belongs to the current set.

#ifdef STATS_SEPARATION
clock_t *timesSep = NULL;
int *nbSearches;
int *sumSizes;
#endif

int searchSeparator(int *S, int n, SET V, Graph g, Separator theSeparator, int nTries, int nFlushes, int depth) {

#ifdef STATS_SEPARATION
    if (timesSep == NULL) {
        timesSep = calloc(10000, sizeof(clock_t));
        nbSearches = calloc(10000, sizeof(int));
        sumSizes = calloc(10000, sizeof(int));
    }
#endif

    if (sep == NULL) sep = newSeparator(g->n, g);
    if (bestSep == NULL) bestSep = newSeparator(g->n, g);

    nbCallsSearchSeparator ++;

    theSeparator->C->n = g->n+1;
    bestEval = -1.0;
    bestIsUninitialized = 1;

    sizeMin = n * ratioMin / 100;
    if (sizeMin == 0) sizeMin = 1;

    sizeMax = n * ratioMax / 100;

    if ((g->n > 300000) && (depth < 6)) nbRunsSeparation = 1; // avoids WA ?

    //saveNbNeighbors(S, n);

    if (modeEvalAtRandom) {
        if (rand()%2 == 1) modeEvalSeparator = EVAL_CARD;
        else modeEvalSeparator = rand() % 4;
    }

    if (0) {
        // test for heur_105: much better when equal to 0
        perForceChoiceNotInC = 0;
        perChooseInCAtRandom = 0;
    }

    clock_t start = clock();

    for (int i = 0; i < nTries; i ++) {

        randomizePriorities(S, n, g);

        initSeparator(V, S, n, g, sep);

        separe(S, n, V, g, sep, nFlushes, modeEvalSeparator);

        if ((stopSearch) || ((clock()-start)/CLOCKS_PER_SEC > maxTimeSeparation/(depth+1))) break;
    }

    if (0 && depth < 3) printf("depth=%d  tSep=%.1fs  total=%.1fs\n", depth,
                          (float) (clock()-start)/CLOCKS_PER_SEC,
                          (float) (clock()-startTime)/CLOCKS_PER_SEC);

#ifdef STATS_SEPARATION
    timesSep[depth] += clock()-start;
    nbSearches[depth] ++;
    sumSizes[depth] += n;
#endif

    if (bestIsUninitialized)
        exit(0);

    if (bestEval > 0)
        copySeparator(bestSep, theSeparator);
    else
        copySeparator(sep, theSeparator);

    //if (verifyTheSeparator && (! verifySeparator(S, n, theSeparator))) exit(0);

    return 1;
}




int chooseFirstVertex(int *S, int n, Graph g) {
    if (1) {
        // Test choose as first vertex one with small degree with a neighbor with max degree
        int t = rand() % nbMaxD;
        for (int *p = S;; p++)
            if (nbNeighborsInS[*p] == maxDegree)
                if (t-- == 0) {
                    int *q = g->lists[*p];
                    int vmin = *q;
                    int min = nbNeighborsInS[*q];
                    q ++;
                    while (q < g->lists[*p]+nbNeighborsInS[*p]) {
                        if (nbNeighborsInS[*q] < min) {
                            min = nbNeighborsInS[*q];
                            vmin = *q;
                        }
                        q ++;
                    }
                    return vmin;
                }
        return NONE;
    }
    else {
        int t = rand() % nbMinD;
        for (int *p = S;; p++)
            if (nbNeighborsInS[*p] == minDegree)
                if (t-- == 0) return *p;
        return NONE;
    }

    if (n <= 10) return S[rand()%n];
    int dmin, vmin;
    vmin = S[rand()%n];
    dmin = nbNeighborsInS[vmin];
    for (int i = 0; i < 5; i ++) {
        int v = S[rand()%n];
        if (nbNeighborsInS[v] < dmin) { dmin = nbNeighborsInS[v]; vmin = v; }
    }
    return vmin;
}




int separe(int *S, int n, SET V, Graph g, Separator s, int nbFlushs, int mode) {

    int firstVertex;

    nbCallsSepare ++;

    Heap A = s->A;
    Heap B = s->B;
    Heap C = s->C;

    for (int i = 0; i < nbFlushs; i ++) {

        if (useMinPartitions) {
            minheapUpdateMin(B, -1, nbNeighborsInB);
            minheapPackMins(B, nbNeighborsInB);
            if (chooseSeedNeighborMaxDegree)
                firstVertex = (i > 0) ? NONE : chooseFirstVertex(S, n, g);
            else
                firstVertex = (i > 0) ? NONE : B->val[rand()%B->nbmin];
        }
        else {
            makeHeap(B);
            firstVertex = (i > 0) ? NONE : chooseFirstVertex(S, n, g);
        }

        heapSetCompare(C, compareCVerticesFlushBtoA);
        makeHeap(C);

        flushBtoA(A, B, C, V, S, n, mode, firstVertex, g);

        if (stopSearch) break;

        if (++ i == nbFlushs) break;

        if (useMinPartitions) {
            minheapUpdateMin(A, -1, nbNeighborsInA);
            minheapPackMins(A, nbNeighborsInA);
        }
        else
            makeHeap(A);

        heapSetCompare(C, compareCVerticesFlushAtoB);
        makeHeap(C);

        flushAtoB(A, B, C, V, S, n, mode, g);

        if (stopSearch) break;
    }

    copySeparator(bestSep, s);

    //if ( ! verifySeparator(S, n, s)) exit(0);

    return 1;
}


double evalSep(int nA, int nB, int nC, int mA, int mB, int mode, int direction) {

    // Should consider here the particular cases (B->n=0, C->n=0,...)

    int min = (nA > nB) ? nB : nA;
    int max = (nA > nB) ? nA : nB;
    int delta = max-min;
    int maxE = (mA > mB) ? mA : mB;
    int minE = (mA > mB) ? mB : mA;

    if ((nA == 0) || (nB == 0)) return 0;


    if (mode == EVAL_ESSAI) {
        // minimize C size while the numbers of constraints in A and B are small and closed
        return (1.0 / (1+nC) / (1+nC) / sqrt(1.0+maxE-minE));
    }

    if (mode == EVAL_TREE_HEIGHT_COMPLETE_GRAPH)
        // approx tree height
        return (1.0 / (nC + sqrt((double) (coeffHeightNbEdges * maxE ))));


    if ((mode == EVAL_0) || (mode == EVAL_SQRT_CARD)) {
        double ra = sqrt((double) (nA));
        double rb = sqrt((double) (nB));
        double rab = sqrt((double) (nA + nB));
        double deltaSqrt = (ra > rb) ? (ra - rb) : (rb - ra);
        if (mode == EVAL_0) return ((rab - deltaSqrt)*(rab - deltaSqrt) / (1+nC) / (1+nC));
        if (nC > 0) return ((rab - deltaSqrt) / nC / nC);
        return ((rab - deltaSqrt) / (1+nC));
    }

    if (mode == EVAL_CARD) return ((double) (min) / (nC+1)); // modified 1/05

    if (mode == EVAL_LEVEL_DENSITY) {
        // evaluate the cost to place vertices of C and the smallest among A and B
        // evaluating the number of vertices of these sets per level
        if (nA < nB) return ((double) (nA+nC) / (nC + sqrt(sqrt((double) 2.0*mB))));
        else return ((double) (nB+nC) / (nC + sqrt(sqrt((double) 2.0*mA))) );
    }

    if (mode == EVAL_MINIMIZE_C) return ((double) 1 / (1+nC));

    if (mode == EVAL_MAX_SOURCE_DENSITY) {
        if (direction == FLUSH_B_A) { // try to isolate a strongly connected kernel in B
            if (nB == 0) return 0;
            if (nB <= sizeMin) return 1.0 / (nC + 1) / (nC + 1);
            return 2.0 * mB / (nB - 1) / (nB - 1) / (nC + 1) / (nC + 1);
        }
        // direction = FLUSH_A_B
        if (nA == 0) return 0;
        if (nA == 1) return 1.0 / (nC + 1);
        return 2.0 * mA / (nA - 1) / (nA - 1) / (nC + 1) / (nC + 1);
    }

    // few vertices
    if (nA+nB+nC <= 7) return ((nC > 0) ? 1.0 / (1+delta) / (1+delta) / nC : 1.0 / (1+delta));

    if (mode == EVAL_MIN_A_B) return min;

    if (mode == EVAL_MAX_C_DENSITY) return ((nC > 0) ? (double) (nbEdgesInS-mA-mB)/nC : nbEdgesInS);

    if (mode == EVAL_MAX_CDENS_MIN_AB_SIZES) return ((nC > 0) ? (double) (nbEdgesInS-mA-mB)/nC/(delta+1) : nbEdgesInS);

    if (mode == EVAL_MAX_AB_SIZES_MIN_AB_EDGES)
        return ((double) min/(1 + ((mA > mB) ? mA : mB)));

    return 0;
}



//
// Choose and remove a vertex from A, B or C (a vertex with few neighbors in B or A)
//

int chooseEltInC(Heap C) {
    if (rand() % 100 < perChooseInCAtRandom) {
        int e = C->val[rand() % (C->n)];
        return e;
    }
    return C->val[0]; // the vertex that has the fewer neighbors in the from set
}



#define CHOICE_AB_MIN_AB 0
#define CHOICE_AB_MIN_AB_MAX_C 1

int choiceModeInAandB = CHOICE_AB_MIN_AB_MAX_C;

int chooseEltInB(Heap B) {

    if (choiceModeInAandB == CHOICE_AB_MIN_AB) {
        if (useMinPartitions) return B->val[rand()%B->nbmin];
        return B->val[0];
    }

    if (choiceModeInAandB == CHOICE_AB_MIN_AB_MAX_C) {
        // Search a vertex in B with minimal number of neighbors in B and max in C.
        int best = B->val[0];
        int nbBMin = nbNeighborsInB[B->val[0]];
        int nbCMax = nbNeighborsInS[B->val[0]] - nbNeighborsInB[B->val[0]];
        int i = 1;
        while (i < B->n) {
            if (nbNeighborsInB[B->val[i]] > nbBMin) break;
            if (nbNeighborsInS[B->val[i]] - nbNeighborsInB[B->val[i]] > nbCMax) {
                best = B->val[i];
                nbCMax = nbNeighborsInS[B->val[i]] - nbNeighborsInB[B->val[i]];
            }
            i++;
        }
        return best;
    }
    return B->val[0];
}



int chooseEltInA(Heap A) {

    if (choiceModeInAandB == CHOICE_AB_MIN_AB) {
        if (useMinPartitions) return A->val[rand()%A->nbmin];
        return A->val[0];
    }

    // (choiceModeInAandB == CHOICE_AB_MIN_AB_MAX_C)
    int i = 1, best = A->val[0];
    int nbAMin = nbNeighborsInA[A->val[0]];
    int nbCMax = nbNeighborsInS[A->val[0]] - nbNeighborsInA[A->val[0]];
    while (i < A->n) {
        if (nbNeighborsInA[A->val[i]] > nbAMin) break;
        if (nbNeighborsInS[A->val[i]] - nbNeighborsInA[A->val[i]] > nbCMax) {
            best = A->val[i];
            nbCMax = nbNeighborsInS[A->val[i]] - nbNeighborsInA[A->val[i]];
        }
        i++;
    }
    return best;

    return A->val[0];
}



int testSourceConnected = 0;
int thresholdConComponents = 13;

void flushBtoA(Heap A, Heap B, Heap C, SET V, int S[], int n, int mode, int seed, Graph g) {
    int e = NONE;
    int bestASize = NONE;

    if (seed != NONE) {
        e = seed;
        goto REMOVE_FROM_B;
    }

    while (1) {
        int sizeB = B->n;

        if (stopSearch) break;

        if ((C->n > 0) && ((B->n == 0) || (rand()%100 >= perForceChoiceNotInC))) { // Should stop if C is empty ?
            e = chooseEltInC(C);
            heapRemove(e, C);
            nbEdgesInA += nbNeighborsInA[e];
        }
        else {
            if (e == NONE) e = B->val[0]; // initial case if no seed is given
            else {
                if ((useMinPartitions == 2) && (B->nbmin <= 0)) minheapSearchAndPackMins(B, nbNeighborsInB);
                e = chooseEltInB(B);
            }

        REMOVE_FROM_B:
            if (useMinPartitions == 1) minheapRemove(e, B, nbNeighborsInB);
            else if (useMinPartitions == 2) minheapJustRemove(e, B, nbNeighborsInB);
            else heapRemove(e, B);
            nbEdgesInB = nbEdgesInB-nbNeighborsInB[e];
            decreaseNbNeighborsInB(e, B, C, V, S, n, g);
        }

        // insert e in A
        heapJustAdd(e, A);
        increaseNbNeighborsInA(e, C, V, S, n, g);

        // move e neighbors which are in B to C
        removeNeighborsFromB(e, A, B, C, V, S, n, g);

        if (A->n >= sizeMax) break;

        if (B->n < sizeMin) break;

        if (A->n >= sizeMin) {
            int nB = B->n, nEB = nbEdgesInB;

            if (testSourceConnected && (B->n > thresholdConComponents) && (B->n < sizeB)) {
                int nbComp = searchConnectedComponentsInHeap(B, nbNeighborsInS, g);
                if (nbComp > 1) {
                    int imax = indexMaxInList(compSizes, nbComp);
                    nB = compSizes[imax];
                    nEB = compNbEdges[imax];
                }
            }

            double e = evalSep(A->n, nB, C->n, nbEdgesInA, nEB, mode, FLUSH_B_A);

            if ((bestIsUninitialized) || (e > bestEval)) {
                // There can be vertices in C with no neighbor in B
                if (1) pourCintoA(A, nbNeighborsInB, nbNeighborsInA, &nbEdgesInA, C, g, 1);
                bestIsUninitialized = 0;
                bestEval = e;
                if (separatorJustMemorizeReceptorSize && (g->n >= separatorJustMemoSizeMin))
                    bestASize = A->n;
                else
                    copySeparator(sep, bestSep);
            }
        }
    }
    if ((bestIsUninitialized) && (A->n < n)) {
        if (1) pourCintoA(A, nbNeighborsInB, nbNeighborsInA, &nbEdgesInA, C, g, 1);
        bestIsUninitialized = 0;
        bestEval = e;
        copySeparator(sep, bestSep);
        return;
    }
    if ((separatorJustMemorizeReceptorSize) && (bestASize != NONE)) {
        copySeparator(sep, bestSep);
        builtBestSeparator(bestSep->A, bestSep->B, bestSep->C, nbNeighborsInA, nbNeighborsInB, bestASize, g);
    }
}




void flushAtoB(Heap A, Heap B, Heap C, SET V, int S[], int n, int mode, Graph g) {
    int e = NONE;
    int bestBSize = NONE;

    while (1) {

        int sizeA = A->n;

        if (stopSearch) break;

        if ((C->n > 0) && ((A->n == 0) || (rand()%100 >= perForceChoiceNotInC)) ) {
            e = chooseEltInC(C);
            heapRemove(e, C);
            nbEdgesInB += nbNeighborsInB[e];
        }
        else {
            if (e == NONE) e = A->val[0];
            else {
                if ((useMinPartitions == 2) && (A->nbmin <= 0)) minheapSearchAndPackMins(A, nbNeighborsInA);
                e = chooseEltInA(A);
            }
            if (useMinPartitions == 1) minheapRemove(e, A, nbNeighborsInA);
            else if (useMinPartitions == 2) minheapJustRemove(e, A, nbNeighborsInA);
            else heapRemove(e, A);
            nbEdgesInA -= nbNeighborsInA[e];
            decreaseNbNeighborsInA(e, A, C, V, S, n, g);
        }

        // insert e in B
        heapJustAdd(e, B);
        increaseNbNeighborsInB(e, C, V, S, n, g);

        // move e neighbors which are in A to C
        removeNeighborsFromA(e, A, B, C, V, S, n, g);

        if (B->n >= sizeMax) break;

        if (A->n < sizeMin) break; // 25/05

        if (B->n >= sizeMin) {
            int nA = A->n, nEA = nbEdgesInA;

            if (testSourceConnected && (A->n > thresholdConComponents) && (A->n < sizeA)) {
                int nbComp = searchConnectedComponentsInHeap(A, nbNeighborsInS, g);
                if (nbComp > 1) {
                    int imax = indexMaxInList(compSizes, nbComp);
                    nA = compSizes[imax];
                    nEA = compNbEdges[imax];
                }
            }

            double e = evalSep(nA, B->n, C->n, nEA, nbEdgesInB, mode, FLUSH_A_B);

            if ((bestIsUninitialized) || (e > bestEval)) {
                if (1) pourCintoA(B, nbNeighborsInA, nbNeighborsInB, &nbEdgesInB, C, g, 1);
                bestIsUninitialized = 0;
                bestEval = e;
                if (separatorJustMemorizeReceptorSize && (g->n >= separatorJustMemoSizeMin)) {
                    bestBSize = B->n;
                }
                else
                    copySeparator(sep, bestSep);
            }
        }
    }
    if ((separatorJustMemorizeReceptorSize) && (bestBSize != NONE)) {
        int nB = B->n;
        copySeparator(sep, bestSep);
        builtBestSeparator(bestSep->B, bestSep->A, bestSep->C, nbNeighborsInB, nbNeighborsInA, bestBSize, g);
    }
}



void builtBestSeparator(Heap dest, Heap src, Heap C, int nbNdest[], int nbNsrc[], int size, Graph g) {
    // get nb neighbors (must not be modified for next calls to flush())
    memcpy(nbNInDestCopy, nbNdest, g->n * sizeof(int));

    // put back last vertices from dest to C
    for (int i = dest->n-1; i >= size; i --) {
        int e = dest->val[i];
        dest->ind[e] = NONE;
        heapJustAdd(e, C);
        for (int *pp = g->lists[e]; *pp != NONE; pp ++)
            nbNInDestCopy[*pp] --;
    }

    // Memorize the number of neighbors in dest and src for C vertices (used for the single branch ordering)
    for (int *p = C->val; p < C->val+C->n; p ++)
        nbNInABCopy[*p] = nbNInDestCopy[*p]+nbNsrc[*p];

    dest->n = size;

    // move vertices of C which have no neighbors in dest, in src
    int i = 0;
    while (i < C->n) {
        int e = C->val[i];
        if (nbNInDestCopy[e] == 0) {
            C->n --;
            C->val[i] = C->val[C->n];
            C->ind[C->val[i]] = i;
            C->ind[e] = NONE;
            heapJustAdd(e, src);
            for (int *pp = g->lists[e]; *pp != NONE; pp ++)
                nbNInABCopy[*pp] ++;
        }
        else i ++;
    }

}





// FlushBtoA() :: e has been moved to A. e neighbors which are in B must be moved to C.
void removeNeighborsFromB(int e, Heap A, Heap B, Heap C, SET V, int *S, int n, Graph g) {
    int nbToThrow = 0;

    for (int *p = g->lists[e]; p-g->lists[e] < nbNeighborsInS[e]; p ++) { //(*p != NONE) : useless, first neighbors are in S
        if (B->ind[*p] != NONE) {

            if (nbNeighborsInB[*p] == 0) verticesToThrow[nbToThrow ++] = *p;

            if (useMinPartitions == 1)
                minheapRemove(*p, B, nbNeighborsInB);
            else if (useMinPartitions == 2)
                minheapJustRemove(*p, B, nbNeighborsInB);
            else
                heapRemove(*p, B);

            heapInsert(*p, C);
            nbEdgesInB -= nbNeighborsInB[*p];

            decreaseNbNeighborsInB(*p, B, C, V, S, n, g);
        }
    }

    if ((nbToThrow == 0) || (C->n == 0) || (B->n <= 1)) return;

    for (int *p = verticesToThrow; p < verticesToThrow+nbToThrow; p ++)
    {
        if ((C->n == 0) || (B->n <= 1)) return;
        if ((C->ind[*p] != NONE) && (nbNeighborsInB[*p] == 0)) {
            heapRemove(*p, C);
            heapJustAdd(*p, A);
            nbEdgesInA += nbNeighborsInA[*p];
            increaseNbNeighborsInA(*p, C, V, S, n, g);
        }
    }
}



void removeNeighborsFromA(int e, Heap A, Heap B, Heap C, SET V, int *S, int n, Graph g) {
    int movesToB = 0;
    for (int *p = g->lists[e]; p-g->lists[e] < nbNeighborsInS[e]; p ++) {
        if (A->ind[*p] != NONE) {

            if (nbNeighborsInA[*p] == 0) movesToB ++;

            if (useMinPartitions == 1)
                minheapRemove(*p, A, nbNeighborsInA);
            else if (useMinPartitions == 2)
                minheapRemove(*p, A, nbNeighborsInA);
            else
                heapRemove(*p, A);

            heapInsert(*p, C);
            nbEdgesInA -= nbNeighborsInA[*p];

            decreaseNbNeighborsInA(*p, A, C, V, S, n, g);
        }
    }

    if ((movesToB == 0) || (C->n == 0) || (A->n <= 1)) return;

    for (int *p = g->lists[e]; p-g->lists[e] < nbNeighborsInS[e]; p ++) {
        if ((C->n == 0) || (A->n <= 1)) return;
        if ((C->ind[*p] != NONE) && (nbNeighborsInA[*p] == 0)) {
            heapRemove(*p, C);
            heapJustAdd(*p, B);
            nbEdgesInB += nbNeighborsInB[*p];
            increaseNbNeighborsInB(*p, C, V, S, n, g);
        }
    }
}



int isSmallerNeighborsInB(int a, int b) {
    if (nbNeighborsInB[a] < nbNeighborsInB[b])
        return 1;
    if (nbNeighborsInB[a] > nbNeighborsInB[b])
        return 0;
    return (priorities[a] > priorities[b]);
}


int isSmallerNeighborsInA(int a, int b) {
    if (nbNeighborsInA[a] < nbNeighborsInA[b])
        return 1;
    if (nbNeighborsInA[a] > nbNeighborsInA[b])
        return 0;
    return (priorities[a] > priorities[b]);
}


int modeChoiceInCflushBA = 0;
int modeChoiceInCflushAB = 0;

// Comparator for C vertices: fewer neighbors in B is better (more neighbors in A also)
// Returns 1 if a is smaller (better) than b
int compareCVerticesFlushBtoA(int a, int b) {

    if (nbNeighborsInB[a] == 0) {
        if ((nbNeighborsInB[b] == 0) && (nbNeighborsInA[a] < nbNeighborsInA[b])) return 0;
        return 1;
    }
    if (nbNeighborsInB[b] == 0) return 0;

    int nbNaB = nbNeighborsInB[a];
    int nbNaA = nbNeighborsInA[a];
    int nbNbB = nbNeighborsInB[b];
    int nbNbA = nbNeighborsInA[b];

    if (modeChoiceInCflushBA == 0) {
        if (nbNeighborsInB[a] < nbNeighborsInB[b])
            return 1;
        if (nbNeighborsInB[a] > nbNeighborsInB[b])
            return 0;
        if (nbNeighborsInA[a] > nbNeighborsInA[b])
            return 1;
        if (nbNeighborsInA[a] < nbNeighborsInA[b])
            return 0;
        return (priorities[a] > priorities[b]);
    }

    if (modeChoiceInCflushBA == 1) {
        // If the difference bteween the numbers of neighbors in B is small but one has many more neighbors in A, then prefer this one
        if (nbNaB < nbNbB) {
            if (nbNaA >= nbNbA) return 1;
            return  (nbNbA - nbNaA < (nbNbB - nbNaB + 2) * (nbNbB - nbNaB + 2) - 1);
        }
        if (nbNaB > nbNbB) {
            if (nbNaA <= nbNbA) return 0;
            return (nbNaA - nbNbA >= (nbNaB - nbNbB + 2) * (nbNaB - nbNbB + 2) - 1);
        }
        if (nbNaA == nbNbA) return (priorities[a] > priorities[b]);
        return (nbNaA > nbNbA);
    }


    if (modeChoiceInCflushBA == 2) {
        if (nbNaB < nbNbB) {
            if (nbNaA >= nbNbA) return 1;
            if ((nbNaA == 0) || (100 * nbNbA) / nbNaA >= 100 * (1 + nbNbB - nbNaB)) return 0;
            return 1;
        }
        if (nbNaB > nbNbB) {
            if (nbNaA <= nbNbA) return 0;
            if ((nbNbA == 0) || (100 * nbNaA) / nbNbA >= 100 * (1 + nbNaB - nbNbB)) return 1;
            return 0;
        }
        if (nbNaA == nbNbA) return (priorities[a] > priorities[b]);
        return (nbNaA > nbNbA);
    }
    return 1;

    if (0) {
        float ea = nbNeighborsInB[a] * nbNeighborsInB[a] / nbNeighborsInA[a];
        float eb = nbNeighborsInB[b] * nbNeighborsInB[b] / nbNeighborsInA[b];

        if (ea < eb) return 1;
        if (ea > eb) return 0;
        return (priorities[a] > priorities[b]);
    }
}



int compareCVerticesFlushAtoB(int a, int b) {

    if (nbNeighborsInA[a] == 0) {
        if ((nbNeighborsInA[b] == 0) && (nbNeighborsInB[a] < nbNeighborsInB[b])) return 0;
        return 1;
    }
    if (nbNeighborsInA[b] == 0) return 0;


    if (modeChoiceInCflushAB == 0) {
        if (nbNeighborsInA[a] < nbNeighborsInA[b])
            return 1;
        if (nbNeighborsInA[a] > nbNeighborsInA[b])
            return 0;
        if (nbNeighborsInB[a] > nbNeighborsInB[b])
            return 1;
        if (nbNeighborsInB[a] < nbNeighborsInB[b])
            return 0;
        return (priorities[a] > priorities[b]);
    }


    if (modeChoiceInCflushAB == 1) {
        if (nbNeighborsInA[a] < nbNeighborsInA[b]) {
            if (nbNeighborsInB[a] >= nbNeighborsInB[b]) return 1;
            if (nbNeighborsInB[b] - nbNeighborsInB[a] >=
                (nbNeighborsInA[b] - nbNeighborsInA[a] + 1) * (nbNeighborsInA[b] - nbNeighborsInA[a] + 1) )
                return 0;
            return 1;
        }
        if (nbNeighborsInA[a] > nbNeighborsInA[b]) {
            if (nbNeighborsInB[a] <= nbNeighborsInB[b]) return 0;
            if (nbNeighborsInB[a] - nbNeighborsInB[b] >=
                (nbNeighborsInA[a] - nbNeighborsInA[b] + 1) * (nbNeighborsInA[a] - nbNeighborsInA[b] + 1) )
                return 1;
            return 0;
        }
        //if (nbNeighborsInB[a] == nbNeighborsInB[b]) {
        if (nbNeighborsInB[a] > nbNeighborsInB[b]) return 1;
        if (nbNeighborsInB[a] < nbNeighborsInB[b]) return 0;
        return (priorities[a] > priorities[b]);
    }

    if (modeChoiceInCflushAB == 2) {
        if (nbNeighborsInA[a] < nbNeighborsInA[b]) {
            if (nbNeighborsInB[a] >= nbNeighborsInB[b]) return 1;
            if ((nbNeighborsInB[a] == 0) || ((100*nbNeighborsInB[b]) / nbNeighborsInB[a] >=
                                             100 * (1 + nbNeighborsInA[b] - nbNeighborsInA[a]) ))
                return 0;
            return 1;
        }
        if (nbNeighborsInA[a] > nbNeighborsInA[b]) {
            if (nbNeighborsInB[a] <= nbNeighborsInB[b]) return 0;
            if ((nbNeighborsInB[b] == 0) || ((100*nbNeighborsInB[a]) / nbNeighborsInB[b] >=
                                             100 * (1 + nbNeighborsInA[a] - nbNeighborsInA[b]) ))
                return 1;
            return 0;
        }
        //if (nbNeighborsInB[a] == nbNeighborsInB[b]) {
        if (nbNeighborsInB[a] > nbNeighborsInB[b]) return 1;
        if (nbNeighborsInB[a] < nbNeighborsInB[b]) return 0;
        return (priorities[a] > priorities[b]);
    }
    return 1;
}



//
// Dust separator: for vertices of C that have no neighbor in A or B
//

void pourCintoA(Heap A, int nbNFrom[], int nbNTo[], int *nbEdges, Heap C, Graph g, int justAdd) {
    for (int i = 0; i < C->n; i++) {
        int e = C->val[i];
        if (nbNFrom[e] == 0) {
            heapRemove(e, C);
            if (justAdd)
                heapJustAdd(e, A);
            else
                heapInsert(e, A);
            *nbEdges = *nbEdges + nbNTo[e];
            for (int *p = g->lists[e]; *p != NONE; p ++) nbNTo[*p]++;
        }
    }
}


//
// C improvement
//

int nbCVerticesWithNoNeighborInAAndB(Heap C, SET V, int *S, int n, Graph g) {
    if (C->n == 0) return 0;
    int nb = 0;
    for (int *p = C->val; p < C->val+C->n; p ++) {
        if ((nbNeighborsInB[*p] == 0) && (nbNeighborsInA[*p] == 0)) nb ++;
    }
    return nb;
}

// Search vertices of C with no neighbor in A nor in B, count the number of neighbors in C,
// finally order them (those with the fewer neighbors in C first) and make a selection
// of independent such vertices. Each of them will form an isolated vertex.
void selectABDisconnectedVertices(Separator sep, SET V, int *S, int n, int date, Graph g) {
    Heap C = sep->C;
    if (C->n == 0) { sep->nbABDV = 0; return ; }
    int nb = 0;
    for (int *p = C->val; p < C->val+C->n; p ++) {
        if ((nbNeighborsInB[*p] == 0) && (nbNeighborsInA[*p] == 0)) {
            dateABDisconnected[*p] = date; // mark *p as an AB Disconnected certex for current separe call
            nbNeighborsABDisconnected[*p] = 0;
            swapValues(C->val+nb, p);
            nb ++;
        }
    }
    if (nb == 0) { sep->nbABDV = 0; return;}

    printf("nb AB Disconnected = %d :: ", nb);
    for (int i = 0; i < nb; i ++) printf("%d (%d,%d)   ", C->val[i], nbNeighborsInB[C->val[i]], nbNeighborsInA[C->val[i]]);
    printf("\n");
    printList(sep->A->val, sep->A->n);
    printList(sep->B->val, sep->B->n);
    printList(C->val, C->n);


    if (nb <= 1) { sep->nbABDV = nb; return ; }

    sep->nbABDV = 1;
    return ;

    // Calculate the degree of each AB Disconnected vertex in the set of AB Disconnected vertices.
    for (int *p = C->val; p < C->val+nb; p ++) {
        for (int *q = g->lists[*p]; q < g->lists[*p]+nbNeighborsInS[*p]; q ++) {
            if (dateABDisconnected[*q] == date)
                nbNeighborsABDisconnected[*q] ++;
        }
    }
    // Order AB Disconnected vertices, the less connected first.
    for (int i = nb; i > 1; i --) {
        int someSwap = 0;
        for (int j = nb-1; j > 0; j --) {
            if (nbNeighborsABDisconnected[C->val[j]] < nbNeighborsABDisconnected[C->val[j - 1]]) {
                swap(C->val, j, j - 1);
                someSwap = 1;
            }
        }
        if (someSwap == 0) break;
    }
    // Select independent AB Disconnected vertices (mark dateABDisconnected[*p] = -nbCallsSepare when selected)
    dateABDisconnected[C->val[0]] = -date;
    int nbS = 1;
    int i = 1;
    while (i < nb) {
        // the ith AB Disconnected vertex is chosen if independent with previous selected vertices
        int v = C->val[i];
        int *q = g->lists[v];
        for ( ; q < g->lists[v]+nbNeighborsInS[v]; q ++)
            if (dateABDisconnected[*q] == -date)
                break;
        if (q == g->lists[v]+nbNeighborsInS[v]) {
            swap(C->val, i, nbS);
            dateABDisconnected[v] = -date;
            nbS ++;
            i ++;
        }
    }
    printf("nbS=%d\n", nbS);

    // verif
    if (0) {
        for (int i = 0; i < nbS - 1; i++)
            for (int j = i + 1; j < nbS; j++)
                if (areNeighbours(C->val[i], C->val[j], g))
                    printf("OUYE\n");
    }
    sep->nbABDV = nbS;
}










//
// Improvement: search moves from A/B to C that generate moves from C to A/B
// Seems useless.


void improveSeparation(Separator s, Graph g, int nSteps) {
    int moves[s->C->n];
    int nbmoves;

    START:
    nbmoves = 0;
    // Moves A->C (neighbors C->B)
    for (int i = 0; i < nSteps; i ++) {
        int u;
        int nbMoves = searchMoveAC(s, g, moves, &u);
        if (u == NONE)
            return;
        makeMove(u, s->A, s->C);
        for (int j = 0; j < nbMoves; j ++) {
            makeMove(moves[j], s->C, s->B);
        }
        nbMoves ++;
    }

    // Moves B->C (neighbors C->A)
    for (int i = 0; i < nSteps; i ++) {
        int u;
        int nbMoves = searchMoveBC(s, g, moves, &u);
        if (u == NONE)