Main.cpp

#include <iostream>
#include <algorithm>
#include <vector>
#include <string>
#include <cmath>
#include <set>
#include "Consts.h"
#include "DOT.h"

// Dense
// Triangle Inequality
// Directed
// Not a DAG
int llgraph[maxN][maxN];

// Lagrange Curve Fitting: 
// p(0,0) p(200, 672) p(400, 840) p(600,1008) p(800,1232) p(1000, 2464)
int cost2time(int cost) {
    //return (9/1e6) * ((long long)cost*cost*cost) - 0.012*(cost*cost) + 5.5 * cost;
    //return 3.5e-11*pow(cost, 5) - 8.17e-8*pow(cost, 4) + 7.35e-5*pow(cost, 3) - 0.0317*pow(cost, 2) + 7.364*cost;
    return cost << 1;
}
// f(x1 + x2) = f(x1) + f(x2)
// f(x) is linear
// Then: f(x) = ax + b
// a*x1 + a*x2 + b = a*x1 + a*x2 + 2b
// When b = 0, that's true

int time2cost(int time) {
    if (time <= 520) return 0;
    return std::log((time - 520) / 80.0) / 0.003;
}

struct State {
    int cost = INF;
    int prevMask = 0;
    int parent = -1;
};

struct Arborescence {
    struct ParetoPoint {
        int s = 0;
        int cost = 0;
    };
    std::vector<std::vector<int>> nodes;
    ParetoPoint paretoPoint;
};

// Longest Path Under Budget (budgeted longest path)
// Maximum path size with total weight < budget
// (0 <= src < maxN) (0 <= d <= maxW)
/*
std::vector<int> budgetedLongestPath(int src, int budget) {
    struct S {
        int budget = 0;
        int length = 0;
        int parent = -1;

        bool operator<(const S& other) const {
            if(this->budget != other.budget) return this->budget < other.budget;
            return this->length > other.length;
        }
    };

    auto domination = [](S a, S b) -> bool {
        return a.length >= b.length && a.budget <= b.budget && (a.length > b.length || a.budget < b.budget);
    };

    // dp invariant: After processing mask, dp[mask][v] contains exactly the pareto optimal elements S(mask, v)
    std::set<S> dp[maxMask+1][maxN]; // Pareto-Frontier https://codeforces.com/blog/entry/149014
    for(int mask = 1; mask <= maxMask; mask++) {
        if(!(mask&(1<<src))) continue;
        for(int i = 0; i < maxN; i++) {
            if(!(mask&(1<<i))) continue;
            std::set<S> S;
            for(int j = 0; j < maxN; j++) {
                if(!(mask&(1<<j))) continue;
                auto nS = dp[mask&(~(1<<j))][j];
                for(auto p : nS) {
                    if(p.budget + llgraph[j][i] > budget) continue;
                    p.budget += llgraph[j][i];
                    p.length++;
                    p.parent = j;
                    bool nondominated = true;
                    for(auto it = S.begin(); it != S.end() && nondominated;) {
                        bool dominated = domination(*it, p);
                        nondominated &= !dominated;
                        bool dominates = domination(p, *it);
                        if(dominates) {
                            auto aux = std::next(it);
                            S.erase(it);
                            it = aux;
                        }
                        else it = std::next(it);
                    }
                    if(nondominated) S.insert(p);
                }
            }
            
            dp[mask][i] = S;
        }
    }

    int end = -1;
    S best = {0, INF};
    for(int i = 0; i < maxN; i++)
        for(int mask = maxMask; mask; mask--)
            for(auto s : dp[mask][i])
                if(s.length > best.length || (s.length == best.length && s.budget < best.budget)) {
                    best = s;
                    end = i;
                }

    std::vector<int> path = {end};
    S node = best;
    for(int tracker = node.parent; tracker != src; tracker = node.parent) {
        path.push_back(tracker);
        //for(int mask = maxMask; mask; mask--)
            
    }

    std::reverse(path.begin(), path.end());
    return path;
}
*/

/*
    Pareto Frontier is the set of trees such that:
    There is no other tree that is: cheaper and reaches at least as many nodes
    Each frontier point answers:
    "If I want exactly k languages, what is the cheapest possible arborescence?"
    or equivalently:
    "If I can afford cost c, what is the maximum coverage?"
*/
Arborescence largestBudgetedArborescence(int src, int d, double starFactor = 0.0) {
    State dp[maxMask+1];
    dp[(1<<src)].cost = 0;
    //Invariant: dp[S] = minimum cost of some directed arborescence rooted at src that reaches exactly the nodes in S
    for(int mask = 1; mask <= maxMask; mask++) {
        if(!(mask&(1<<src)) || (dp[mask].cost != INF && dp[mask].cost > d)) continue;
        for(int i = 0; i < maxN; i++) {
            if(!(mask&(1<<i))) continue;
            for(int j = 0; j < maxN; j++) {
                if(mask&(1<<j)) continue;
                int cost = dp[mask].cost + llgraph[i][j]*(i == src ? 1.0-starFactor : 1);
                int nxMask = mask|(1<<j);
                if(dp[nxMask].cost > cost) {
                    dp[nxMask].cost = cost; 
                    dp[nxMask].parent = i;
                    dp[nxMask].prevMask = mask;
                }
            }
        }
    }

    State best;
    int bestS = 0;
    for(int mask = 1; mask <= maxMask; mask++) {
        if(dp[mask].cost > d) continue;
        int blen = __builtin_popcount(bestS);
        int curlen = __builtin_popcount(mask);
        if(blen < curlen || (blen == curlen && best.cost > dp[mask].cost)) {
            best = dp[mask];
            bestS = mask;
        }
    }

    std::vector<std::vector<int>> arborescence(maxN, std::vector<int>(maxN, 0));
    State tracker = best;
    int curmask = bestS;
    while(curmask != (1 << src)) {
        int node = __builtin_ctz(curmask ^ tracker.prevMask);
        arborescence[tracker.parent][node] = llgraph[tracker.parent][node]*(tracker.parent == src ? 1.0-starFactor : 1);
        printf("[%d|%s] = %d|%d => [%d|%s]\n", tracker.parent, langs[tracker.parent], llgraph[tracker.parent][node], arborescence[tracker.parent][node], node, langs[node]);
        curmask &= ~(1<<node);
        tracker = dp[curmask];
    }
    return {arborescence, {__builtin_popcount(bestS), dp[bestS].cost}};
}

Arborescence minCostFixedSizeArborescence(int src, int s, double starFactor = 0.0, int blacklistMask = 0, int whitelistMask = maxMask) {
    State dp[maxMask+1];
    dp[(1<<src)].cost = 0;
    //Invariant: dp[S] = minimum cost of some directed arborescence rooted at src that reaches exactly the nodes in S
    for(int mask = 1; mask <= maxMask; mask++) {
        if(!(mask&(1<<src)) || __builtin_popcount(mask) > s || (mask&blacklistMask) || (mask&(~whitelistMask))) continue;
        for(int i = 0; i < maxN; i++) {
            if(!(mask&(1<<i))) continue;
            for(int j = 0; j < maxN; j++) {
                if(mask&(1<<j)) continue;
                int cost = dp[mask].cost + llgraph[i][j]*(i == src ? 1.0-starFactor : 1);
                int nxMask = mask|(1<<j);
                if(dp[nxMask].cost > cost) {
                    dp[nxMask].cost = cost; 
                    dp[nxMask].parent = i;
                    dp[nxMask].prevMask = mask;
                }
            }
        }
    }

    State best;
    int bestS = 0;
    for(int mask = 1; mask <= maxMask; mask++) {
        int curlen = __builtin_popcount(mask);
        if(curlen != s || (mask&blacklistMask) || (mask&(~whitelistMask))) continue;
        if(best.cost > dp[mask].cost) {
            best = dp[mask];
            bestS = mask;
        }
    }

    printf("p(%d) = %d ~ %dh\n", s, dp[bestS].cost, cost2time(dp[bestS].cost));
    std::vector<std::vector<int>> arborescence(maxN, std::vector<int>(maxN, 0));
    State tracker = best;
    int curmask = bestS;
    while(curmask != (1 << src)) {
        int node = __builtin_ctz(curmask ^ tracker.prevMask);
        arborescence[tracker.parent][node] = llgraph[tracker.parent][node]*(tracker.parent == src ? 1.0-starFactor : 1);
        printf("[%d|%s] = %d|%d => [%d|%s]\n", tracker.parent, langs[tracker.parent], llgraph[tracker.parent][node], arborescence[tracker.parent][node], node, langs[node]);
        curmask &= ~(1<<node);
        tracker = dp[curmask];
    }
    return {arborescence, {s, dp[bestS].cost}};
}

// Steiner Path with Cardinality Constraint
// Shortest path that contains L, U, and ≥ s nodes
std::vector<int> steinerPath(int src, int end, int s, int blacklistMask = 0, int whitelistMask = maxMask) {
    State dp[maxMask+1][maxN];
    dp[(1<<src)][src].cost = 0;
    for(int mask = 1; mask <= maxMask; mask++) {
        if(!(mask&(1<<src)) || __builtin_popcount(mask) > s+2 || (mask&blacklistMask) || (mask&(~whitelistMask))) continue;
        for(int i = 0; i < maxN; i++) {
            if(!(mask&(1<<i)) || i == end) continue;
            for(int j = 0; j < maxN; j++) {
                if(mask&(1<<j)) continue;
                int cost = dp[mask][i].cost + llgraph[i][j];
                int nxMask = mask|(1<<j);
                if(dp[nxMask][j].cost > cost) {
                    dp[nxMask][j].cost = cost;
                    dp[nxMask][j].parent = i;
                    dp[nxMask][j].prevMask = mask;
                }
            }
        }
    }

    State best;
    int bestS = 0;
    for(int mask = 1; mask <= maxMask; mask++) {
        if(!(mask&(1<<src)) || !(mask&(1<<end)) || __builtin_popcount(mask) != s+2 || (mask&blacklistMask) || (mask&(~whitelistMask))) continue;
        if(best.cost > dp[mask][end].cost) {
            best = dp[mask][end];
            bestS = mask;
        }
    }

    std::vector<int> path;
    printf("([%d|%s] => [%d|%s]) p(%d) = %d ~ %dh\n", src, langs[src], end, langs[end], s+2, best.cost, cost2time(best.cost));
    State tracker = best;
    int curmask = bestS;
    while(curmask != (1<<src)) {
        int node = __builtin_ctz(curmask ^ tracker.prevMask);
        path.push_back(node);
        curmask = tracker.prevMask;
        tracker = dp[tracker.prevMask][tracker.parent];
    }
    path.push_back(src);
    std::reverse(path.begin(), path.end());
    for(int i = 0; i < path.size()-1; i++)
        printf("[%d|%s] = %d => ", path[i], langs[path[i]], llgraph[path[i]][path[i+1]]);
    printf("[%d|%s]\n", path.back(), langs[path.back()]);
    return path;
}

int getEasiestLanguage() {
    int lang = -1, minD = INF;
    for(int i = 0; i < maxN; i++) {
        int d = 0;
        for(int j = 0; j < maxN; j++)
            d += llgraph[j][i];
        if(minD > d) {
            minD = d;
            lang = i;
        }
    }

    printf("%s is the easiest language with ~%dh to learn\n", langs[lang], cost2time(minD/maxN));
    return lang;
}

int getHardestLanguage() {
    int lang = -1, maxD = 0;
    for(int i = 0; i < maxN; i++) {
        int d = 0;
        for(int j = 0; j < maxN; j++)
            d += llgraph[j][i];
        if(maxD < d) {
            maxD = d;
            lang = i;
        }
    }

    printf("%s is the hardest language with ~%dh to learn\n", langs[lang], cost2time(maxD/maxN));
    return lang;
}

void buildLinearGraph() {
    for(int i = 0; i < maxN; i++)
        for(int j = 0; j < maxN; j++) {
            double base = (1000-similarity[i][j]);
            double support = 1.0 + 0.5*(0.1*(popularity[j]/1000.0) + 0.9*(presence[j]/1000.0));
            llgraph[i][j] = base / support;
        }
    DOT::Save(llgraph);
}

void buildGraph() {
    for(int i = 0; i < maxN; i++)
        for(int j = 0; j < maxN; j++) {
            const int scriptDist = 20*(script[j]-script[i]);
            const int morphDist = 25*(morphComplex[j]-morphComplex[i]);
            const int tonalDist = 35*(toneLevel[j]-toneLevel[i]);
            const int woDist = 20*(wordOrder[j] != wordOrder[i]);
            const double simGap = 1000.0 - similarity[i][j];
            const double typology = scriptDist + morphDist + tonalDist + woDist;
            const double beta = 0.0015; // [0.001, 0.01]
            const double base = simGap * exp(beta * (typology-80));
            const double alpha = 0.0012;
            const double simf = 1400 * (1 - exp(-alpha * base));
            const double transferFactor = (similarity[i][j] / 1000.0) + 0.1;
            const double support = 1.0 + transferFactor * 0.6 * (0.1*(popularity[j]/1000.0) + 0.9*(presence[j]/1000.0));
            llgraph[i][j] = simf / support;
        }
    DOT::Save(llgraph);
}

void header() {
    for(int i = 0; i < maxN; i++)
        printf(" %d) %s\n", i+1, langs[i]);
}

int main() {
    header();
    buildGraph();

    puts("Largest Budgeted");
    largestBudgetedArborescence(LANGUAGE::PORTUGUESE, 3000, 0.1);

    std::vector<Arborescence::ParetoPoint> paretoFrontier(maxN+1);
    llgraph[PORTUGUESE][ENGLISH] = 70;
    puts("Min Cost Fixed Size");
    //for(int s = 1; s <= maxN; s++)
    //     paretoFrontier[s] = minCostFixedSizeArborescence(LANGUAGE::PORTUGUESE, s, 0.0, 0, (1<<PORTUGUESE)|(1<<CHINESE)|(1<<JAPANESE)|(1<<ENGLISH)|(1<<GERMAN)|(1<<RUSSIAN)|(1<<HEBREW)).paretoPoint;
    
    for(int s = 0; s <= maxN; s++) {
        printf("pp(%d, %d) | pp(%d, %d)\n", paretoFrontier[s].s, paretoFrontier[s].cost, paretoFrontier[s].s, cost2time(paretoFrontier[s].cost));
    }

    for(int s = 0; s <= maxN; s++)
        printf("%d, ", cost2time(paretoFrontier[s].cost));
    putchar('\n');

    puts("Steiner Path");
    steinerPath(LANGUAGE::PORTUGUESE, LANGUAGE::CHINESE, 2, 0, (1<<PORTUGUESE)|(1<<CHINESE)|(1<<KOREAN)|(1<<JAPANESE)|(1<<ENGLISH)|(1<<GERMAN)|(1<<RUSSIAN)|(1<<HEBREW));
    steinerPath(LANGUAGE::JAPANESE, LANGUAGE::CHINESE, 1);
    steinerPath(LANGUAGE::JAPANESE, LANGUAGE::CHINESE, 0);
    getEasiestLanguage();
    getHardestLanguage();
    return 0;
}