API tools faq
paste
Advertisement
Guest User

Untitled

a guest
May 17th, 2019
271
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 23.72 KB | None | 0 0
raw download clone embed print report
  1. #pragma GCC optimize("Ofast,unroll-loops,inline,omit-frame-pointer,no-asynchronous-unwind-tables,fast-math")
  2. #pragma GCC target("tune=native,avx,avx2,fma")
  3. #define WITHOUT_LI 1
  4. #ifdef CP_GEN
  5. #define INCLUDE #include
  6. INCLUDE <algorithm>
  7. INCLUDE <chrono>
  8. INCLUDE <cinttypes>
  9. INCLUDE <cstring>
  10. INCLUDE <iostream>
  11. INCLUDE <memory>
  12. INCLUDE <queue>
  13. INCLUDE <set>
  14. INCLUDE <stdio.h>
  15. INCLUDE <tuple>
  16. INCLUDE <vector>
  17. INCLUDE <cmath>
  18. #else
  19. #include <algorithm>
  20. #include <stdio.h>
  21. #include <vector>
  22. #include <queue>
  23. #include <cinttypes>
  24. #include <chrono>
  25. #include <tuple>
  26. #include <iostream>
  27. #include <memory>
  28. #include <cstring>
  29. #include <set>
  30. #include <cmath>
  31. #endif
  32.  
  33. using namespace std;
  34.  
  35. #include "macros.h"
  36. #include "types.h"
  37. #include "util.h"
  38.  
  39. #include "pdqsort.h"
  40. #include "random.h"
  41.  
  42. #ifdef CP_GEN
  43. INCLUDE <windows.h>
  44. #endif
  45.  
  46. ld getTime() {
  47. #ifndef CP_GEN
  48.   return 1.0 * clock() / CLOCKS_PER_SEC;
  49. #else
  50.     HANDLE hProcess;
  51.     hProcess = GetCurrentProcess();
  52.     FILETIME CreationTime, ExitTime, KernelTime, UserTime;
  53.     GetProcessTimes(hProcess, &CreationTime, &ExitTime, &KernelTime, &UserTime);
  54.     uint64_t kernel = ((uint64_t)KernelTime.dwHighDateTime << 32) | KernelTime.dwLowDateTime;
  55.     uint64_t user = ((uint64_t) UserTime.dwHighDateTime << 32) | UserTime.dwLowDateTime;
  56.     return (kernel + user) / 1e7;
  57. #endif
  58. }
  59.  
  60. using vd = vector<ld>;
  61. using vvd = vector<vd>;
  62. const ld EPS = 1e-9;
  63.  
  64. struct LPSolver {
  65.   int m, n;
  66.   vi N, B;
  67.   vvd D;
  68.  
  69.   LPSolver(const vvd &A, const vd &b, const vd &c) :
  70.     m(b.size()), n(c.size()), N(n + 1), B(m), D(m + 2, vd(n + 2)) {
  71.     for (int i = 0; i < m; i++)
  72.       for (int j = 0; j < n; j++) D[i][j] = A[i][j];
  73.     for (int i = 0; i < m; i++)
  74.       { B[i] = n + i; D[i][n] = -1; D[i][n + 1] = b[i]; }
  75.     for (int j = 0; j < n; j++)
  76.       { N[j] = j; D[m][j] = -c[j]; }
  77.     N[n] = -1; D[m + 1][n] = 1;
  78.   }
  79.  
  80.   void Pivot(int r, int s) {
  81.     static tuple<int, ld> L[100000], R[100000];
  82.     int nl = 0, nr = 0;
  83.     ld inv = 1.0 / D[r][s];
  84.     { for (int i = 0; i < m + 2; i++) if (i != r && (ld)abs((ld)D[i][s]) > EPS) {
  85.           L[nl++] = mt(i,D[i][s]);
  86.         }
  87.       for (int j = 0; j < n + 2; j++) if (j != s && (ld)abs((ld)D[r][j]) > EPS) {
  88.           R[nr++] = mt(j,D[r][j]);
  89.         }
  90.       FOR(j,nr) FOR(i,nl) {
  91.           D[get<0>(L[i])][get<0>(R[j])] -= get<1>(L[i]) * get<1>(R[j]) * inv;
  92.         }
  93.     }
  94.     // the previous block is a (usually) more efficient version of~:
  95.     // for (int i = 0; i < m + 2; i++) if (i != r)
  96.     //   for (int j = 0; j < n + 2; j++) if (j != s)
  97.     // D[i][j] -= D[r][j] * D[i][s] * inv;
  98.     for (int j = 0; j < n + 2; j++) if (j != s) D[r][j] *= inv;
  99.     for (int i = 0; i < m + 2; i++) if (i != r) D[i][s] *= -inv;
  100.     D[r][s] = inv;
  101.     swap(B[r], N[s]);
  102.   }
  103.  
  104.   bool Simplex(int phase) {
  105.     int x = phase == 1 ? m + 1 : m;
  106.     while (true) {
  107.       int s = -1;
  108.       for (int j = 0; j <= n; j++) {
  109.         if (phase == 2 && N[j] == -1) continue;
  110.         if (s == -1 || D[x][j] < D[x][s] ||
  111.             (D[x][j] == D[x][s] && N[j] < N[s])) s = j;
  112.       }
  113.       if (D[x][s] > -EPS) return true;
  114.       int r = -1;
  115.       for (int i = 0; i < m; i++) {
  116.         if (D[i][s] < EPS) continue;
  117.         if (r == -1 || D[i][n + 1] / D[i][s] < D[r][n + 1] / D[r][s] ||
  118.             ((D[i][n + 1] / D[i][s]) == (D[r][n + 1] / D[r][s]) &&
  119.              B[i] < B[r])) r = i;
  120.       }
  121.       if (r == -1) return false;
  122.       Pivot(r, s);
  123.     }
  124.   }
  125.  
  126.   ld Solve(vd &x) {
  127.     int r = 0;
  128.     for (int i = 1; i < m; i++) if (D[i][n + 1] < D[r][n + 1]) r = i;
  129.     if (D[r][n + 1] < -EPS) {
  130.       Pivot(r, n);
  131.       if (!Simplex(1) || D[m + 1][n + 1] < -EPS)
  132.         return -numeric_limits<ld>::infinity();
  133.       for (int i = 0; i < m; i++) if (B[i] == -1) {
  134.         int s = -1;
  135.         for (int j = 0; j <= n; j++)
  136.           if (s == -1 || D[i][j] < D[i][s] ||
  137.               (D[i][j] == D[i][s] && N[j] < N[s])) s = j;
  138.         Pivot(i, s);
  139.       }
  140.     }
  141.     if (!Simplex(2)) return numeric_limits<ld>::infinity();
  142.     x = vd(n);
  143.     for (int i = 0; i < m; i++) if (B[i] < n) x[B[i]] = D[i][n + 1];
  144.     return D[m][n + 1];
  145.   }
  146. };
  147.  
  148. const int MAXN = 2000;
  149. const int MAXT = 1000;
  150. const int MAXX = 101;
  151. const int MAXP = 7;
  152.  
  153. const int infty = 1<<28;
  154.  
  155. int n;
  156. int sx, sy;
  157.  
  158. struct location {
  159.   int x, y, p;
  160.   int duration;
  161.   int from, to;
  162.  
  163.   int score;
  164.   int sdist;
  165.  
  166.   void read(){
  167.     cin>>x>>y>>duration>>p>>from>>to;
  168.     to -= duration;
  169.     init();
  170.   }
  171.  
  172.   void init() {
  173.     score = duration * p * (p+5);
  174.     sdist = dist(sx,sy);
  175.   }
  176.  
  177.   inline int dist(location const& o) const {
  178.     return abs(x-o.x) + abs(y-o.y);
  179.   }
  180.  
  181.   inline int dist(int x_, int y_) const {
  182.     return abs(x-x_) + abs(y-y_);
  183.   }
  184. };
  185.  
  186. struct MaxFlow {
  187.   struct Edge {
  188.     li from, to, cap, flow, index;
  189.     Edge() = default;
  190.     Edge(li from, li to, li cap, li flow, li index) :
  191.       from(from), to(to), cap(cap), flow(flow), index(index) {}
  192.   };
  193.  
  194.   li N;
  195.   int ng[2+2*MAXN];
  196.   Edge G[2+2*MAXN][2+2*MAXN];
  197.   li excess[2+2*MAXN];
  198.   li dist[2+2*MAXN];
  199.   li active[2+2*MAXN];
  200.   li count[4+4*MAXN];
  201.  
  202.   queue<li> Q;
  203.  
  204.   MaxFlow() { }
  205.  
  206.   void reset() {
  207.     N = 0;
  208.     Q = queue<li>();
  209.   }
  210.  
  211.   li addNode() {
  212.     ng[N] = 0;
  213.     excess[N] = 0;
  214.     dist[N] = 0;
  215.     active[N] = 0;
  216.     count[2*N] = 0;
  217.     count[2*N+1] = 0;
  218.     N += 1;
  219.     return N-1;
  220.   }
  221.  
  222.   li addEdge(li from, li to, li cap, li flow = 0) {
  223.     li ix = ng[from];
  224.     G[from][ng[from]++] = Edge(from, to, cap, flow, ng[to] + (from == to ? 1 : 0));
  225.     G[to][ng[to]++] = Edge(to, from, 0, -flow, ix);
  226.     return ix;
  227.   }
  228.  
  229.   void enqueue(li v) {
  230.     if (!active[v] && excess[v] > 0) {
  231.       active[v] = true;
  232.       Q.push(v);
  233.     }
  234.   }
  235.  
  236.   void push(Edge &e) {
  237.     li amt = min(excess[e.from], e.cap - e.flow);
  238.     if(dist[e.from] <= dist[e.to] || amt == 0) return;
  239.     e.flow += amt;
  240.     G[e.to][e.index].flow -= amt;
  241.     excess[e.to] += amt;
  242.     excess[e.from] -= amt;
  243.     enqueue(e.to);
  244.   }
  245.  
  246.   void gap(li k) {
  247.     FOR(v, N) {
  248.       if(dist[v] < k) continue;
  249.       count[dist[v]]--;
  250.       dist[v] = max(dist[v], N+1);
  251.       count[dist[v]]++;
  252.       enqueue(v);
  253.     }
  254.   }
  255.  
  256.   void relabel(li v) {
  257.     count[dist[v]]--;
  258.     dist[v] = 2*N;
  259.     FOR(i,ng[v]) if(G[v][i].cap - G[v][i].flow > 0) {
  260.         dist[v] = min(dist[v], dist[G[v][i].to] + 1);
  261.       }
  262.     count[dist[v]]++;
  263.     enqueue(v);
  264.   }
  265.  
  266.   void discharge(li v) {
  267.     for(li i = 0; excess[v] > 0 && i < (li)ng[v]; ++i) push(G[v][i]);
  268.     if(excess[v] > 0) {
  269.       if(count[dist[v]] == 1) {
  270.         gap(dist[v]);
  271.       } else {
  272.         relabel(v);
  273.       }
  274.     }
  275.   }
  276.  
  277.   li flow(li s, li t) {
  278.     count[0] = N-1;
  279.     count[N] = 1;
  280.     dist[s] = N;
  281.     active[s] = active[t] = true;
  282.     FOR(i,ng[s]) {
  283.       excess[s] += G[s][i].cap;
  284.       push(G[s][i]);
  285.     }
  286.     while (!Q.empty()) {
  287.       li v = Q.front();
  288.       Q.pop();
  289.       active[v] = false;
  290.       discharge(v);
  291.     }
  292.     li totflow = 0;
  293.     FOR(i,ng[s]) totflow += G[s][i].flow;
  294.     return totflow;
  295.   }
  296. };
  297.  
  298. const li INF = (1<<28);
  299.  
  300. li ha[MAXN+16][MAXN+16];
  301. int u[MAXN+16+1], v[MAXN+16+1], p[MAXN+16+1], way[MAXN+16+1];
  302. int minv[MAXN+16+1], used[MAXN+16+1];
  303.  
  304. // Returns the minimum cost maximum matching
  305. // !!! !!!
  306. // !!! Hangs when n > m !!!
  307. // !!! Transpose the matrix if needed !!!
  308. // !!! !!!
  309. vi hungarian(li n, li m) {
  310.   memset(u,0,(n+1)*sizeof(int));
  311.   memset(v,0,(m+1)*sizeof(int));
  312.   memset(p,0,(m+1)*sizeof(int));
  313.   memset(way,0,(m+1)*sizeof(int));
  314.   FORU(i,1,n) {
  315.     p[0] = i;
  316.     li j0 = 0;
  317.     memset(minv,127,(m+1)*sizeof(int));
  318.     memset(used,0,(m+1)*sizeof(int));
  319.     do {
  320.       used[j0] = true;
  321.       li i0 = p[j0], delta = INF, j1 = 0;
  322.       FORU(j,1,m) if(!used[j]) {
  323.         li cur = ha[i0][j] - u[i0] - v[j];
  324.         if(cur < minv[j]) { minv[j] = cur; way[j] = j0; }
  325.         if(minv[j] < delta) { delta = minv[j]; j1 = j; }
  326.       }
  327.       FORU(j,0,m) {
  328.         if(used[j]) { u[p[j]] += delta; v[j] -= delta; }
  329.         else { minv[j] -= delta; }
  330.       }
  331.       j0 = j1;
  332.     } while(p[j0] != 0);
  333.     do {
  334.       li j1 = way[j0];
  335.       p[j0] = p[j1];
  336.       j0 = j1;
  337.     } while(j0);
  338.   }
  339.   return vi(p, p+m+1);
  340. }
  341.  
  342.  
  343. location locs[MAXN];
  344. int locsDists[MAXN][MAXN];
  345.  
  346. void read_input(){
  347.   cin >> n;
  348.   n -= 1;
  349.   cin>>sx>>sy;
  350.   { int d; cin>>d>>d>>d>>d; }
  351.   FOR(i,n) locs[i].read();
  352.   FOR(i,n) {
  353.     FOR(j,n) locsDists[i][j] = locs[i].dist(locs[j]);
  354.     pdqsort(locsDists[i], locsDists[i]+n);
  355.   }
  356.  
  357. }
  358.  
  359. const ld PHASE1_TIME = 10.0;
  360. const ld PHASE2_TIME = 12.0;
  361. const ld MAX_TIME    = 14.5;
  362.  
  363. const ld MIN_DIST = 1.0/100.0;
  364. const ld MAX_DIST = 1.0/5.0;
  365.  
  366. int inflow[MAXN], outflow[MAXN];
  367. int intime[MAXN], outtime[MAXN];
  368. array<int, 7> weight[MAXN];
  369.  
  370. int N[MAXN];
  371. int nc;
  372. int C[MAXN], IC[MAXN];
  373.  
  374. int nq;
  375. int Q[MAXN];
  376.  
  377. int flow_A[MAXN], flow_B[MAXN];
  378. int flow_nes;
  379. tuple<int,int,int> flow_ES[MAXN*MAXN];
  380. MaxFlow F;
  381.  
  382. ld flowTime = 0;
  383. ld fullFlowTime = 0;
  384. ld initTime = 0;
  385. ld fullInitTime = 0;
  386. ld hungarianTime = 0;
  387. ld fullHungarianTime = 0;
  388. ld sortTime = 0;
  389. ld removeTime = 0;
  390. ld calcTime = 0;
  391.  
  392. int ngraph[MAXN], nrgraph[MAXN];
  393. array<int,7> graph[MAXN], rgraph[MAXN];
  394.  
  395. void solve() {
  396.   li goal = 0;
  397.   FOR(i,n) goal += locs[i].p;
  398.  
  399.   li  cur = 0;
  400.   int niter = 0;
  401.  
  402.   // Initialization
  403.   auto init = [&](int maxDist) {
  404.     ld t0 = getTime();
  405.  
  406.     nc = 0;
  407.     memset(N, 0, n * sizeof(int));
  408.     FOR(i,n) FOR(j0,ngraph[i]) N[graph[i][j0]]++;
  409.     nq = 0;
  410.     FOR(i,n) if(!N[i]) Q[nq++] = i;
  411.     while(nq) {
  412.       swap(Q[random(nq)], Q[nq-1]);
  413.       int i = Q[--nq];
  414.       C[nc++] = i;
  415.       FOR(j0,ngraph[i]) {
  416.         int j = graph[i][j0];
  417.         N[j]--;
  418.         if(!N[j]) Q[nq++] = j;
  419.       }
  420.     }
  421.  
  422.     FOR(i0, n) {
  423.       int i = C[i0];
  424.       intime[i] = locs[i].from + locs[i].duration;
  425.       FOR(j0, nrgraph[i]) intime[i] = max(intime[i], intime[rgraph[i][j0]] + locs[i].dist(locs[rgraph[i][j0]]) + locs[i].duration);
  426.     }
  427.     FORD(i0,n-1,0) {
  428.       int i = C[i0];
  429.       outtime[i] = locs[i].to;
  430.       FOR(j0, ngraph[i]) outtime[i] = min(outtime[i], outtime[graph[i][j0]] - locs[i].dist(locs[graph[i][j0]]) - locs[i].duration);
  431.     }
  432.     vi time(n);
  433.     if(random(2)) {
  434.       FOR(i,n) time[i] = intime[i] - locs[i].duration;
  435.     }else{
  436.       FOR(i,n) time[i] = outtime[i];
  437.     }
  438.  
  439.     F.reset();
  440.     int S = F.addNode(), T = F.addNode();
  441.     FOR(i,n) {
  442.       flow_A[i] = F.addNode();
  443.       F.addEdge(S, flow_A[i], locs[i].p);
  444.     }
  445.     FOR(i,n) {
  446.       flow_B[i] = F.addNode();
  447.       F.addEdge(flow_B[i], T, locs[i].p);
  448.     }
  449.     flow_nes = 0;
  450.     FOR(i,n) FOR(j,n) if(i!=j && time[i] + locs[i].duration + locs[i].dist(locs[j]) <= time[j] && locs[i].dist(locs[j]) <= locsDists[i][maxDist]) {
  451.       flow_ES[flow_nes++] = mt(i,j,F.addEdge(flow_A[i], flow_B[j], min(locs[i].p, locs[j].p)));
  452.     }
  453.  
  454.     ld t1 = getTime();
  455.     F.flow(S, T);
  456.     initTime += getTime() - t1;
  457.  
  458.     cur = 0;
  459.     FOR(i,n) {
  460.       ngraph[i] = 0;
  461.       nrgraph[i] = 0;
  462.       inflow[i] = 0;
  463.       outflow[i] = 0;
  464.     }
  465.     FOR(k, flow_nes) {
  466.       int i,j,ix; tie(i,j,ix) = flow_ES[k];
  467.       int f = F.G[flow_A[i]][ix].flow;
  468.       if(f) {
  469.         cur += f;
  470.         weight[i][ngraph[i]] = f;
  471.         graph[i][ngraph[i]++] = j;
  472.         rgraph[j][nrgraph[j]++] = i;
  473.         outflow[i] += f;
  474.         inflow[j] += f;
  475.       }
  476.     }
  477.  
  478.     fullInitTime += getTime() - t0;
  479.   };
  480.  
  481.   init(n*MIN_DIST);
  482.  
  483.   while(1) {
  484.     niter++;
  485.     double t = getTime();
  486.     int phase = 1;
  487.     ld done = (ld)t/(ld)PHASE1_TIME;
  488.     if(t > PHASE1_TIME) {
  489.       done = (ld)(t-PHASE1_TIME)/(ld)(PHASE2_TIME-PHASE1_TIME);
  490.       phase = 2;
  491.     }
  492.     if(t > PHASE2_TIME) {
  493.       done = (ld)(t-PHASE2_TIME)/(ld)(MAX_TIME-PHASE2_TIME);
  494.       phase = 3;
  495.     }
  496.     if(t > MAX_TIME) break;
  497.     int maxDist = min<int>(phase == 1 ? (n * MIN_DIST) + pow(done, 1.5) * (n*MAX_DIST-n*MIN_DIST) : n*MAX_DIST, n-1);
  498.  
  499.     if(phase <= 2 && niter % 1200 == 1199) {
  500.       init(maxDist);
  501.       continue;
  502.     }
  503.  
  504.     ld t2 = getTime();
  505.     // compute a random topological sort
  506.     nc = 0;
  507.     int ty = random(2);
  508.     if(ty == 0) {
  509.       memset(N, 0, n * sizeof(int));
  510.       FOR(i,n) FOR(j0, ngraph[i]) N[graph[i][j0]]++;
  511.       nq = 0;
  512.       FOR(i,n) if(!N[i]) Q[nq++] = i;
  513.       while(nq) {
  514.         swap(Q[random(nq)], Q[nq-1]);
  515.         int i = Q[--nq];
  516.         IC[i] = nc;
  517.         C[nc++] = i;
  518.         FOR(j0, ngraph[i]) {
  519.           int j = graph[i][j0];
  520.           N[j]--;
  521.           if(!N[j]) Q[nq++] = j;
  522.         }
  523.       }
  524.     }else if(ty == 1){
  525.       memset(N, 0, n * sizeof(int));
  526.       FOR(i,n) FOR(j0, nrgraph[i]) N[rgraph[i][j0]]++;
  527.       nq = 0;
  528.       FOR(i,n) if(!N[i]) Q[nq++] = i;
  529.       while(nq) {
  530.         swap(Q[random(nq)], Q[nq-1]);
  531.         int i = Q[--nq];
  532.         IC[i] = n-1-nc;
  533.         C[n-1-(nc++)] = i;
  534.         FOR(j0, nrgraph[i]) {
  535.           int j = rgraph[i][j0];
  536.           N[j]--;
  537.           if(!N[j]) Q[nq++] = j;
  538.         }
  539.       }
  540.  
  541.     }
  542.     sortTime += getTime() - t2;
  543.  
  544.     ld t4 = getTime();
  545.     // compute min/max times
  546.     FOR(i0,n) {
  547.       int i = C[i0];
  548.       intime[i] = locs[i].from + locs[i].duration;
  549.       FOR(j0, nrgraph[i]) intime[i] = max(intime[i], intime[rgraph[i][j0]] + locs[i].dist(locs[rgraph[i][j0]]) + locs[i].duration);
  550.     }
  551.     FORD(i0,n-1,0) {
  552.       int i = C[i0];
  553.       outtime[i] = locs[i].to;
  554.       FOR(j0, ngraph[i]) outtime[i] = min(outtime[i], outtime[graph[i][j0]] - locs[i].dist(locs[graph[i][j0]]) - locs[i].duration);
  555.     }
  556.     calcTime += getTime() - t4;
  557.  
  558.     if(phase <= 2 || (phase == 3 && done < 0.5)) {
  559.       { int i = random(n);
  560.         if(ngraph[i] != 1) goto l_NoImp;
  561.         int j0 = random(ngraph[i]);
  562.         int j = graph[i][j0];
  563.         if(nrgraph[j] != 1) goto l_NoImp;
  564.         if(inflow[i] > locs[j].p) goto l_NoImp;
  565.         if(outflow[j] > locs[i].p) goto l_NoImp;
  566.         if(locs[j].dist(locs[i]) > locsDists[j][maxDist]) goto l_NoImp;
  567.         int tj = locs[j].from + locs[j].duration;
  568.         li mx0 = 0, mx1 = 0;
  569.         li mxg0 = 0, mxg1 = 0;
  570.         FOR(k0, nrgraph[i]) {
  571.           mx0 = max(mx0, locs[i].dist(locs[rgraph[i][k0]]));
  572.           mx1 = max(mx1, locs[j].dist(locs[rgraph[i][k0]]));
  573.           if(locs[rgraph[i][k0]].dist(locs[j]) > locsDists[rgraph[i][k0]][maxDist]) goto l_NoImp;
  574.           tj = max(tj, intime[rgraph[i][k0]] + locs[j].dist(locs[rgraph[i][k0]]) + locs[j].duration);
  575.         }
  576.         if(tj > locs[j].to) goto l_NoImp;
  577.         int ti = max(locs[i].from + locs[i].duration, tj + locs[i].dist(locs[j]) + locs[i].duration);
  578.         if(ti > locs[i].to) goto l_NoImp;
  579.         FOR(k0, ngraph[j]) {
  580.           mxg0 = max(mxg0, locs[j].dist(locs[graph[j][k0]]));
  581.           mxg1 = max(mxg1, locs[i].dist(locs[graph[j][k0]]));
  582.           if(locs[graph[j][k0]].dist(locs[i]) > locsDists[i][maxDist]) goto l_NoImp;
  583.           if(ti + locs[i].dist(locs[graph[j][k0]]) > outtime[graph[j][k0]]) goto l_NoImp;
  584.         }
  585.         if(max(mx0,mxg0) < max(mx1,mxg1)) goto l_NoImp;
  586.         // cout << "OK " << i << ":" << nrgraph[i] << " " << ngraph[i] << " " << j << ":" << nrgraph[j] << " " << ngraph[j] << endl;
  587.  
  588.         swap(nrgraph[i], nrgraph[j]);
  589.         swap(rgraph[i], rgraph[j]);
  590.         swap(ngraph[i], ngraph[j]);
  591.         swap(graph[i], graph[j]);
  592.         swap(inflow[i], inflow[j]);
  593.         swap(outflow[i], outflow[j]);
  594.         swap(weight[i], weight[j]);
  595.         // weight[j][i] = weight[i][j];
  596.         // weight[i][j] = 0;
  597.         rgraph[i][0] = j;
  598.         graph[j][0] = i;
  599.         FOR(k0, ngraph[i]) {
  600.           int k = graph[i][k0];
  601.           *find(begin(rgraph[k]), begin(rgraph[k]) + nrgraph[k], j) = i;
  602.         }
  603.         FOR(k0, nrgraph[j]) {
  604.           int k = rgraph[j][k0];
  605.           *find(begin(graph[k]), begin(graph[k]) + ngraph[k], i) = j;
  606.         }
  607.  
  608.         continue;
  609.       }
  610.     l_NoImp:;
  611.     }
  612.  
  613.     ld t3 = getTime();
  614.     // remove edges between the two parts
  615.     int th = random(n+1);
  616.     FOR(i,n) if(IC[i]<th) {
  617.       FOR(j0, ngraph[i]) while(j0 < ngraph[i] && IC[graph[i][j0]] >= th) {
  618.         int j = graph[i][j0];
  619.         cur -= weight[i][j0];
  620.         outflow[i] -= weight[i][j0];
  621.         inflow[j] -= weight[i][j0];
  622.         swap(graph[i][j0], graph[i][ngraph[i]-1]);
  623.         swap(weight[i][j0], weight[i][ngraph[i]-1]);
  624.         ngraph[i]--;
  625.       }
  626.     }
  627.     FOR(i,n) if(IC[i]>=th) {
  628.       FOR(j0, nrgraph[i]) while(j0 < nrgraph[i] && IC[rgraph[i][j0]] < th) {
  629.         swap(rgraph[i][j0], rgraph[i][nrgraph[i]-1]);
  630.         nrgraph[i]--;
  631.       }
  632.     }
  633.     removeTime += getTime() - t3;
  634.  
  635.     if(phase >= 2) {
  636.       ld t0 = getTime();
  637.       vi as, bs;
  638.       FOR(i0,th) FOR(k, locs[C[i0]].p - outflow[C[i0]]) as.pb(C[i0]);
  639.       FORU(i0,th,n-1) FOR(k, locs[C[i0]].p - inflow[C[i0]]) bs.pb(C[i0]);
  640.       int na = as.size(), nb = bs.size();
  641.       if(na <= nb) {
  642.         FOR(i0,na) FOR(j0,nb) {
  643.           int i = as[i0], j = bs[j0];
  644.           if(locs[i].dist(locs[j]) <= locsDists[i][maxDist] && intime[i] + locs[i].dist(locs[j]) <= outtime[j]) {
  645.             ha[i0+1][j0+1] = (phase == 2 ? 4.0 * (1 - done) * locs[i].dist(locs[j]) : 0) +
  646.               ((1000 - intime[i]) + + outtime[j]) +
  647.               ((400 - locs[i].sdist) + (-locs[j].sdist));
  648.           } else {
  649.             ha[i0+1][j0+1] = 1e6;
  650.           }
  651.         }
  652.  
  653.         ld t1 = getTime();
  654.         vi p = hungarian(na, nb);
  655.         hungarianTime += getTime() - t1;
  656.         FORU(j0,1,nb) if(p[j0]) {
  657.           int i = as[p[j0]-1];
  658.           int j = bs[j0-1];
  659.           if(intime[i] + locs[i].dist(locs[j]) <= outtime[j]) {
  660.             int j0 = 0;
  661.             while(j0 < ngraph[i] && graph[i][j0] != j) j0++;
  662.             if(j0 == ngraph[i]) {
  663.               weight[i][ngraph[i]] = 0;
  664.               graph[i][ngraph[i]++] = j;
  665.               rgraph[j][nrgraph[j]++] = i;
  666.             }
  667.             weight[i][j0] += 1;
  668.             outflow[i] += 1;
  669.             inflow[j] += 1;
  670.             cur += 1;
  671.           }
  672.         }
  673.       }else{
  674.         FOR(i0,na) FOR(j0,nb) {
  675.           int i = as[i0], j = bs[j0];
  676.           if(locs[i].dist(locs[j]) <= locsDists[i][maxDist] && intime[i] + locs[i].dist(locs[j]) <= outtime[j]) {
  677.             ha[j0+1][i0+1] = (phase == 2 ? 4.0 * (1 - done) * locs[i].dist(locs[j]) : 0) +
  678.               ((1000 - intime[i]) + + outtime[j]) +
  679.               ((400 - locs[i].sdist) + (-locs[j].sdist));
  680.           } else {
  681.             ha[j0+1][i0+1] = 1e6;
  682.           }
  683.         }
  684.  
  685.         ld t0 = getTime();
  686.         vi p = hungarian(nb, na);
  687.         hungarianTime += getTime() - t0;
  688.         FORU(i0,1,na) if(p[i0]) {
  689.           int i = as[i0-1];
  690.           int j = bs[p[i0]-1];
  691.           if(intime[i] + locs[i].dist(locs[j]) <= outtime[j]) {
  692.             int j0 = 0;
  693.             while(j0 < ngraph[i] && graph[i][j0] != j) j0++;
  694.             if(j0 == ngraph[i]) {
  695.               weight[i][ngraph[i]] = 0;
  696.               graph[i][ngraph[i]++] = j;
  697.               rgraph[j][nrgraph[j]++] = i;
  698.             }
  699.             weight[i][j0] += 1;
  700.             outflow[i] += 1;
  701.             inflow[j] += 1;
  702.             cur += 1;
  703.           }
  704.         }
  705.       }
  706.  
  707.       fullHungarianTime += getTime() - t0;
  708.  
  709.     }else{
  710.       ld t0 = getTime();
  711.       vi as, bs;
  712.       FOR(i0,th) if(outflow[C[i0]] < locs[C[i0]].p) as.pb(C[i0]);
  713.       FORU(i0,th,n-1) if(inflow[C[i0]] < locs[C[i0]].p) bs.pb(C[i0]);
  714.  
  715.       F.reset();
  716.       int S = F.addNode(), T = F.addNode();
  717.       for(int i : as) {
  718.         flow_A[i] = F.addNode();
  719.         F.addEdge(S, flow_A[i], locs[i].p - outflow[i]);
  720.       }
  721.       for(int i : bs) {
  722.         flow_A[i] = F.addNode();
  723.         F.addEdge(flow_A[i], T, locs[i].p - inflow[i]);
  724.       }
  725.       flow_nes = 0;
  726.       for(int i : as) for(int j : bs) if(locs[i].dist(locs[j]) <= locsDists[i][maxDist] && intime[i] + locs[i].dist(locs[j]) <= outtime[j]) {
  727.             flow_ES[flow_nes++] = mt(i,j,F.addEdge(flow_A[i], flow_A[j], min(locs[i].p, locs[j].p)));
  728.           }
  729.       ld t1 = getTime();
  730.       int f = F.flow(S, T);
  731.       flowTime += getTime() - t1;
  732.       (void)f;
  733.  
  734.       // add the new edges
  735.       FOR(k, flow_nes) {
  736.         int i,j,ix; tie(i,j,ix) = flow_ES[k];
  737.         li f = F.G[flow_A[i]][ix].flow;
  738.         if(f) {
  739.           weight[i][ngraph[i]] = f;
  740.           cur += f;
  741.           graph[i][ngraph[i]++] = j;
  742.           rgraph[j][nrgraph[j]++] = i;
  743.           outflow[i] += f;
  744.           inflow[j] += f;
  745.         }
  746.       }
  747.       fullFlowTime += getTime() - t0;
  748.     }
  749.  
  750.     if(niter % 100 == 0) {
  751.       cerr << "niter = " << niter << ", cur = " << goal-cur << ", done = " << done << endl;
  752.     }
  753.   }
  754.  
  755.   // Initialize the LPSolver with a feasible solution
  756.   vi O;
  757.   { vi N(n, 0);
  758.     FOR(i,n) FOR(j0, ngraph[i]) N[graph[i][j0]]++;
  759.     vi Q;
  760.     FOR(i,n) if(!N[i]) Q.pb(i);
  761.     while(!Q.empty()) {
  762.       swap(Q[random(Q.size())], Q.back());
  763.       int i = Q.back(); Q.pop_back();
  764.       O.pb(i);
  765.       FOR(j0, ngraph[i]) {
  766.         int j = graph[i][j0];
  767.         N[j]--;
  768.         if(!N[j]) Q.pb(j);
  769.       }
  770.     }
  771.     for(int i : O) {
  772.       intime[i] = locs[i].from;
  773.       FOR(j0, nrgraph[i]) intime[i] = max(intime[i], intime[rgraph[i][j0]] + locs[i].dist(locs[rgraph[i][j0]]) + locs[rgraph[i][j0]].duration);
  774.     }
  775.   }
  776.  
  777.   // Run the LPSolver
  778.   { vvd A;
  779.     vd B;
  780.     vd C(n);
  781.     FOR(i,n) C[i] = outflow[i] - inflow[i];
  782.     FOR(i,n) {
  783.       A.eb(n, 0);
  784.       A.back()[i] = 1;
  785.       B.eb(locs[i].to - intime[i]);
  786.     }
  787.     FOR(i,n) {
  788.       set<int> seen;
  789.       FOR(j0, ngraph[i]) {
  790.         int j = graph[i][j0];
  791.         if(seen.count(j)) continue;
  792.         seen.insert(j);
  793.         A.eb(n, 0);
  794.         A.back()[i] = 1;
  795.         A.back()[j] = -1;
  796.         B.eb(intime[j] - intime[i] - locs[i].duration - locs[i].dist(locs[j]));
  797.       }
  798.     }
  799.     LPSolver S(A,B,C);
  800.     vd x;
  801.     ld co = S.Solve(x);
  802.     (void)co;
  803.     FOR(i,n) intime[i] = intime[i]+x[i];
  804.   }
  805.  
  806.   li cost = 0;
  807.   FOR(i,n) cost += (2 * locs[i].p - inflow[i] - outflow[i]) * locs[i].sdist;
  808.   FOR(i,n) cost += intime[i] * inflow[i] - (intime[i] + locs[i].duration) * outflow[i] + locs[i].p * locs[i].duration;
  809.   li baseSc = 0;
  810.   FOR(i,n) if(locs[i].p) baseSc += locs[i].score;
  811.   cerr << "baseSc = " << baseSc << endl;
  812.   cerr << "nw cost = " << (goal-cur) * 240 << endl;
  813.   cerr << "route cost = " << cost << endl;
  814.   cerr << "total = " << baseSc - (goal-cur) * 240 - cost << endl;
  815.  
  816.   cerr << "flow time = " << flowTime << endl;
  817.   cerr << "full flow time = " << fullFlowTime << endl;
  818.   cerr << "init time = " << initTime << endl;
  819.   cerr << "full init time = " << fullInitTime << endl;
  820.   cerr << "hungarian time = " << hungarianTime << endl;
  821.   cerr << "full hungarian time = " << fullHungarianTime << endl;
  822.   cerr << "sort time = " << sortTime << endl;
  823.   cerr << "remove time = " << removeTime << endl;
  824.   cerr << "calc time = " << calcTime << endl;
  825.  
  826. #ifdef CP_GEN
  827.   { vvi routes;
  828.     vvi routesAt(n);
  829.     for(int i : O) {
  830.       FOR(j0, nrgraph[i]) {
  831.         int i0 = 0; while(graph[rgraph[i][j0]][i0] != i) i0++;
  832.         FOR(k, weight[rgraph[i][j0]][i0]) {
  833.           int r = routesAt[rgraph[i][j0]].back();
  834.           routesAt[rgraph[i][j0]].pop_back();
  835.           routes[r].pb(i);
  836.           routesAt[i].pb(r);
  837.         }
  838.       }
  839.       FOR(k, locs[i].p - inflow[i]) {
  840.         int r = routes.size();
  841.         routes.eb();
  842.         routes[r].pb(i);
  843.         routesAt[i].pb(r);
  844.       }
  845.     }
  846.     for(auto const& R : routes) {
  847.       cout << "start " << intime[R[0]] - locs[R[0]].sdist << " 1\n";
  848.       for(int i : R) {
  849.         cout << "arrive " << intime[i] << " " << 2+i << '\n';
  850.         cout << "work " << intime[i] << " " << intime[i]+locs[i].duration << " " << 2+i << '\n';
  851.       }
  852.       cout << "arrive " << intime[R.back()]+locs[R.back()].duration+locs[R.back()].sdist << " 1\n";
  853.       cout << "end\n";
  854.     }
  855.     cout << flush;
  856.   }
  857. #endif
  858. }
  859.  
  860. int main() {
  861.   random_reset();
  862.   read_input();
  863.   cerr << "n = " << n << endl;
  864.   solve();
  865.  
  866.   return 0;
  867. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement
Public Pastes
Advertisement
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!