#include "lib.h"
#include <stdlib.h>

#define DBG 0

static
ch input[] = {
#embed "day8_input.txt"
};

static
ch test[] = {
#embed "day8_test.txt"
};


typedef num point __attribute__((vector_size(4*sizeof(num))));

alignas(64) static
point positions[1000];

alignas(64) static
unsigned short circuits[1000];

alignas(1024) static
struct { num min_dist; unsigned short min_target; } targets[1000][1024];

alignas(64) static
unsigned short next_targets[1000];

alignas(64) static
num next_dists[1000];

alignas(64) static
unsigned short next_indices[1000];

static
int targetcompar(const void *a, const void *b) {
    const num ad = ((typeof(**targets)*)a)->min_dist;
    const num bd = ((typeof(**targets)*)b)->min_dist;
    if (ad > bd) return 1;
    else if (ad < bd) return -1;
    else return 0;
}

static
num do_part1(size_t file_len, ch file[file_len], bool istest) {
    const num iters = istest ? 10 : 1000;
    const num pointcount = istest ? 20 : 1000;
    num iter = 0;
#if DBG
    pnl();
#endif

    // read in positions from file
    { ch *s = file;
    for (num i = 0; i < pointcount; i++) {
        unsigned x = grabnum(s, &s);
        s++;
        unsigned y = grabnum(s, &s);
        s++;
        unsigned z = grabnum(s, &s);
        s++;
        positions[i] = (point){ x,y,z };
    } }
    // set up initial circuits
    for (num i = 0; i < pointcount; i += 4) {
        circuits[i+0] = i+0;
        circuits[i+1] = i+1;
        circuits[i+2] = i+2;
        circuits[i+3] = i+3;
    }
    // find distances between all points
    for (num i = 0; i < pointcount; i++) {
        point pi = positions[i];
        for (num j = 0; j < pointcount; j += 4) {
            point pj0 = positions[j+0];
            point pj1 = positions[j+1];
            point pj2 = positions[j+2];
            point pj3 = positions[j+3];

            point dvp0 = (pi - pj0) * (pi - pj0);
            point dvp1 = (pi - pj1) * (pi - pj1);
            point dvp2 = (pi - pj2) * (pi - pj2);
            point dvp3 = (pi - pj3) * (pi - pj3);

            num d0 = dvp0[0] + dvp0[1] + dvp0[2] + dvp0[3];
            num d1 = dvp1[0] + dvp1[1] + dvp1[2] + dvp1[3];
            num d2 = dvp2[0] + dvp2[1] + dvp2[2] + dvp2[3];
            num d3 = dvp3[0] + dvp3[1] + dvp3[2] + dvp3[3];

            targets[i][j+0] = (typeof(**targets)) { d0, j+0 };
            targets[i][j+1] = (typeof(**targets)) { d1, j+1 };
            targets[i][j+2] = (typeof(**targets)) { d2, j+2 };
            targets[i][j+3] = (typeof(**targets)) { d3, j+3 };
        }
        // sort targets. slow.
        qsort(targets[i], pointcount, sizeof(**targets), targetcompar);
    }
    // set up initial next targets
    for (num i = 0; i < pointcount; i++) {
        next_targets[i] = targets[i][1].min_target;
        next_dists[i]   = targets[i][1].min_dist;
        next_indices[i] = 1;
    }

    // FIXED JOIN ITERATIONS
    do {
        // find min dist
        num min_dist = -1;
        num min_i = 0;
        for (num i = 0; i < pointcount; i++)
            if (next_dists[i] < min_dist)
                min_dist = next_dists[min_i = i];
        num min_j = next_targets[min_i];
        // -- join the pair --
#if DBG
        point ipos = positions[min_i];
        point jpos = positions[min_j];
        printd(ipos[0]); print(","); printd(ipos[1]); print(","); printd(ipos[2]);
        print("\t->\t");
        printd(jpos[0]); print(","); printd(jpos[1]); print(","); printd(jpos[2]);
        pnl();
#endif
        // update circuits
        unsigned old_circuit = circuits[min_j];
        unsigned new_circuit = circuits[min_i];
        if (old_circuit != new_circuit) {
            for (num c = 0; c < pointcount; c++)
                if (circuits[c] == old_circuit)
                    circuits[c] = new_circuit;
        }
        // update targets
        unsigned niidx = ++next_indices[min_i];
        unsigned njidx = ++next_indices[min_j];
        next_targets[min_i] = targets[min_i][niidx].min_target;
        next_targets[min_j] = targets[min_j][njidx].min_target;
        next_dists[min_i] = targets[min_i][niidx].min_dist;
        next_dists[min_j] = targets[min_j][njidx].min_dist;

        iter++;
    } while(iter < iters);

    // reuse next_targets to store circuit counts
    for (num i = 0; i < pointcount; i++) {
        next_targets[i] = 0;
    }
    for (num i = 0; i < pointcount; i++) {
        next_targets[circuits[i]]++;
    }
    num circ0 = 0, circ1 = 0, circ2 = 0;
    for (num i = 0; i < pointcount; i++) {
        if (next_targets[i] > circ0) {
            circ2 = circ1;
            circ1 = circ0;
            circ0 = next_targets[i];
        } else if (next_targets[i] > circ1) {
            circ2 = circ1;
            circ1 = next_targets[i];
        } else if (next_targets[i] > circ2) {
            circ2 = next_targets[i];
        }
    }

    return circ0 * circ1 * circ2;
}

static
num do_part2(size_t file_len, ch file[file_len], bool istest) {
    const num pointcount = istest ? 20 : 1000;
#if DBG
    pnl();
#endif

    // read in positions from file
    { ch *s = file;
    for (num i = 0; i < pointcount; i++) {
        unsigned x = grabnum(s, &s);
        s++;
        unsigned y = grabnum(s, &s);
        s++;
        unsigned z = grabnum(s, &s);
        s++;
        positions[i] = (point){ x,y,z };
    } }
    // set up initial circuits
    for (num i = 0; i < pointcount; i += 4) {
        circuits[i+0] = i+0;
        circuits[i+1] = i+1;
        circuits[i+2] = i+2;
        circuits[i+3] = i+3;
    }
    // find distances between all points
    for (num i = 0; i < pointcount; i++) {
        point pi = positions[i];
        for (num j = 0; j < pointcount; j += 4) {
            point pj0 = positions[j+0];
            point pj1 = positions[j+1];
            point pj2 = positions[j+2];
            point pj3 = positions[j+3];

            point dvp0 = (pi - pj0) * (pi - pj0);
            point dvp1 = (pi - pj1) * (pi - pj1);
            point dvp2 = (pi - pj2) * (pi - pj2);
            point dvp3 = (pi - pj3) * (pi - pj3);

            num d0 = dvp0[0] + dvp0[1] + dvp0[2] + dvp0[3];
            num d1 = dvp1[0] + dvp1[1] + dvp1[2] + dvp1[3];
            num d2 = dvp2[0] + dvp2[1] + dvp2[2] + dvp2[3];
            num d3 = dvp3[0] + dvp3[1] + dvp3[2] + dvp3[3];

            targets[i][j+0] = (typeof(**targets)) { d0, j+0 };
            targets[i][j+1] = (typeof(**targets)) { d1, j+1 };
            targets[i][j+2] = (typeof(**targets)) { d2, j+2 };
            targets[i][j+3] = (typeof(**targets)) { d3, j+3 };
        }
        // sort targets. slow.
        qsort(targets[i], pointcount, sizeof(**targets), targetcompar);
    }
    // set up initial next targets
    for (num i = 0; i < pointcount; i++) {
        next_targets[i] = targets[i][1].min_target;
        next_dists[i]   = targets[i][1].min_dist;
        next_indices[i] = 1;
    }

    unsigned last_i = 0, last_j = 0;
    // main iterations
    do {
        // find min dist
        num min_dist = -1;
        num min_i = 0;
        for (num i = 0; i < pointcount; i++)
            if (next_dists[i] < min_dist)
                min_dist = next_dists[min_i = i];
        num min_j = next_targets[min_i];
        last_i = min_i;
        last_j = min_j;
        // -- join the pair --
#if DBG
        point ipos = positions[min_i];
        point jpos = positions[min_j];
        printd(ipos[0]); print(","); printd(ipos[1]); print(","); printd(ipos[2]);
        print("\t->\t");
        printd(jpos[0]); print(","); printd(jpos[1]); print(","); printd(jpos[2]);
        pnl();
#endif
        // update circuits
        unsigned old_circuit = circuits[min_j];
        unsigned new_circuit = circuits[min_i];
        if (old_circuit != new_circuit) {
            bool done = true;
            for (num c = 0; c < pointcount; c++) {
                if (circuits[c] == old_circuit)
                    circuits[c] = new_circuit;
                if (circuits[c] != new_circuit)
                    done = false;
            }
            if (done) break;
        }
        // update targets
        unsigned niidx = ++next_indices[min_i];
        unsigned njidx = ++next_indices[min_j];
        next_targets[min_i] = targets[min_i][niidx].min_target;
        next_targets[min_j] = targets[min_j][njidx].min_target;
        next_dists[min_i] = targets[min_i][niidx].min_dist;
        next_dists[min_j] = targets[min_j][njidx].min_dist;
    } while(true);

    return positions[last_i][0] * positions[last_j][0];
}

#define RUN_TEST1 1
#define RUN_PART1 1
#define RUN_TEST2 1
#define RUN_PART2 1

#define TEST1_EXPECT 40
#define TEST2_EXPECT 25272

void run() {
#if RUN_TEST1
    print("PART 1 TEST: ");
    if (num v = do_part1(countof(test), test, true); v != TEST1_EXPECT) {
        print("FAILED (got ");
        printd(v);
        print(", expected " xstr(TEST1_EXPECT) ")\n");
    } else {
        print("PASSED\n");
    }
#endif

#if RUN_PART1
    print("PART 1 RESULT: ");
    printd(do_part1(countof(input), input, false));
    print("\n");
#endif

#if RUN_TEST2
    print("PART 2 TEST: ");
    if (num v = do_part2(countof(test), test, true); v != TEST2_EXPECT) {
        print("FAILED (got ");
        printd(v);
        print(", expected " xstr(TEST2_EXPECT) ")\n");
    } else {
        print("PASSED\n");
    }
#endif

#if RUN_PART2
    print("PART 2 RESULT: ");
    printd(do_part2(countof(input), input, false));
    print("\n");
#endif

    exit_group(0);
}