Gerrymander North Carolina The challenge

How well can you gerrymander North Carolina into 13 voting districts?

In this challenge, you use the following files to draw different maps for Republicans and Democrats.

File 1: NC_2016_presidential_election.csv

File 2: NC_county_adjacency.csv

File 1 gives county-level voting data from the 2016 presidential election (based on this and that), while File 2 gives adjacency information between counties (based on this). For File 2, the first entry of each row is adjacent to the other entries in that row. For example, Alexander is adjacent to Alexander, Caldwell, Catawba, Iredell and Wilkes. (Note that every county is adjacent to itself.)

Admissibility

Your code should be written in a way that can gerrymander arbitrary maps. In particular, your code must take each of the following as inputs:


A set \$C\$ of counties and a set \$P\$ of parties.
Weights \$(v_{cp})\$ that count voters in county \$c\in C\$ who vote for party \$p\in P\$. (E.g., File 1.)
A connected graph \$G\$ with vertex set \$C\$. (E.g., File 2.)
A number of districts \$k\$.
A party \$p^\star\in P\$ to gerrymander in favor of.


The output must be an assignment of counties to districts such that the counties in each district induce a connected subgraph of \$G\$. For example, for the counties in Files 1 and 2 (listed in alphabetical order), the following is one of 13! sequences that specify the partition of counties into districts pictured below:

3,5,5,10,5,5,13,11,10,8,9,12,1,12,13,13,3,9,3,12,13,12,9,8,13,10,13,13,4,4,8,11,11,5,7,9,11,12,11,13,4,11,7,12,12,11,1,13,4,12,7,13,3,13,9,12,12,11,12,2,12,1,1,7,8,11,13,3,13,13,8,13,3,11,9,3,10,10,3,1,9,8,10,1,5,5,12,12,13,10,11,6,11,11,5,8,5,7,5,12




(This map is loosely based on the current map of U.S. congressional districts in North Carolina.)

Scoring

Your score is determined by how your code performs for various problem instances.

Given a district assignment, let \$d\$ denote the number of districts for which party \$p^\star\$ receives the plurality of votes cast, and let \$a\$ and \$b\$ denote the minimum and maximum number of voters in a district, respectively. Then the score for the district assignment is given by

$$ \text{assignment score} = d - \frac{b}{a}. $$

For example, if \$p^\star\$ is the Republican party, then the district assignment described above has \$d=10\$, \$a=247739\$, and \$b=527337\$, leading to a score of about \$7.87\$.

To compute your score, gerrymander North Carolina (as specified in Files 1 and 2) into \$k=13\$ districts in favor of the Republican party, and then in favor of the Democratic party. Let \$S_1\$ denote the sum of these two assignment scores. Next, gerrymander 10 different modifications of North Carolina, each obtained by removing one of the following counties from Files 1 and 2:

Cleveland,Davidson,Hyde,Johnston,Madison,Montgomery,Pender,Scotland,Warren,Wilkes


(Since North Carolina is biconnected, each of these are acceptable problem instances.) For each modification, gerrymander into \$k=13\$ districts in favor of the Republican party, and then in favor of the Democratic party. Let \$S_2\$ denote the sum of these 20 assignment scores. Then the score of your submission is given by

$$ \text{submission score} = S_1+\frac{S_2}{10}. $$

In addition to your code, you are encouraged to provide an illustration of your two maps that gerrymander all of North Carolina for Republicans and for Democrats. (DRA2020 might be a helpful resource.)
C++, score = 33.526203, all test cases in 2min34s
(of course, the score can be improved by running the program for more time; please do not do that to claim a victory or something)
That score is usually 11 won districts for the Republicans and from 9 to 11 for the Democrats.
In this revision, the score decreased significantly, because I found some really large bugs in the input parsing. My program also verifies its solution now.
//#define _GLIBCXX_DEBUG
#include &lt;iostream&gt;
#include &lt;streambuf&gt;
#include &lt;iomanip&gt;
#include &lt;fstream&gt;
#include &lt;cstring&gt;
#include &lt;bitset&gt;
#include &lt;cassert&gt;
#include &lt;cstdio&gt;
#include &lt;vector&gt;
#include &lt;algorithm&gt;
#include &lt;cmath&gt;
#include &lt;climits&gt;
#include &lt;random&gt;
#include &lt;set&gt;
#include &lt;map&gt;
#include &lt;deque&gt;
#include &lt;string&gt;

constexpr uint64_t rotl(uint64_t x, char k)
{
    return (x&lt;&lt;k) | (x&gt;&gt;(64-k)); //typo of the day: replace the &gt;&gt; with &lt;&lt;
}
uint64_t xs128()
{
    static uint64_t s_0 = 0, s_1 = 1;
    const uint64_t s0 = s_0;
    uint64_t s1 = s_1;
    uint64_t res = s0 + s1;
    s1 ^= s0;
    s_0 = rotl(s0, 24) ^ s1 ^ (s1 &lt;&lt; 16);
    s_1 = rotl(s1, 37);
    return res &gt;&gt; 4;
}
std::vector&lt;std::string&gt; getCsvLine(std::istream&amp; str)
{
    https://stackoverflow.com/a/1120224 - saved some time writing a parser...
    std::vector&lt;std::string&gt; result;
    std::string line; std::getline(str,line);
    std::stringstream lineStream(line);
    std::string cell;
    while(std::getline(lineStream,cell, ','))
        result.push_back(cell);
    return result;
}
struct county
{
    int good = 0;
    int sum = 0;
    std::vector&lt;int&gt; others;
    std::string name;
};
struct district
{
    int good = 0;
    int sum = 0;
    std::vector&lt;int&gt; others;
    void operator += (county&amp; rhs)
    {
        good += rhs.good; sum += rhs.sum;
        for(int i = 0; i &lt; others.size(); i++)
            others[i] += rhs.others[i];
    }
    void operator -= (county&amp; rhs)
    {
        good -= rhs.good; sum -= rhs.sum;
        for(int i = 0; i &lt; others.size(); i++)
            others[i] -= rhs.others[i];
    }
};
//seems like exactly 100 districts exist, below CHAR_MAX
std::vector&lt;std::vector&lt;char&gt;&gt; graph;
std::vector&lt;char&gt; curdistrict;
std::vector&lt;county&gt; counties;
std::vector&lt;district&gt; districts;
std::vector&lt;char&gt; marks;
bool pluralityCheck(const district&amp; d)
{
    return d.good &gt; *std::max_element(d.others.begin(), d.others.end());
}
int n, k, excluded;
double getFitness(std::vector&lt;char&gt;&amp; dmap)
{
    std::vector&lt;district&gt; districts(k);
    for(district&amp; d : districts)
        d.others.resize(counties[excluded==0?1:0].others.size());
    for(int i = 0; i &lt; n; i++)
        if(i != excluded) districts[dmap[i]] += counties[i];
    double score = 0;
    int min = INT_MAX, max = INT_MIN;
    for(int i = 0; i &lt; k; i++)
        score += pluralityCheck(districts[i]),
        min = std::min(min, districts[i].sum), max = std::max(max, districts[i].sum);
    score -= double(max) / min;
    return score;
}
void dfs(char v, char color)
{
    marks[v] = color;
    for(int to : graph[v])
        if(marks[to] != color &amp;&amp; curdistrict[to] == curdistrict[v])
            dfs(to, color);
}
bool dfsFor(char at, char target, char color)
{
    marks[at] = color;
    if(at == target) return true;
    if(at == excluded) return false;
    bool ans = false;
    for(int to : graph[at])
        if(marks[to] != color &amp;&amp; curdistrict[to] == curdistrict[at])
            ans |= dfsFor(to, target, color);
    return ans;
}
bool isCutpoint(int v)
{
    static int call = 0;
    char color = call + 1;
    call = (call + 1) % 255;
    if(color == 1) for(char&amp; el : marks) el = 0; //reset every 256 calls
    marks[v] = color;
    int calls = 0;
    for(char to : graph[v])
    {
        if(curdistrict[to] != curdistrict[v]) continue;
        if(marks[to] == color) continue;
        calls++;
        if(calls == 2) break;
        dfs(to, color);
    }
    return calls != 1; //0 -&gt; kills a district, &gt;=2 -&gt; a cutpoint -&gt; splits a district
}
std::map&lt;std::string, int&gt; dnamemap;
double solve(int n, int k, int target, int excluded)
{
    ::n = n; ::k = k; ::excluded = excluded;
    marks = std::vector&lt;char&gt;(n);
    curdistrict = std::vector&lt;char&gt;(n, -1); graph = std::vector&lt;std::vector&lt;char&gt;&gt;(n);
    districts = std::vector&lt;district&gt;(k); counties = std::vector&lt;county&gt;(n);
    dnamemap = std::map&lt;std::string, int&gt;();
    std::ifstream data(&quot;NC_2016_presidential_election.csv&quot;);
    getCsvLine(data); //skip header
    for(int i = 0; i &lt; n; i++)
    {
        std::vector&lt;std::string&gt; row = getCsvLine(data);
        dnamemap[row[0]] = i;
        counties[i].name = row[0];
        if(i == excluded) continue;
        for(int j = 1; j &lt; row.size(); j++)
        {
            int v = std::stoi(row[j]);
            if(j - 1 == target)
                counties[i].good = v;
            else counties[i].others.push_back(v);
            counties[i].sum += v;
        }

    }
    data = std::ifstream(&quot;NC_county2.csv&quot;); //no header
    //NC_county2.csv is obtained by cat NC_county_adjacency.csv | sort | tr -d '\r'
    //the default counties file caused not one but TWO horrible parsing bugs
    //they silently increased the score.
    //There were \r\n line endings instead of \n, and my program mishandled them.
    //and *the worst*. It seems to be in alphabetic order, but it actually simply isn't.
    for(int i = 0; i &lt; n; i++)
    {
        std::vector&lt;std::string&gt; row = getCsvLine(data);
        if(i == excluded) continue;
        for(int j = 1; j &lt; row.size(); j++)
        {
            if(row[j] == row[0]) continue;
            int to = dnamemap[row[j]];
            if(to == excluded) continue;
            graph[i].push_back(to);
        }
    }
    for(district&amp; d : districts)
        d.others.resize(counties[excluded==0?1:0].others.size());
    int di = 0;
    for(int i = 0; di &lt; k; i++) if(i != excluded) curdistrict[i] = di++;
    for(int it = n; it --&gt; 0;)
        for(int i = 0; i &lt; n; i++)
        {
            if(i == excluded || curdistrict[i] != -1) continue;
            for(int to : graph[i])
            {
                if(to == excluded || curdistrict[to] == -1) continue;
                curdistrict[i] = curdistrict[to];
                break;
            }
        }
    if(excluded &gt;= 0) curdistrict[excluded] = -1;
    for(int i = 0; i &lt; n; i++) printf(&quot;%d &quot;, curdistrict[i]);
    printf(&quot;%d\n&quot;, excluded);
    for(int i = 0; i &lt; n; i++) if(excluded != i) districts[curdistrict[i]] += counties[i];
    double temp = 3;
    int its = 2e7;
    std::vector&lt;char&gt; bestplacement = curdistrict;
    std::vector&lt;district&gt; bestdistricts = districts;
    int lastimprovement = 0;
    double bestscore = getFitness(curdistrict);
    int perblock = 2e4; //max number of iterations without improvement
    double delta = 1.01 * temp / its;
    double curscore = bestscore;
    double withoutimprovement = 0;
    for(int it = 0; it &lt; its; it++)
    {
        temp -= delta;
        if(it % 1048576 == 0)
            printf(&quot;it=%d, cs=%lf       \r&quot;, it, curscore);
        int i = xs128() % n;
        while(i == excluded) i = xs128() % n;
        if(isCutpoint(i)) continue; //&lt;- around 50% of program time
        std::vector&lt;int&gt; nearby;
        bool ok = false;
        for(int to : graph[i])
            if(curdistrict[to] != curdistrict[i] &amp;&amp; curdistrict[to] != -1) ok = true;
        if(!ok) continue;
        int nd = -1;
        while(nd == curdistrict[i] || nd == -1)
            nd = curdistrict[graph[i][xs128() % graph[i].size()]];
        //calculate new min size, new max size etc.
        //can be done via std::set and other tree-based stuff, but it's much simpler to use a normal arrays (for k=13)
        //*might* use a segment tree to optimize later, but it's not spending a lot of time here
        char oldmaj1 = pluralityCheck(districts[curdistrict[i]]);
        char newmaj1 = pluralityCheck(districts[nd]);
        int maxs = INT_MIN, mins = INT_MAX;
        for(int i = 0; i &lt; k; i++)
            maxs = std::max(maxs, districts[i].sum), mins = std::min(mins, districts[i].sum);
        double pen1 = double(maxs) / mins;
        districts[curdistrict[i]] -= counties[i];
        districts[nd] += counties[i];
        char oldmaj2 = pluralityCheck(districts[curdistrict[i]]);
        char newmaj2 = pluralityCheck(districts[nd]);
        maxs = INT_MIN, mins = INT_MAX;
        for(int i = 0; i &lt; k; i++)
            maxs = std::max(maxs, districts[i].sum), mins = std::min(mins, districts[i].sum);
        double pen2 = double(maxs) / mins;
        double delta = pen1 - pen2 + (oldmaj2 + newmaj2 - oldmaj1 - newmaj1);
        if(delta &lt;= 0 &amp;&amp; (temp &lt;= 0 || ldexpf(std::exp(delta / temp), 60) &lt; xs128()))
        {
            districts[nd] -= counties[i], districts[curdistrict[i]] += counties[i];
        }
        else
        {
            curdistrict[i] = nd;
            curscore += delta;
            if(curscore &gt; bestscore)
            {
                bestscore = curscore; lastimprovement = it;
                bestplacement = curdistrict; withoutimprovement = 0;
                bestdistricts = districts;
            }
        }
        if(curscore &lt; bestscore) withoutimprovement += bestscore - curscore;
        if(it - lastimprovement &gt; perblock || withoutimprovement &gt; 3e4)
        {
            //restart at best score so far
            withoutimprovement = 0;
            lastimprovement = it;
            curdistrict = bestplacement;
            districts = bestdistricts;
            curscore = bestscore;
        }
    }
    for(int i = 0; i &lt; n; i++)
        printf(&quot;%d &quot;, curdistrict[i]);
    printf(&quot;\n&quot;);
    for(district&amp; d : districts)
    {
        printf(&quot;%d total, maj=%d, %d good; rest:&quot;, d.sum, (int)pluralityCheck(d), d.good);
        for(int el : d.others) printf(&quot; %d&quot;, el);
        printf(&quot;\n&quot;);
    }
    for(int i = 0; i &lt; n; i++) //verify validness
    for(int j = i+1; j &lt; n; j++)
    {
        if(i == excluded || j == excluded) continue;
        if(curdistrict[i] != curdistrict[j]) continue;
        for(char&amp; el : marks) el = 0;
        bool good = dfsFor(i, j, 1);
        if(!good) printf(&quot;invalid: i=%d, j=%d, d=%d!\n&quot;, i, j, curdistrict[i]);
    }
    double score = getFitness(curdistrict);
    printf(&quot;%lf\n&quot;, score);
    return score;
}
int main()//int64_t argc, char*argv[])
{
    //randomly choose and recolor vertices
    //do not recolor a vertex if it is a cutpoint (more known as an articulation point or cut vertex) or alone
    //to check, DFS over same-colored vertices from it. If != 1 calls were done, it's alone or a cutpoint
    //initial state: all counties belong to district 0, but first 13 belong to district i (dumb but why bother)
    //for each district maintain the voters for both parties and total number
    //oh, and don't merely randomly recolor: use Simulated Annealing (TM) (it actually improves a lot!)
    setbuf(stdout, 0);
    int k = 13, n = 100;
    std::vector&lt;int&gt; toExclude { -1, 22, 28, 47, 50, 57, 61, 70, 82, 92, 97 }; //the final score computing
    //std::vector&lt;int&gt; toExclude { -1 };
    double score = 0;
    for(int el : toExclude)
    {
        double s = 0;
        s += solve(n, k, 0, el);
        s += solve(n, k, 1, el);
        if(el == -1) score += s;
        else score += s / 10;
    }
    printf(&quot;%lf\n&quot;, score);
}

Uses (attempts to use) simulated annealing to improve from a dumb assignment by making random changes if they do not ruin things (that is, destroy or split districts). Restarts from the best known assignment when either 20000 iterations passed without an improvement or the accumulated difference between the best score and the current score since the last improvement exceeds 30000.
Can be compiled with clang++ hax.cpp -Ofast -march=native -flto -no-pie -o a.out. Assumes files NC_2016_presidential_election.csv and NC_county2.csv are in the current directory, where NC_county2.csv is the NC_county_adjacency.csv, but with its lines actually sorted in alphabetical order, and with Unix (\n) line endings.
Quoting from the question, &quot;you are encouraged to provide an illustration&quot;. I am definitely encouraged (I would really like to know whether or not the code has no bugs), but I have no idea how to create an illustration. I also have no idea how to export the output into a format readable by the linked website.
Assignment information on the normal map for the Democrats, then for the Republicans:


begin assignment for target=0, excluded=-1
Alamance,4
Alexander,8
Alleghany,8
Anson,12
Ashe,8
Avery,5
Beaufort,6
Bertie,6
Bladen,9
Brunswick,9
Buncombe,5
Burke,10
Cabarrus,8
Caldwell,8
Camden,7
Carteret,9
Caswell,0
Catawba,10
Chatham,12
Cherokee,10
Chowan,7
Clay,10
Cleveland,10
Columbus,2
Craven,9
Cumberland,2
Currituck,7
Dare,7
Davidson,12
Davie,8
Duplin,6
Durham,12
Edgecombe,6
Forsyth,3
Franklin,7
Gaston,10
Gates,7
Graham,10
Granville,0
Greene,6
Guilford,4
Halifax,7
Harnett,12
Haywood,10
Henderson,10
Hertford,7
Hoke,2
Hyde,7
Iredell,8
Jackson,10
Johnston,12
Jones,6
Lee,12
Lenoir,6
Lincoln,11
Macon,10
Madison,10
Martin,6
McDowell,5
Mecklenburg,11
Mitchell,5
Montgomery,12
Moore,12
Nash,7
New Hanover,9
Northampton,7
Onslow,9
Orange,0
Pamlico,9
Pasquotank,7
Pender,9
Perquimans,7
Person,0
Pitt,6
Polk,10
Randolph,4
Richmond,12
Robeson,2
Rockingham,0
Rowan,8
Rutherford,10
Sampson,9
Scotland,2
Stanly,8
Stokes,0
Surry,3
Swain,10
Transylvania,10
Tyrrell,7
Union,8
Vance,7
Wake,1
Warren,7
Washington,6
Watauga,5
Wayne,9
Wilkes,8
Wilson,6
Yadkin,8
Yancey,10
end assignment for target=0, excluded=-1
begin assignment for target=1, excluded=-1
Alamance,0
Alexander,1
Alleghany,11
Anson,7
Ashe,11
Avery,10
Beaufort,6
Bertie,6
Bladen,8
Brunswick,12
Buncombe,10
Burke,10
Cabarrus,7
Caldwell,11
Camden,9
Carteret,6
Caswell,0
Catawba,1
Chatham,2
Cherokee,4
Chowan,9
Clay,4
Cleveland,4
Columbus,12
Craven,9
Cumberland,8
Currituck,6
Dare,6
Davidson,0
Davie,7
Duplin,9
Durham,0
Edgecombe,9
Forsyth,1
Franklin,6
Gaston,4
Gates,9
Graham,10
Granville,6
Greene,6
Guilford,11
Halifax,6
Harnett,8
Haywood,10
Henderson,4
Hertford,9
Hoke,2
Hyde,6
Iredell,7
Jackson,4
Johnston,8
Jones,9
Lee,2
Lenoir,9
Lincoln,10
Macon,4
Madison,10
Martin,9
McDowell,10
Mecklenburg,5
Mitchell,10
Montgomery,2
Moore,2
Nash,9
New Hanover,12
Northampton,9
Onslow,12
Orange,2
Pamlico,9
Pasquotank,6
Pender,12
Perquimans,6
Person,0
Pitt,6
Polk,4
Randolph,0
Richmond,7
Robeson,8
Rockingham,0
Rowan,2
Rutherford,4
Sampson,12
Scotland,2
Stanly,2
Stokes,11
Surry,11
Swain,10
Transylvania,4
Tyrrell,6
Union,7
Vance,6
Wake,3
Warren,9
Washington,6
Watauga,11
Wayne,9
Wilkes,1
Wilson,9
Yadkin,1
Yancey,10
end assignment for target=1, excluded=-1