Advertisement
Guest User

TCO16R2 MM by Psyho

a guest
May 25th, 2016
636
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 37.67 KB | None | 0 0
  1. // Author: Psyho
  2. // Blog: http://psyho.gg/
  3. // Twitter: https://twitter.com/fakepsyho
  4.  
  5. const double TIME_MULT = 1.0;
  6. const bool FULL_EVALUATE = false;
  7.  
  8. #include <bits/stdc++.h>
  9. // #include <time.h>
  10.  
  11. using namespace std;
  12.  
  13. #undef assert
  14. #define assert(x) ;
  15.  
  16. #define INLINE   inline __attribute__ ((always_inline))
  17. #define NOINLINE __attribute__ ((noinline))
  18.  
  19. #define ALIGNED __attribute__ ((aligned(16)))
  20.  
  21. #define likely(x)   __builtin_expect(!!(x),1)
  22. #define unlikely(x) __builtin_expect(!!(x),0)
  23.  
  24. #define SSELOAD(a)     _mm_load_si128((__m128i*)&a)
  25. #define SSESTORE(a, b) _mm_store_si128((__m128i*)&a, b)
  26.  
  27. #define FOR(i,a,b)  for(int i=(a);i<(b);++i)
  28. #define REP(i,a)    FOR(i,0,a)
  29. #define ZERO(m)     memset(m,0,sizeof(m))
  30. #define ALL(x)      x.begin(),x.end()
  31. #define PB          push_back
  32. #define S           size()
  33. #define byte        unsigned char
  34. #define LL          long long
  35. #define ULL         unsigned long long
  36. #define LD          long double
  37. #define MP          make_pair
  38. #define X           first
  39. #define Y           second
  40. #define VC          vector
  41. #define PII         pair<int, int>
  42. #define PDD         pair<double, double>
  43. #define VI          VC<int>
  44. #define VVI         VC<VI>
  45. #define VVVI        VC<VVI>
  46. #define VPII        VC<PII>
  47. #define VVPII       VC<VPII>
  48. #define VVVPII      VC<VVPII>
  49. #define VD          VC<double>
  50. #define VVD         VC<VD>
  51. #define VVVD        VC<VVD>
  52. #define VPDD        VC<PDD>
  53. #define VVPDD       VC<VPDD>
  54. #define VVVPDD      VC<VVPDD>
  55. #define VS          VC<string>
  56. #define VVS         VC<VS>
  57. #define VVVS        VC<VVS>
  58. #define DB(a)       cerr << #a << ": " << (a) << endl;
  59.  
  60. template<class A, class B> ostream& operator<<(ostream &os, pair<A,B> p) {os << "(" << p.X << "," << p.Y << ")"; return os;}
  61. template<class A, class B, class C> ostream& operator<<(ostream &os, tuple<A,B,C> p) {os << "(" << get<0>(p) << "," << get<1>(p) << "," << get<2>(p) << ")"; return os;}
  62. template<class T> ostream& operator<<(ostream &os, VC<T> v) {os << "{"; REP(i, v.S) {if (i) os << ", "; os << v[i];} os << "}"; return os;}
  63. template<class T> ostream& operator<<(ostream &os, set<T> s) {VS vs(ALL(s)); return os << vs;}
  64. template<class T> string i2s(T x) {ostringstream o; o << x; return o.str();}
  65. VS splt(string s, char c = ' ') {VS all; int p = 0, np; while (np = s.find(c, p), np >= 0) {all.PB(s.substr(p, np - p)); p = np + 1;} all.PB(s.substr(p)); return all;}
  66.  
  67. double getTime() {
  68.     ULL timelo, timehi;
  69.     __asm__ volatile ("rdtsc" : "=a" (timelo), "=d" (timehi));
  70.     return ((timehi << 32) + timelo) / 2.5e9;
  71.     // timespec tv;
  72.     // clock_gettime(CLOCK_MONOTONIC, &tv);
  73.     // return tv.tv_sec + 1e-9 * tv.tv_nsec;
  74. }
  75.  
  76. struct RNG {
  77.     unsigned int x = 123456789;
  78.     unsigned int y = 362436069;
  79.     unsigned int z = 521288629;
  80.    
  81.     unsigned int rand() {
  82.         x ^= x << 16;
  83.         x ^= x >> 5;
  84.         x ^= x << 1;
  85.  
  86.         unsigned int t = x;
  87.         x = y;
  88.         y = z;
  89.         z = t ^ x ^ y;
  90.  
  91.         return z;      
  92.     }
  93.    
  94.     INLINE int next(int x) {
  95.         return rand() % x;
  96.     }
  97.    
  98.     INLINE int nextAnd(int x) {
  99.         return rand() & x;
  100.     }
  101.    
  102.     INLINE int next(int a, int b) {
  103.         return a + (rand() % (b - a));
  104.     }
  105.    
  106.     INLINE double nextDouble() {
  107.         return (rand() + 0.5) * (1.0 / 4294967296.0);
  108.     }
  109.    
  110.     template<class T>
  111.     INLINE void shuffle(VC<T> &v, int end = -1) {
  112.         if (end == -1) end = v.S;
  113.         REP(i, end) swap(v[i], v[next(i, v.S)]);
  114.     }
  115. };
  116.  
  117. RNG rng;
  118.  
  119. template <int MAXV>
  120. class MCMF {
  121.     struct Edge {
  122.         int target;
  123.         int flow;
  124.         int cost;
  125.         int reverse;
  126.        
  127.         Edge(int target, int flow, int cost, int reverse) {
  128.             this->target = target;
  129.             this->flow = flow;
  130.             this->cost = cost;
  131.             this->reverse = reverse;
  132.         }
  133.     };
  134.  
  135.     VC<Edge> edges[MAXV];
  136.     int prev[MAXV];
  137.     int dist[MAXV];
  138.     int queueExist[MAXV];
  139.     int queueData[MAXV + 1];
  140.     int queueStart;
  141.     int queueEnd;
  142.    
  143.     int source;
  144.     int sink;
  145.  
  146. public:
  147.     MCMF(int source, int sink) {
  148.         this->source = source;
  149.         this->sink = sink;
  150.     }
  151.        
  152.     void addEdge(int a, int b, int flow, int cost = 0) {
  153.         int sa = edges[a].S;
  154.         int sb = edges[b].S;
  155.         edges[a].PB(Edge(b, flow, cost, sb));
  156.         edges[b].PB(Edge(a, 0, -cost, sa));
  157.     }
  158.  
  159.     int checkFlow(int a, int b) {
  160.         for (Edge &edge : edges[b]) if (edge.target == a) return edge.flow;
  161.         return 0;
  162.     }
  163.  
  164. private:
  165.     bool queueDone() {
  166.         return queueStart == queueEnd;
  167.     }
  168.  
  169.     void queueClear() {
  170.         ZERO(queueExist);
  171.         queueStart = 0;
  172.         queueEnd = 0;
  173.     }
  174.  
  175.     int queuePop() {
  176.         int rv = queueData[queueStart];
  177.         queueExist[rv] = false;
  178.         queueStart = (queueStart + 1) % (MAXV + 1);
  179.         return rv;
  180.     }
  181.  
  182.     void queueAdd(int v) {
  183.         if (queueExist[v]) return;
  184.         queueExist[v] = true;
  185.         queueData[queueEnd] = v;
  186.         queueEnd = (queueEnd + 1) % (MAXV + 1);
  187.     }
  188.    
  189. public:
  190.     LL run() {
  191.         LL totalCost = 0;
  192.         int totalFlow = 0;
  193.  
  194.         while (true) {
  195.             memset(dist, 0x3F, sizeof(dist));
  196.             dist[source] = 0;
  197.             prev[sink] = -1;
  198.             prev[source] = -1;
  199.             queueClear();
  200.             queueAdd(source);
  201.             while (!queueDone()) {
  202.                 int v = queuePop();
  203.                 REP(i, edges[v].S) if (edges[v][i].flow > 0 && dist[edges[v][i].target] > dist[v] + edges[v][i].cost) {
  204.                     dist[edges[v][i].target] = dist[v] + edges[v][i].cost;
  205.                     prev[edges[v][i].target] = edges[v][i].reverse;
  206.                     queueAdd(edges[v][i].target);
  207.                 }
  208.             }
  209.            
  210.             if (prev[sink] == -1 || dist[sink] > (1 << 20)) break;
  211.            
  212.             int maxFlow = 1 << 30;
  213.             int v = sink;
  214.             while (v != source) {
  215.                 int pv = edges[v][prev[v]].target;
  216.                 maxFlow = min(maxFlow, edges[pv][edges[v][prev[v]].reverse].flow);
  217.                 v = pv;
  218.             }
  219.            
  220.             assert(maxFlow > 0);
  221.            
  222.             totalFlow += maxFlow;
  223.            
  224.             v = sink;
  225.             while (v != source) {
  226.                 int pv = edges[v][prev[v]].target;
  227.                 edges[pv][edges[v][prev[v]].reverse].flow -= maxFlow;
  228.                 edges[v][prev[v]].flow += maxFlow;
  229.                 totalCost += edges[pv][edges[v][prev[v]].reverse].cost * maxFlow;
  230.                 v = pv;
  231.             }
  232.            
  233.         }
  234.        
  235.         return totalCost;
  236.     }
  237.  
  238. };
  239.  
  240. const int MAX_STARS = 2000;
  241. const int MAX_UFOS = 20;
  242. const int MAX_SHIPS = 10;
  243. const int MAX_NEIGHBOURS = 30;
  244. const double UFO_MULT = 0.001;
  245.  
  246. PII PP[MAX_STARS];
  247. double OPD[MAX_STARS][MAX_STARS+4];
  248. double PD[MAX_STARS][MAX_STARS+4];
  249. int PN[MAX_STARS][MAX_STARS];
  250. int PO[MAX_STARS][MAX_STARS];
  251. int PC[MAX_STARS][MAX_NEIGHBOURS];
  252. PII PPOS[MAX_STARS];
  253. int PNO;
  254.  
  255. int SP[MAX_SHIPS];
  256. int SNO;
  257.  
  258. VI UFOS[MAX_UFOS];
  259. int UNO;
  260. int UPROBTURN = -1000;
  261. VI UPROB[MAX_UFOS];
  262. int UFOUSED[MAX_UFOS];
  263.  
  264. double UFODIST[MAX_UFOS];
  265.  
  266. bool USEREALUFOW = false;
  267. int REALUFOW[MAX_UFOS];
  268.  
  269. double PEXPCOST;
  270.  
  271. int XVSNO = 1;
  272. int XVS[MAX_STARS];
  273. int vs[MAX_STARS];
  274. int vscopy[MAX_STARS];
  275. int totvs = 0;
  276.  
  277. int turn = -1;
  278.  
  279. const int MCMST_SIMS = 25;
  280. const int UFO_DIST_SIZE = 32768;
  281. // const double TOLERANCE = 100;
  282.  
  283. const double TIME_LIMIT = 19.5 * TIME_MULT;
  284. const double MAX_MC_TIME = 10.0 * TIME_MULT;
  285. const double MAX_MCMST_TIME = 15.0 * TIME_MULT;
  286.  
  287. int NEXTMOVE[MAX_SHIPS];
  288.  
  289. double totalTime = 0;
  290.  
  291. int TPATH[MAX_STARS+10];
  292. int PATH[MAX_SHIPS][MAX_STARS+10];
  293. int PATHLEN[MAX_SHIPS];
  294. int PATHTURN[MAX_SHIPS];
  295. int BPATH[MAX_SHIPS][MAX_STARS+10];
  296. int BPATHLEN[MAX_SHIPS];
  297. double BPATHSUM[MAX_SHIPS][MAX_STARS+10];
  298.  
  299. bool FINAL_PATH = false;
  300.  
  301. int hops = 0;
  302. int detours = 0;
  303. int hopped[MAX_SHIPS];
  304.  
  305. VI bestMatching(VVI &w) {
  306.     assert(w.S == SNO);
  307.     REP(i, SNO) assert(w[i].S == UNO);
  308.    
  309.     MCMF<42> mcmf(40, 41);
  310.     REP(i, SNO) REP(j, UNO) if (w[i][j] >= 0) mcmf.addEdge(i, 20 + j, 1, w[i][j]);
  311.     REP(i, SNO) mcmf.addEdge(40, i, 1, 0);
  312.     REP(i, UNO) mcmf.addEdge(20 + i, 41, 1, 0);
  313.     mcmf.run();
  314.    
  315.     VI rv(SNO, -1);
  316.     REP(i, SNO) REP(j, UNO) if (mcmf.checkFlow(i, 20 + j)) rv[i] = j;
  317.     return rv;
  318. }
  319.  
  320. struct Strategy {
  321.     int hopSize;
  322.     bool transport;
  323.     VI forceUfo;
  324.    
  325.     Strategy() {
  326.         hopSize = 0;
  327.         transport = false;
  328.     }
  329. };
  330.  
  331. Strategy str;
  332.  
  333.  
  334. double calcSolution(int s) {
  335.     double rv = 0;
  336.     REP(i, PATHLEN[s]-1) rv += PD[PATH[s][i]][PATH[s][i+1]];
  337.     return rv;
  338. }
  339.  
  340. double calcSolution() {
  341.     double rv = 0;
  342.     REP(s, SNO) REP(i, PATHLEN[s]-1) rv += PD[PATH[s][i]][PATH[s][i+1]];
  343.     return rv;
  344. }
  345.  
  346. INLINE double pathDist(int s, int pos, int p) {
  347.     return PD[PATH[s][pos]][p];
  348. }
  349.  
  350. double greedy(VI pos, VI planets) {
  351.     double rv = 0;
  352.     REP(s, SNO) {
  353.         PATHLEN[s] = 1;
  354.         PATH[s][0] = pos[0];
  355.     }
  356.    
  357.     while (planets.S) {
  358.         double bv = 1e30;
  359.         int bs = -1;
  360.         int bp = -1;
  361.         int bi = -1;
  362.         REP(s, SNO) FOR(p, 1, PATHLEN[s]+1) REP(i, planets.S) {
  363.             double av = 0;
  364.             av += PD[PATH[s][p-1]][planets[i]];
  365.             if (p < PATHLEN[s]) av += PD[PATH[s][p]][planets[i]];
  366.             if (av < bv) {
  367.                 bv = av;
  368.                 bs = s;
  369.                 bp = p;
  370.                 bi = i;
  371.             }
  372.         }
  373.         PATHLEN[bs]++;
  374.         for (int i = PATHLEN[bs]; i > bp; i--) PATH[bs][i] = PATH[bs][i-1];
  375.         PATH[bs][bp] = planets[bi];
  376.         rv += bv;
  377.         planets[bi] = planets.back();
  378.         planets.pop_back();
  379.     }
  380.     return rv;
  381. }
  382.  
  383. double quickLB(VI &vs, VI &ships) {
  384.     double rv = 0;
  385.     double nearShip = 1e30;
  386.     REP(i, PNO) if (!vs[i]) {
  387.         int found = 0;
  388.         double d1 = 0;
  389.         double d2 = 0;
  390.         FOR(j, 1, PNO) {
  391.             int p = PN[i][j];
  392.             if (!vs[p]) {
  393.                 rv += PD[i][p];
  394.                 if (++found >= 2) break;
  395.             }
  396.         }
  397.         REP(j, ships.S) nearShip = min(nearShip, PD[i][ships[j]]);
  398.     }
  399.     rv /= 2;
  400.     return rv + (nearShip > 1e10 ? 0 : nearShip);
  401. }
  402.  
  403. VI genDistribution(int no, int x, int size) {
  404.     VI rv(size);
  405.     REP(i, size) {
  406.         rv[i] = no;
  407.         REP(j, x) rv[i] = min(rv[i], rng.next(1, no));
  408.     }
  409.     return rv;
  410. }
  411.  
  412. struct State {
  413.     VI vs;
  414.     VI ships;
  415.     VI hopped;
  416.     VVI ufo;
  417.     VI ufow;
  418.     VVI ufoprob;
  419.     int offset;
  420.     int turnsLeft;
  421. };
  422.  
  423. int FUP[MAX_STARS];
  424. int FUR[MAX_STARS];
  425.  
  426. void FUinit(int n) {
  427.     REP(i, n) FUP[i] = i;
  428.     REP(i, n) FUR[i] = 0;
  429. }
  430.  
  431.  
  432. void FUunion(int x, int y) {
  433.     int xr = FUP[x];
  434.     int yr = FUP[y];
  435.     if (xr == yr) return;
  436.     if (FUR[xr] < FUR[yr]) {
  437.         FUP[xr] = yr;
  438.     } else if (FUR[xr] > FUR[yr]) {
  439.         FUP[yr] = xr;
  440.     } else {
  441.         FUP[yr] = xr;
  442.         FUR[xr]++;
  443.     }
  444. }
  445.  
  446. int FUfind(int x) {
  447.     return FUP[x] == x ? x : FUP[x] = FUfind(FUP[x]);
  448. }
  449.  
  450. double calcMST() {
  451.     VI v;
  452.     REP(i, PNO) if (!vs[i]) v.PB(i);
  453.     if (v.S <= 1) return 0;
  454.    
  455.     FUinit(v.S);
  456.    
  457.     priority_queue<pair<double,PII>> pq;
  458.     REP(i, v.S) REP(j, i) pq.push(MP(-PD[v[i]][v[j]], MP(i, j)));
  459.    
  460.     double rv = 0;
  461.     int edgesAdded = 1;
  462.     while (edgesAdded < v.S) {
  463.         auto p = pq.top(); pq.pop();
  464.         if (FUfind(p.Y.X) == FUfind(p.Y.Y)) continue;
  465.         rv += -p.X;
  466.         edgesAdded++;
  467.         FUunion(p.Y.X, p.Y.Y);
  468.     }
  469.     return rv;
  470. }
  471.  
  472. double calcMST(VI &ships) {
  473.     VI v;
  474.     REP(i, PNO) if (!vs[i]) v.PB(i);
  475.     if (v.S == 0) return 0;
  476.    
  477.     FUinit(v.S+1);
  478.    
  479.     priority_queue<pair<double,PII>> pq;
  480.     REP(i, v.S) REP(j, i) pq.push(MP(-PD[v[i]][v[j]], MP(i, j)));
  481.     REP(i, v.S) {
  482.         double minDist = 1e30;
  483.         REP(j, ships.S) minDist = min(minDist, PD[v[i]][ships[j]]);
  484.         pq.push(MP(-minDist, MP(i, v.S)));
  485.     }
  486.    
  487.     double rv = 0;
  488.     int edgesAdded = 0;
  489.     while (edgesAdded < v.S) {
  490.         auto p = pq.top(); pq.pop();
  491.         if (FUfind(p.Y.X) == FUfind(p.Y.Y)) continue;
  492.         rv += -p.X;
  493.         edgesAdded++;
  494.         FUunion(p.Y.X, p.Y.Y);
  495.     }
  496.     return rv;
  497. }
  498.  
  499. double calcMST(VI &vs, VI &ships) {
  500.     VI v;
  501.     REP(i, PNO) if (!vs[i]) v.PB(i);
  502.     if (v.S == 0) return 0;
  503.    
  504.     FUinit(v.S+1);
  505.    
  506.     priority_queue<pair<double,PII>> pq;
  507.     REP(i, v.S) REP(j, i) pq.push(MP(-PD[v[i]][v[j]], MP(i, j)));
  508.     REP(i, v.S) {
  509.         double minDist = 1e30;
  510.         REP(j, ships.S) minDist = min(minDist, PD[v[i]][ships[j]]);
  511.         pq.push(MP(-minDist, MP(i, v.S)));
  512.     }
  513.    
  514.     double rv = 0;
  515.     int edgesAdded = 0;
  516.     while (edgesAdded < v.S) {
  517.         auto p = pq.top(); pq.pop();
  518.         if (FUfind(p.Y.X) == FUfind(p.Y.Y)) continue;
  519.         rv += -p.X;
  520.         edgesAdded++;
  521.         FUunion(p.Y.X, p.Y.Y);
  522.     }
  523.     return rv;
  524. }
  525.  
  526. const int MCMST_MAX_TURNS = 300;
  527. const int MCMST_MAX_PLANETS = 300;
  528. #define PDII pair<double, PII>
  529. int MCMSTEDGESNO = 0;
  530. int MCMSTXEDGESNO = 0;
  531. PDII MCMSTEDGES[MCMST_MAX_PLANETS * MCMST_MAX_PLANETS + 4];
  532. PDII MCMSTXEDGES[MCMST_MAX_PLANETS];
  533. double calcMCMST(VI final, VI ufo, VI ufoused, int turns, int sims) {
  534.     VC<RNG> r(UNO); REP(u, ufo.S) REP(i, turns) r[u].rand();
  535.    
  536.     double rv = 0;
  537.     double ud = 0;
  538.     REP(i, PNO) vscopy[i] = vs[i];
  539.    
  540.     VI ov;
  541.     REP(i, PNO) if (!vs[i]) ov.PB(i);
  542.    
  543.     MCMSTEDGESNO = 0;
  544.     REP(i, ov.S) REP(j, i) MCMSTEDGES[MCMSTEDGESNO++] = MP(PD[ov[i]][ov[j]], MP(i, j));
  545.     sort(MCMSTEDGES, MCMSTEDGES + MCMSTEDGESNO);   
  546.     MCMSTEDGES[MCMSTEDGESNO].X = 1e30;
  547.    
  548.     REP(sim, sims) {
  549.         REP(i, PNO) vs[i] = vscopy[i];
  550.         VI pos = final;
  551.         REP(i, UNO) if (ufoused[i]) {
  552.             int x = ufo[i];
  553.             REP(t, turns) {
  554.                 int y = PN[x][UPROB[i][r[i].nextAnd(UFO_DIST_SIZE-1)]];
  555.                 // ud += PD[x][y];
  556.                 x = y;
  557.                 vs[x] = 1;
  558.             }
  559.             pos.PB(x);
  560.         }
  561.         assert(pos.S == SNO);
  562.        
  563.         FUinit(ov.S+1);    
  564.        
  565.         MCMSTXEDGESNO = 0;
  566.         REP(i, ov.S) if (!vs[ov[i]]) {
  567.             double minDist = 1e30;
  568.             REP(j, pos.S) minDist = min(minDist, PD[ov[i]][pos[j]]);
  569.             MCMSTXEDGES[MCMSTXEDGESNO++] = MP(minDist, MP(i, ov.S));
  570.         }
  571.         sort(MCMSTXEDGES, MCMSTXEDGES + MCMSTXEDGESNO);
  572.         MCMSTXEDGES[MCMSTXEDGESNO].X = 1e30;
  573.        
  574.         int EPOS = 0;
  575.         int XEPOS = 0;
  576.        
  577.         int edgesAdded = 0;
  578.         while (edgesAdded < MCMSTXEDGESNO) {
  579.             if (MCMSTEDGES[EPOS].X < MCMSTXEDGES[XEPOS].X) {
  580.                 auto &p = MCMSTEDGES[EPOS++];
  581.                 if (vs[ov[p.Y.X]] || vs[ov[p.Y.Y]]) continue;
  582.                 if (FUfind(p.Y.X) == FUfind(p.Y.Y)) continue;
  583.                 rv += p.X;
  584.                 edgesAdded++;
  585.                 FUunion(p.Y.X, p.Y.Y);
  586.             } else {
  587.                 auto &p = MCMSTXEDGES[XEPOS++];
  588.                 if (vs[ov[p.Y.X]]) continue;
  589.                 if (FUfind(p.Y.X) == FUfind(p.Y.Y)) continue;
  590.                 rv += p.X;
  591.                 edgesAdded++;
  592.                 FUunion(p.Y.X, p.Y.Y);
  593.             }
  594.         }
  595.     }
  596.     REP(i, PNO) vs[i] = vscopy[i];
  597.     // rv += ud * UFO_MULT;
  598.     rv /= sims;
  599.     return rv;
  600. }
  601.  
  602.  
  603.  
  604. double monteCarlo(State &state, Strategy &str) {
  605.     VI vs = state.vs;
  606.    
  607.     int planetsLeft = 0;
  608.     REP(i, PNO) if (!vs[i]) planetsLeft++;
  609.     if (planetsLeft == 0) return 0;
  610.    
  611.     int opl = planetsLeft;
  612.    
  613.     int hopSize = str.hopSize;
  614.    
  615.     double rv = 0;
  616.     VI ufos = state.ufo[0];
  617.     VI next = state.ufo[0];
  618.     VI ships = state.ships;
  619.     VI usedUfo(UNO, -1);
  620.     VI hopped = state.hopped;
  621.    
  622.     int turn = 0;
  623.     // First Turn
  624.     VVI w(SNO, VI(UNO, -1));
  625.     REP(s, SNO) REP(u, UNO) {
  626.         double minDist = 1e30;
  627.         if (ships[s] == state.ufo[0][u]) minDist = 0;
  628.         minDist = min(minDist, PD[ships[s]][state.ufo[1][u]]);
  629.         minDist = min(minDist, PD[ships[s]][state.ufo[2][u]]);
  630.         if (minDist > hopSize && !str.forceUfo[u]) continue;
  631.         if (hopped[s] && !str.forceUfo[u] && minDist > 0) continue;
  632.         w[s][u] = (int)round(minDist) + (str.forceUfo[u] ? 0 : 2000);
  633.     }
  634.     VI matching = bestMatching(w);
  635.     REP(s, SNO) if (matching[s] != -1) {
  636.         int u = matching[s];
  637.         int p = ships[s] == state.ufo[0][u] || PD[ships[s]][state.ufo[1][u]] < PD[ships[s]][state.ufo[2][u]] ? state.ufo[1][u] : state.ufo[2][u];
  638.         rv += PD[ships[s]][p] * (ships[s] == state.ufo[0][u] ? UFO_MULT : 1.0);
  639.         hopped[s] = 1;
  640.         ships[s] = p;
  641.         planetsLeft -= vs[p] == 0;
  642.         vs[p] = 1;
  643.     }
  644.     ufos = state.ufo[1];
  645.     turn++;
  646.    
  647.    
  648.    
  649.     //Rest
  650.     while (state.turnsLeft - turn > planetsLeft && planetsLeft) {
  651.         if (turn >= 2) {
  652.             REP(i, UNO) next[i] = PN[ufos[i]][state.ufoprob[i][turn]];
  653.         } else if (turn == 1) {
  654.             next = state.ufo[2];
  655.         }
  656.        
  657.         int firstShip = 0;
  658.         REP(u, UNO) {
  659.             double bv = 1e30;
  660.             int bs = -1;
  661.             FOR(s, firstShip, SNO) {
  662.                 if (ships[s] != ufos[u]) {
  663.                     if (hopped[s] || PD[ships[s]][next[u]] > hopSize) continue;
  664.                     if (PD[ships[s]][next[u]] < bv) {
  665.                         bv = PD[ships[s]][next[u]];
  666.                         bs = s;
  667.                     }
  668.                 } else {
  669.                     bs = s;
  670.                     break;
  671.                 }
  672.             }
  673.            
  674.             if (bs == -1) continue;
  675.            
  676.             hopped[bs] = 1;
  677.             rv += PD[ships[bs]][next[u]] * (ships[bs] == ufos[u] ? UFO_MULT : 1.0);
  678.             planetsLeft -= vs[next[u]] == 0;
  679.             vs[next[u]] = 1;
  680.             ships[bs] = next[u];
  681.             swap(ships[firstShip], ships[bs]);
  682.             swap(hopped[firstShip], hopped[bs]);
  683.             firstShip++;
  684.         }
  685.        
  686.         REP(i, UNO) ufos[i] = next[i];
  687.         turn++;
  688.     }
  689.    
  690.     return rv + quickLB(vs, ships);
  691. }
  692.  
  693. VD monteCarlo(State &state, VC<Strategy> &str, VVI beam) {
  694.     VD rv(str.S);
  695.     VI left; REP(i, str.S) left.PB(i);
  696.     int totalSims = 0;
  697.     for (VI &beamLevel : beam) {
  698.         REP (i, beamLevel[0]) {
  699.             REP(u, UNO) rng.shuffle(state.ufoprob[u], state.turnsLeft + 5);
  700.             for (int j : left) rv[j] += monteCarlo(state, str[j]);
  701.         }
  702.         totalSims += beamLevel[0];
  703.         VC<pair<double, int>> vp;
  704.         for (int j : left) vp.PB(MP(rv[j], j));
  705.         sort(ALL(vp));
  706.        
  707.         left.clear();
  708.         REP(i, beamLevel[1]) left.PB(vp[i].Y);
  709.         FOR(i, beamLevel[1], vp.S) rv[vp[i].Y] /= totalSims;
  710.         REP(i, rv.S) rv[i] += 10000;
  711.     }
  712.    
  713.     return rv;
  714. }
  715.  
  716. double monteCarlo(State &state, Strategy &str, VVVI &prob, int from = -1, int to = -1) {
  717.     if (to == -1) from = 0, to = prob.S;
  718.     double rv = 0;
  719.     FOR(i, from, to) {
  720.         REP(u, UNO) state.ufoprob[u] = prob[i][u];
  721.         rv += monteCarlo(state, str);
  722.     }
  723.     return rv / (to - from);
  724. }
  725.  
  726. INLINE void updatePPOS(int s, int a, int b) {
  727.     FOR(i, a, b+1) PPOS[PATH[s][i]] = MP(s, i);
  728. }
  729.  
  730. bool performDetourStarted = false;
  731. double finalPathGain = 0;
  732.  
  733. void showFinalStats() {
  734.     cerr << "[Final Stats]" << endl;
  735.     cerr << "TestTime = " << totalTime << endl;
  736.     if (abs(finalPathGain) <= 1e-6) finalPathGain = 0;
  737.     DB(finalPathGain);
  738. }
  739.  
  740.  
  741.  
  742. class StarTraveller {public:
  743.  
  744. void passUfoW(VI &ufow) {
  745.     USEREALUFOW = true;
  746.     REP(i, ufow.S) REALUFOW[i] = ufow[i];
  747. }
  748.  
  749. int init(VI &stars) {
  750.     // REP(i, 123456) rng.rand();
  751.    
  752.     auto startTime = getTime();
  753.    
  754.     PNO = stars.S / 2;
  755.     REP(i, PNO) PP[i] = MP(stars[i*2], stars[i*2+1]);
  756.     REP(i, PNO) REP(j, PNO) PD[i][j] = sqrt((PP[i].X-PP[j].X)*(PP[i].X-PP[j].X) + (PP[i].Y-PP[j].Y)*(PP[i].Y-PP[j].Y));
  757.    
  758.     pair<double, int> vp[MAX_STARS];
  759.     REP(i, PNO) {
  760.         REP(j, PNO) vp[j] = MP(PD[i][j], j);
  761.         sort(vp, vp + PNO);
  762.         REP(j, PNO) PO[i][vp[j].Y] = j;
  763.         REP(j, PNO) PN[i][j] = vp[j].Y;
  764.     }
  765.    
  766.     memset(PATH, -1, sizeof(PATH));
  767.     memset(NEXTMOVE, -1, sizeof(NEXTMOVE));
  768.    
  769.     totalTime += getTime() - startTime;
  770.     DB(totalTime);
  771.     return 0;
  772. }
  773.  
  774. VI makeMoves(VI ufos, VI ships) {
  775.     auto startTime = getTime();
  776.    
  777.     turn++;
  778.     int turnsLeft = 4 * PNO - turn;
  779.    
  780.     if (FINAL_PATH) {
  781.         bool done = true;
  782.         REP(s, SNO) {
  783.             int bp = PATHTURN[s]+1<BPATHLEN[s] ? BPATH[s][PATHTURN[s]+1] : ships[s];
  784.             bool useufo = false;
  785.             bool shortcut = false;
  786.             double bv = PD[ships[s]][bp];
  787.             REP(u, UNO) {
  788.                 if (ufos[u*3] == ships[s]) {
  789.                     double av1 = BPATHLEN[s]-PATHTURN[s]-1<=turnsLeft ? UFO_MULT * PD[ships[s]][ufos[u*3+1]] + PD[ufos[u*3+1]][BPATH[s][PATHTURN[s]+1]] : 1e30;
  790.                     double av2 = BPATHLEN[s]-PATHTURN[s]+0<=turnsLeft ? UFO_MULT * (PD[ships[s]][ufos[u*3+1]] + PD[ufos[u*3+1]][ufos[u*3+2]]) + PD[ufos[u*3+2]][BPATH[s][PATHTURN[s]+1]] : 1e30;
  791.                     double av = min(av1, av2);
  792.                     if (av < bv) {
  793.                         if (av1 < av2) finalPathGain += bv - av;
  794.                         bv = av;
  795.                         shortcut = true;
  796.                         useufo = true;
  797.                         bp = ufos[u*3+1];
  798.                     }
  799.                 } else {
  800.                     double av = BPATHLEN[s]-PATHTURN[s]+0<=turnsLeft ? UFO_MULT * PD[ufos[u*3+1]][ufos[u*3+2]] + PD[ships[s]][ufos[u*3+1]] + PD[ufos[u*3+2]][BPATH[s][PATHTURN[s]+1]] : 1e30;
  801.                     if (av < bv) {
  802.                         bv = av;
  803.                         shortcut = true;
  804.                         useufo = false;
  805.                         bp = ufos[u*3+1];
  806.                     }
  807.                 }
  808.             }
  809.            
  810.             if (!shortcut && ufos.S && BPATHLEN[s]-PATHTURN[s]+1<=turnsLeft && PATHTURN[s]+1<BPATHLEN[s]) {
  811.                 if (rng.nextDouble() < PD[ships[s]][bp] * (turnsLeft-(BPATHLEN[s]-PATHTURN[s])) / (BPATHSUM[s][BPATHLEN[s]]-BPATHSUM[s][PATHTURN[s]]+1e-3)) {
  812.                     done = false;
  813.                     continue;
  814.                 }
  815.             }
  816.            
  817.             PATHTURN[s] += bp == BPATH[s][PATHTURN[s]+1];
  818.             done &= PATHTURN[s] >= BPATHLEN[s]-1;
  819.             finalPathGain -= PD[ships[s]][bp] * (useufo ? UFO_MULT : 1.0);
  820.             ships[s] = bp;
  821.         }
  822.        
  823.         totalTime += getTime() - startTime;
  824.         if (done) showFinalStats();
  825.         return ships;
  826.     }
  827.    
  828.     if (turn) for (int x : ships) vs[x] = 1;
  829.    
  830.    
  831.     SNO = ships.S;
  832.     UNO = ufos.S / 3;
  833.    
  834.     if (turn == 0)
  835.         REP(u, UNO) UFODIST[u] += PO[ufos[u*3]][ufos[u*3+1]];
  836.     REP(u, UNO) UFODIST[u] += PO[ufos[u*3+1]][ufos[u*3+2]];
  837.     if (turn == 0) {
  838.         str.transport = false;
  839.         str.hopSize = 0;
  840.         str.forceUfo = VI(UNO, 0);
  841.     }
  842.    
  843.    
  844.     bool allHopped = true;
  845.     REP(i, SNO) allHopped &= hopped[i] > 0;
  846.     VI ufoFollowed(UNO, 0);
  847.     VI shipFollowing(SNO, -1);
  848.     int ufosUsed = 0;
  849.     REP(u, UNO) REP(s, SNO) if (ufoFollowed[u] == 0 && shipFollowing[s] == -1 && ships[s] == ufos[3*u]) {
  850.         ufoFollowed[u] = 1;
  851.         shipFollowing[s] = u;
  852.         ufosUsed++;
  853.         break;
  854.     }
  855.    
  856.            
  857.     int planetsVisited = 0;
  858.     for (int x : vs) planetsVisited += x;
  859.     int planetsLeft = PNO - planetsVisited;
  860.    
  861.     bool MCHopCheck = false;
  862.     bool MCForceUfoCheck = false;
  863.     bool MCDisableUfoCheck = false;
  864.    
  865.     if (ufosUsed < UNO && (turn == 2 || turn == 25 || turn >= 75 && turn % 75 == 0)) {
  866.         MCHopCheck = !allHopped;
  867.         MCForceUfoCheck = true;
  868.         // MCDisableUfoCheck = true;
  869.     }
  870.    
  871.     bool hopSizeChange = false;
  872.     if ((MCHopCheck || MCForceUfoCheck || MCDisableUfoCheck)) DB(totalTime);
  873.     if ((MCHopCheck || MCForceUfoCheck || MCDisableUfoCheck) && totalTime < MAX_MC_TIME) {
  874.         State state;
  875.         state.turnsLeft = turnsLeft;
  876.         state.vs = VI(vs, vs + PNO);
  877.         state.ships = ships;
  878.         state.hopped = VI(hopped, hopped+SNO);
  879.        
  880.         VC<pair<int, int>> vp(UNO); REP(i, UNO) vp[i] = MP(max(10, min(PNO / 12, (int)round((turn + 2) * PNO / UFODIST[i]))), i);
  881.         sort(ALL(vp));
  882.         VI origUfo(UNO);
  883.         state.ufo = VVI(3, VI(UNO));
  884.         state.ufow = VI(UNO);
  885.         REP(j, 3) REP(i, UNO) state.ufo[j][i] = ufos[3*vp[i].Y+j];
  886.         REP(i, UNO) state.ufow[i] = vp[i].X;
  887.         REP(i, UNO) origUfo[i] = vp[i].Y;
  888.        
  889.         // if (turn <= 100 || turn - UPROBTURN > 600) {
  890.             REP(i, UNO) UPROB[i] = genDistribution(PNO, state.ufow[i], UFO_DIST_SIZE);
  891.             UPROBTURN = turn;
  892.         // }
  893.         state.ufoprob = VVI(UNO);
  894.         REP(i, UNO) state.ufoprob[i] = UPROB[i];
  895.        
  896.         const int MAX_SIMS = 80;
  897.        
  898.         auto mcstart = getTime();
  899.         DB(state.ufow);
  900.        
  901.         Strategy bestStr = str;
  902.         double bv = 1e30;
  903.        
  904.         str.forceUfo = VI(UNO, 0);
  905.         if (MCHopCheck) {
  906.             hopSizeChange = true;
  907.             VI hopSizes;
  908.             if (PNO < 1200)
  909.                 hopSizes = {0, 10, 20, 35, 50, 80, 125, 200, 300, 500, 750, 1500};
  910.             else
  911.                 hopSizes = {0, 20, 50, 125, 300, 750, 1500};
  912.             VC<Strategy> strs(hopSizes.S);
  913.             REP(i, hopSizes.S) {
  914.                 strs[i] = str;
  915.                 strs[i].transport = false;
  916.                 strs[i].hopSize = hopSizes[i];
  917.             }
  918.             VVI mcBeam;
  919.             if (FULL_EVALUATE) {
  920.                 mcBeam = {{50, 0}};
  921.             } else {
  922.                 mcBeam = {{20, 3}, {30, 0}};
  923.             }
  924.            
  925.             VD av = monteCarlo(state, strs, mcBeam);
  926.             REP(i, hopSizes.S)
  927.                 cerr << "Hop: " << hopSizes[i] << " Exp: " << av[i] << " Time Taken: " << getTime() - mcstart << endl;
  928.                
  929.             VC<pair<double, int>> vp;
  930.             REP(i, av.S) vp.PB(MP(av[i], i));
  931.             sort(ALL(vp));
  932.             bv = vp[0].X;
  933.             bestStr = strs[vp[0].Y];
  934.            
  935.             str = bestStr;
  936.             DB(str.hopSize);
  937.         }
  938.        
  939.         if (MCForceUfoCheck) {
  940.             const int SIM1_NO = 10;
  941.             const int SIM2_NO = 20;
  942.             VVVI prob(SIM2_NO, VVI(UNO));
  943.             REP(i, prob.S) REP(j, UNO) {
  944.                 prob[i][j] = state.ufoprob[j];
  945.                 rng.shuffle(prob[i][j], turnsLeft + 5);
  946.             }
  947.             double bv1 = monteCarlo(state, str, prob, 0, SIM1_NO);
  948.             double bv2 = monteCarlo(state, str, prob, SIM1_NO, SIM2_NO);
  949.             double bv = (bv1 + bv2) / 2;
  950.             Strategy mcstr = str;
  951.             bool noForce = true;
  952.             REP(i, UNO) if ((!noForce || !ufoFollowed[i]) && !str.forceUfo[i]) {
  953.                 mcstr.forceUfo[i] = 1;
  954.                 double av1 = monteCarlo(state, mcstr, prob, 0, SIM1_NO);
  955.                 if (av1 > bv1 && !FULL_EVALUATE) {
  956.                     mcstr.forceUfo[i] = 0;
  957.                     continue;
  958.                 }
  959.                 double av2 = monteCarlo(state, mcstr, prob, SIM1_NO, SIM2_NO);
  960.                 double av = (av1 + av2) / 2;
  961.                 if (av < bv - 100) {
  962.                     bv = av;
  963.                     bv1 = av1;
  964.                     str.forceUfo[origUfo[i]] = 1;
  965.                     noForce = false;
  966.                 } else {
  967.                     mcstr.forceUfo[i] = 0;
  968.                 }
  969.                 cerr << "ForceUfo: " << i << " Exp: " << av << " Time Taken: " << getTime() - mcstart << endl;
  970.             }
  971.         }      
  972.        
  973.         totalTime += getTime() - startTime;
  974.         startTime = getTime();
  975.     }
  976.    
  977.     REP(u, UNO) REP(s, SNO) if (ufoFollowed[u] == 0 && shipFollowing[s] == -1 && ships[s] == ufos[3*u+1]) {
  978.         ufoFollowed[u] = 1;
  979.         shipFollowing[s] = u;
  980.         ufosUsed++;
  981.         break;
  982.     }
  983.     REP(u, UNO) REP(s, SNO) if (ufoFollowed[u] == 0 && shipFollowing[s] == -1 && ships[s] == ufos[3*u+2]) {
  984.         ufoFollowed[u] = 1;
  985.         shipFollowing[s] = u;
  986.         ufosUsed++;
  987.         break;
  988.     }
  989.    
  990.     VI rv = VI(SNO, -1);
  991.     if (UNO == 0 || planetsLeft >= turnsLeft - 2) {
  992. #ifndef LOCAL
  993.         double oldTime = totalTime;
  994.         totalTime = 20.0 - Server::getRemainingTime(0) * 1e-3;
  995.         cerr << "Old: " << oldTime << " vs New: " << totalTime << endl;
  996. #endif
  997.  
  998.         VI notVisited;
  999.         REP(p, PNO) if (!vs[p]) notVisited.PB(p);
  1000.         REP(i, notVisited.S) swap(notVisited[i], notVisited[rng.next(i, notVisited.S)]);
  1001.        
  1002.         DB(greedy(ships, notVisited)); 
  1003.        
  1004.         REP(s, SNO) {
  1005.             PATHLEN[s] = 1;
  1006.             PATH[s][0] = ships[s];
  1007.         }
  1008.         REP(i, notVisited.S) PATH[i%SNO][PATHLEN[i%SNO]++] = notVisited[i];
  1009.         REP(i, SNO) PATH[i][PATHLEN[i]] = PNO;
  1010.         REP(s, SNO) {
  1011.             BPATHLEN[s] = PATHLEN[s];
  1012.             REP(i, PATHLEN[s]+1) BPATH[s][i] = PATH[s][i];
  1013.         }
  1014.        
  1015.         REP(i, PNO) REP(j, PNO) OPD[i][j] = PD[i][j];
  1016.        
  1017.         REP(s, SNO) vs[ships[s]] = 0;
  1018.         REP(i, PNO) {
  1019.             int pos = 0;
  1020.             REP(j, PNO) if (PN[i][j] != i && !vs[PN[i][j]]) {
  1021.                 PC[i][pos++] = PN[i][j];
  1022.                 if (pos == MAX_NEIGHBOURS) break;
  1023.             }
  1024.         }
  1025.         REP(s, SNO) vs[ships[s]] = 1;
  1026.        
  1027.         REP(s, SNO) updatePPOS(s, 0, PATHLEN[s]-1);
  1028.        
  1029.         REP(u, UNO) {
  1030.             int p = ufos[3*u];
  1031.             if (!vs[p]) continue;
  1032.             REP(i, PNO) {
  1033.                 double d = OPD[p][i];
  1034.                 d = min(d, OPD[ufos[3*u+1]][i]);
  1035.                 d = min(d, OPD[ufos[3*u+2]][i]);
  1036.                 PD[p][i] = PD[i][p] = d;
  1037.             }
  1038.         }
  1039.        
  1040.         cerr << "[Mode] TSP at: " << totalTime << endl;
  1041.         DB(notVisited.S);
  1042.        
  1043.         double bv = calcSolution();
  1044.         double xv = bv;
  1045.         VD shipCost(SNO);
  1046.         REP(i, SNO) shipCost[i] = calcSolution(i);
  1047.        
  1048.         int steps = 0;
  1049.         int steps2 = 0;
  1050.         int acc = 0;
  1051.         int bestacc = 0;
  1052.         double timePassed = 0.0;
  1053.         double temp = 0.0;
  1054.        
  1055.         int NEIGHBOURS = min(MAX_NEIGHBOURS, (int)notVisited.S);
  1056.         while (true) {
  1057.             if ((steps & 1023) == 0) {
  1058.                 if (notVisited.S == 1) break;
  1059.                 timePassed = (getTime() - startTime) / (TIME_LIMIT - totalTime);
  1060.                 temp = (1.0 - pow(timePassed, 0.15)) * 160;//(notVisited.S < 50 ? 320 : notVisited.S < 100 ? 240 : 160);
  1061.                 if (timePassed >= 1.0) break;
  1062.             }
  1063.             steps++;
  1064.            
  1065.            
  1066.             int s0 = rng.next(SNO);
  1067.             if (PATHLEN[s0] == 1) continue;
  1068.             int p0 = rng.next(1, PATHLEN[s0]);
  1069.             int mode = rng.rand() & 1;
  1070.            
  1071.             int p1, s1;
  1072.             if (mode == 0) {
  1073.                 PII p = PPOS[PC[PATH[s0][p0]][rng.next(NEIGHBOURS)]];
  1074.                 s1 = p.X;
  1075.                 p1 = p.Y;
  1076.             } else {
  1077.                 p1 = 0;
  1078.                 s1 = rng.next(SNO);
  1079.             }
  1080.             int p0b, p1b;
  1081.             int type;
  1082.            
  1083.             double diff = 0;
  1084.             if (s0 == s1) {
  1085.                 if (p1 == 0) p1 = rng.next(1, PATHLEN[s0]);
  1086.                 type = rng.next(2) == 0;
  1087.                 if (type == 0) {
  1088.                     if (p0 == p1) continue;
  1089.                     if (p0 > p1) swap(p0, p1);
  1090.                     diff -= pathDist(s0, p0-1, PATH[s0][p0]);
  1091.                     diff -= pathDist(s0, p1+1, PATH[s0][p1]);
  1092.                     diff += pathDist(s0, p0-1, PATH[s0][p1]);
  1093.                     diff += pathDist(s0, p1+1, PATH[s0][p0]);
  1094.                 } else {
  1095.                     if (abs(p0-p1)<=2) continue;
  1096.                     diff -= pathDist(s0, p0-1, PATH[s0][p0]);
  1097.                     diff -= pathDist(s0, p0+1, PATH[s0][p0]);
  1098.                     diff += pathDist(s0, p0+1, PATH[s0][p0-1]);
  1099.                     diff -= pathDist(s0, p1, PATH[s0][p1-1]);
  1100.                     diff += pathDist(s0, p1-1, PATH[s0][p0]);
  1101.                     diff += pathDist(s0, p1, PATH[s0][p0]);
  1102.                 }
  1103.             } else {
  1104.                 type = rng.next(2) == 0;
  1105.                 if (type == 0) {
  1106.                     if (PATHLEN[s1] == 1) continue;
  1107.                     if (p1 == 0) p1 = rng.next(1, PATHLEN[s1]);
  1108.                     p0b = rng.next(1, PATHLEN[s0]);
  1109.                     p1b = rng.next(1, PATHLEN[s1]);
  1110.                     if (p0 > p0b) swap(p0, p0b);
  1111.                     if (p1 > p1b) swap(p1, p1b);
  1112.                     // int move = (p0b-p0)-(p1b-p1);
  1113.                     // if (PATHLEN[s0]-move>turnsLeft || PATHLEN[s1]+move>turnsLeft) continue;
  1114.                     diff -= pathDist(s0, p0-1, PATH[s0][p0]);
  1115.                     diff -= pathDist(s0, p0b+1, PATH[s0][p0b]);
  1116.                     diff += pathDist(s0, p0-1, PATH[s1][p1]);
  1117.                     diff += pathDist(s0, p0b+1, PATH[s1][p1b]);
  1118.                     diff -= pathDist(s1, p1-1, PATH[s1][p1]);
  1119.                     diff -= pathDist(s1, p1b+1, PATH[s1][p1b]);
  1120.                     diff += pathDist(s1, p1-1, PATH[s0][p0]);
  1121.                     diff += pathDist(s1, p1b+1, PATH[s0][p0b]);
  1122.                 } else {
  1123.                     p1 = p1 == 0 ? rng.next(1, PATHLEN[s1]+1) : p1 + 1;
  1124.                     if (rng.next(2)) {
  1125.                         p0b = rng.next(1, PATHLEN[s0]);
  1126.                         if (p0 > p0b) swap(p0, p0b);
  1127.                     } else {
  1128.                         p0b = p0;
  1129.                     }
  1130.                     p1b = p1-1;
  1131.                    
  1132.                     // if (PATHLEN[s1]+(p0b-p0)+1>turnsLeft) continue;
  1133.                     diff -= pathDist(s0, p0-1, PATH[s0][p0]);
  1134.                     diff -= pathDist(s0, p0b+1, PATH[s0][p0b]);
  1135.                     diff += pathDist(s0, p0b+1, PATH[s0][p0-1]);
  1136.                     diff -= pathDist(s1, p1, PATH[s1][p1-1]);
  1137.                     diff += pathDist(s1, p1-1, PATH[s0][p0]);
  1138.                     diff += pathDist(s1, p1, PATH[s0][p0b]);
  1139.                 }
  1140.             }
  1141.            
  1142.             if (diff <= rng.nextDouble() * temp) {
  1143.                 if (s0 == s1) {
  1144.                     if (type == 0) {
  1145.                         reverse(PATH[s0]+p0, PATH[s0]+p1+1);
  1146.                         updatePPOS(s0, p0, p1);
  1147.                     } else {
  1148.                         int tmp = PATH[s0][p0];
  1149.                         if (p0 < p1) {
  1150.                             FOR(i, p0, p1-1) PATH[s0][i] = PATH[s0][i+1];
  1151.                             PATH[s0][p1-1] = tmp;
  1152.                         } else {
  1153.                             for (int i = p0; i > p1; i--) PATH[s0][i] = PATH[s0][i-1];
  1154.                             PATH[s0][p1] = tmp;
  1155.                         }
  1156.                         updatePPOS(s0, p0, p1);
  1157.                     }
  1158.                 } else {
  1159.                     if (p0b-p0>p1b-p1) {
  1160.                         swap(s0, s1);
  1161.                         swap(p0, p1);
  1162.                         swap(p0b, p1b);
  1163.                     }
  1164.                     int move = (p1b-p1)-(p0b-p0);
  1165.                     assert(move >= 0);
  1166.                     FOR(i, p0, p0b+1) TPATH[i] = PATH[s0][i];
  1167.                     PATHLEN[s0] += move;
  1168.                     PATHLEN[s1] -= move;
  1169.                     assert(PATHLEN[s0] <= turnsLeft+1);
  1170.                     assert(PATHLEN[s1] <= turnsLeft+1);
  1171.                     for (int i = PATHLEN[s0]-1; i >= p0b+1+move; i--) PATH[s0][i] = PATH[s0][i-move]; //TODO: fix
  1172.                     FOR(i, p1, p1b+1) PATH[s0][i-p1+p0] = PATH[s1][i];
  1173.                     FOR(i, p0, p0b+1) PATH[s1][i-p0+p1] = TPATH[i];
  1174.                     FOR(i, p1b+1-move, PATHLEN[s1]) PATH[s1][i] = PATH[s1][i+move];
  1175.                     updatePPOS(s0, p0, PATHLEN[s0]-1);
  1176.                     updatePPOS(s1, 1, PATHLEN[s1]-1);
  1177.                 }
  1178.                 PATH[s0][PATHLEN[s0]] = PNO;
  1179.                 PATH[s1][PATHLEN[s1]] = PNO;
  1180.                 bv += diff;
  1181.                 acc++;
  1182.                 if (bv < xv) {
  1183.                     xv = bv;
  1184.                     bestacc++;
  1185.                     REP(s, SNO) {
  1186.                         BPATHLEN[s] = PATHLEN[s];
  1187.                         REP(i, PATHLEN[s]) BPATH[s][i] = PATH[s][i];
  1188.                     }
  1189.                 }
  1190.             }
  1191.         }
  1192.         REP(i, PNO) REP(j, PNO) PD[i][j] = OPD[i][j];
  1193.        
  1194.         DB(calcSolution());
  1195.         DB(xv);
  1196.         DB(steps);
  1197.         DB(acc);
  1198.         DB(bestacc);
  1199.         DB(detours);
  1200.         FINAL_PATH = true;
  1201.         REP(s, SNO) {
  1202.             BPATHSUM[s][0] = 0;
  1203.             FOR(i, 1, BPATHLEN[s]) {
  1204.                 BPATHSUM[s][i] = BPATHSUM[s][i-1] + PD[BPATH[s][i-1]][BPATH[s][i]];
  1205.             }
  1206.         }
  1207.         finalPathGain = xv;
  1208.         totalTime += getTime() - startTime;
  1209.         DB(totalTime);
  1210.         return makeMoves(ufos, ships);
  1211.     } else {
  1212.         //Wait for UFO
  1213.         bool performDetour = false;
  1214.         double mcmstEstimate = 0;
  1215.         VI upos, fpos;
  1216.         if (planetsLeft + 1200 > turnsLeft && planetsLeft < 500) {
  1217.             if (!performDetourStarted) {
  1218.                 REP(i, UNO) UPROB[i] = genDistribution(PNO, max(12, min(PNO / 12, (int)round((turn + 2) * PNO / UFODIST[i]))), UFO_DIST_SIZE);
  1219.                 performDetourStarted = true;
  1220.                 cerr << "[Mode] Detour at: " << totalTime << endl;
  1221.             }
  1222.             UPROBTURN = turn;
  1223.             performDetour = true;
  1224.             ZERO(XVS);
  1225.             XVSNO = 1;
  1226.            
  1227.            
  1228.            
  1229.             REP(u, UNO) upos.PB(ufos[3*u+2]);
  1230.             REP(s, SNO) if (shipFollowing[s] == -1) fpos.PB(ships[s]);
  1231.            
  1232.             if (totalTime < MAX_MCMST_TIME && (planetsLeft < MCMST_MAX_PLANETS && turnsLeft - planetsLeft < MCMST_MAX_TURNS)) {
  1233.                 mcmstEstimate = calcMCMST(fpos, upos, ufoFollowed, turnsLeft - planetsLeft, MCMST_SIMS);
  1234.             }
  1235.             PEXPCOST = max(50.0, 90 - (turnsLeft - planetsLeft) * 0.20);
  1236.         }
  1237.        
  1238.         VI usedUfo(UNO);
  1239.         rv = VI(SNO, -1);
  1240.         int ufosUsed = 0;
  1241.  
  1242.         if (hopSizeChange) {
  1243.             DB(VI(hopped, hopped+SNO));
  1244.             DB(str.hopSize);
  1245.             VVI w(SNO, VI(UNO, -1));
  1246.             REP(s, SNO) REP(u, UNO) {
  1247.                 double minDist = 1e30;
  1248.                 if (ships[s] == ufos[3*u]) minDist = 0;
  1249.                 minDist = min(minDist, PD[ships[s]][ufos[3*u+1]]);
  1250.                 minDist = min(minDist, PD[ships[s]][ufos[3*u+2]]);
  1251.                 if (minDist > str.hopSize && !str.forceUfo[u]) continue;
  1252.                 if (hopped[s] && !str.forceUfo[u] && minDist > 0) continue;
  1253.                 w[s][u] = (int)round(minDist) + (str.forceUfo[u] ? 0 : 2000);
  1254.             }
  1255.             VI matching = bestMatching(w);
  1256.             REP(s, SNO) if (matching[s] != -1) {
  1257.                 int u = matching[s];
  1258.                 int p = ufos[3*u+1];
  1259.                 if (ships[s] != ufos[3*u] && PD[ships[s]][ufos[3*u+2]] < PD[ships[s]][ufos[3*u+1]]) {
  1260.                     // cerr << "Saved: " << PD[ships[s]][ufos[3*u+1]] - PD[ships[s]][ufos[3*u+2]] << endl;
  1261.                     p = ships[s];
  1262.                     NEXTMOVE[s] = u;
  1263.                 }
  1264.                 hopped[s] = 1;
  1265.                 rv[s] = p;
  1266.                 vs[p] = 1;
  1267.             }
  1268.             DB(VI(hopped, hopped+SNO));
  1269.             DB(matching);
  1270.         } else {
  1271.             if (mcmstEstimate && turn % 10 == 0) {
  1272.                 DB(ufoFollowed);
  1273.                 DB(fpos);
  1274.                 cerr << "MCMST Calculate Jump" << endl;
  1275.                 REP(u, UNO) if (!ufoFollowed[u]) {
  1276.                     int s = rng.next(SNO);
  1277.                     // REP(i, SNO) if (PD[ships[i]][ufos[3*u+1]] < PD[ships[s]][ufos[3*u+1]]) s = i;
  1278.                     if (rv[s] != -1) continue;
  1279.                     if (shipFollowing[s] == -1) {
  1280.                         bool found = false;
  1281.                         REP(i, fpos.S) if (fpos[i] == ships[s]) {
  1282.                             fpos.erase(fpos.begin() + i);
  1283.                             found = true;
  1284.                             break;
  1285.                         }
  1286.                         assert(found);
  1287.                     } else {
  1288.                         ufoFollowed[shipFollowing[s]] = 0;
  1289.                     }
  1290.                     ufoFollowed[u] = 1;
  1291.                     double estimate = calcMCMST(fpos, upos, ufoFollowed, turnsLeft - planetsLeft, MCMST_SIMS);
  1292.                     ufoFollowed[u] = 0;
  1293.                     if (shipFollowing[s] == -1) {
  1294.                         fpos.PB(ships[s]);
  1295.                     } else {
  1296.                         ufoFollowed[shipFollowing[s]] = 1;
  1297.                     }
  1298.                     double dist = min(PD[ships[s]][ufos[3*u+1]], PD[ships[s]][ufos[3*u+2]]);
  1299.                     if (estimate + dist < mcmstEstimate - 200) {
  1300.                         ufosUsed += usedUfo[u] == 0;
  1301.                         usedUfo[u] = 1;
  1302.                         rv[s] = PD[ships[s]][ufos[3*u+1]] < PD[ships[s]][ufos[3*u+2]] ? ufos[3*u+1] : ufos[3*u+2];
  1303.                         vs[rv[s]] = 1;
  1304.                         hopped[s] = 1;
  1305.                         if (shipFollowing[s] == -1) {
  1306.                             bool found = false;
  1307.                             REP(i, fpos.S) if (fpos[i] == ships[s]) {
  1308.                                 fpos.erase(fpos.begin() + i);
  1309.                                 found = true;
  1310.                                 break;
  1311.                             }
  1312.                             assert(found);
  1313.                         } else {
  1314.                             ufoFollowed[shipFollowing[s]] = 0;
  1315.                         }
  1316.                         ufoFollowed[u] = 1;
  1317.                         shipFollowing[s] = u;
  1318.                         cerr << "Performed Jump: " << mcmstEstimate << " -> " << estimate << " + " << dist << " = " << estimate + dist << endl;
  1319.                     }
  1320.                 }
  1321.                 // cerr << "Finished Jumping" << endl;
  1322.             }
  1323.            
  1324.             while (true) {
  1325.                 double bv = -1e30;
  1326.                 int bs = -1;
  1327.                 int bp = -1;
  1328.                 int bu = -1;
  1329.                 int bt = -1;
  1330.                 REP(s, SNO) if (rv[s] == -1) {
  1331.                     REP(u, UNO) if (usedUfo[u] == 0 || false && performDetour && SNO > UNO) {
  1332.                         int p1 = ufos[u*3+0];
  1333.                         int p2 = ufos[u*3+1];
  1334.                         int p3 = ufos[u*3+2];
  1335.                         double av = usedUfo[u] * -1e10;
  1336.                         bool sameplanet = p1 == ships[s] || p2 == ships[s] || p3 == ships[s];
  1337.                         // if (!sameplanet && usedUfo[u]) continue;
  1338.                         int tp = p1 == ships[s] ? p2 : ships[s];
  1339.                         if (!sameplanet) {
  1340.                             double d1 = PD[ships[s]][ufos[u*3+1]];
  1341.                             double d2 = PD[ships[s]][ufos[u*3+2]];
  1342.                             if (min(d1, d2) > (hopped[s] ? 0 : str.hopSize)) continue;
  1343.                             if (d1 < d2) {
  1344.                                 tp = p2;
  1345.                                 av += d1 * -1e6;
  1346.                             } else {
  1347.                                 tp = -1;
  1348.                                 av += d2 * -1e6;
  1349.                             }
  1350.                         }
  1351.                         av += UFODIST[u] / (turn+2) * 1e3;
  1352.                         av += (p1 == ships[s]);
  1353.                         if (av > bv) {
  1354.                             bv = av;
  1355.                             bu = u;
  1356.                             bs = s;
  1357.                             bp = tp;
  1358.                             bt = 0;
  1359.                         }
  1360.                     }
  1361.                    
  1362.                     if (NEXTMOVE[s] != -1) {
  1363.                         bu = NEXTMOVE[s];
  1364.                         bs = s;
  1365.                         bt = 0;
  1366.                         if (ships[bs] == ufos[3*bu] || PD[ships[bs]][ufos[3*bu+1]] < PD[ships[bs]][ufos[3*bu+2]]) {
  1367.                             bp = ufos[3*bu+1];
  1368.                             NEXTMOVE[s] = -1;
  1369.                         } else {
  1370.                             bp = ships[bs];
  1371.                         }
  1372.                         break;
  1373.                     }
  1374.                 }
  1375.                
  1376.                 if (bu == -1) break;
  1377.                
  1378.                 if (performDetour && ufos[bu*3] == ships[bs] && vs[ufos[bu*3+1]] && vs[ufos[bu*3+2]]) {
  1379.                     bv = 1e30;
  1380.                     int xp = -1;
  1381.                     REP(i, PNO) if (!vs[i]) {
  1382.                         double av = min(PD[i][ships[bs]], min(PD[i][ufos[bu*3+1]], PD[i][ufos[bu*3+2]]));
  1383.                         if (av < bv) {
  1384.                             bv = av;
  1385.                             xp = i;
  1386.                         }
  1387.                     }
  1388.                    
  1389.                     if (xp != -1) {
  1390.                         if (mcmstEstimate && ufoFollowed[bu]) {
  1391.                             ufoFollowed[bu] = 0;
  1392.                             fpos.PB(xp);
  1393.                             double estimate = calcMCMST(fpos, upos, ufoFollowed, turnsLeft - planetsLeft, MCMST_SIMS);
  1394.                             fpos.pop_back();
  1395.                             ufoFollowed[bu] = 1;
  1396.                             DB(totalTime);
  1397.                             DB(turn);
  1398.                             DB(bu);
  1399.                             DB(PEXPCOST);
  1400.                             DB(bv);
  1401.                             DB(estimate);
  1402.                             DB(mcmstEstimate);
  1403.                             if (estimate + bv > mcmstEstimate - (15 + (turnsLeft - planetsLeft) * 0.5)) {
  1404.                                 xp = -1;
  1405.                             }
  1406.                         } else if (bv < PEXPCOST) {
  1407.                             int turns = turnsLeft - planetsLeft;
  1408.                             const int sims = 100;
  1409.                             int vis1 = 0;
  1410.                             int visx = 0;
  1411.                             FOR(sim, 1, sims+1) {
  1412.                                 XVSNO++;
  1413.                                 int planet = ufos[bu*3+2];
  1414.                                 int no = 0;
  1415.                                 REP(i, turns) {
  1416.                                     planet = PN[planet][UPROB[bu][rng.nextAnd(UFO_DIST_SIZE-1)]];
  1417.                                     if (!vs[planet] && XVS[planet] < XVSNO) {
  1418.                                         XVS[planet] = XVSNO;
  1419.                                         no++;
  1420.                                     }
  1421.                                 }
  1422.                                 vis1 += no > 0;
  1423.                                 visx += no;
  1424.                                 if (sim == 20 && (vis1 >= 20 || vis1 <= 4) || sim == 50 && (vis1 >= 48 || vis1 <= 20)) {
  1425.                                     vis1 = vis1 > sim / 2 ? sims : 0;
  1426.                                     break;
  1427.                                 }
  1428.                             }
  1429.                             if (vis1 > 0.9 * sims || visx > 1.4 * sims) xp = -1;
  1430.                         } else {
  1431.                             xp = -1;
  1432.                         }
  1433.                     }
  1434.                     if (xp != -1 && PD[xp][ships[bs]] <= min(PD[xp][ufos[bu*3+1]], PD[xp][ufos[bu*3+2]])) {
  1435.                         bp = xp;
  1436.                         bt = 1;
  1437.                         DB(bv);
  1438.                     }
  1439.                 }
  1440.                
  1441.                 // detours += bt == 1;
  1442.                 bp = bt == 0 && ufos[bu*3] == ships[bs] && vs[ufos[bu*3+1]] && ufos[bu*3+2] == ships[bs] ? ships[bs] : bp;
  1443.                 if (bp != -1) {
  1444.                     hopped[bs] = 1;
  1445.                     rv[bs] = bp;
  1446.                     vs[bp] = 1;
  1447.                 }
  1448.                 ufosUsed += usedUfo[bu] == 0;
  1449.                 usedUfo[bu]++;
  1450.                 if (ufosUsed == UNO && !performDetour) break;
  1451.             }
  1452.         }
  1453.         REP(i, SNO) if (rv[i] == -1) rv[i] = ships[i];
  1454.        
  1455.         if (false && turn % 1000 == 0) {
  1456.             DB(turn);
  1457.             DB(hops);
  1458.             REP(i, UNO)
  1459.                 cerr << i << ' ' << usedUfo[i] << ' ' << UFODIST[i]*1.0/(turn+2) << endl;
  1460.         }
  1461.     }
  1462.    
  1463.     bool allVS = true;
  1464.     REP(i, PNO) if (!vs[i]) {
  1465.         allVS = false;
  1466.         break;
  1467.     }
  1468.    
  1469.     totalTime += getTime() - startTime;
  1470.     if (allVS) showFinalStats();
  1471.     return rv;
  1472. }
  1473.  
  1474. };
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement