Advertisement
Guest User

Untitled

a guest
Feb 23rd, 2019
104
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 12.10 KB | None | 0 0
  1. #ifdef LOCAL
  2. //#define _GLIBCXX_DEBUG
  3. // #pragma comment(linker, "/STACK:102400000,102400000") // 在 Windows 上有效
  4. #endif
  5. #pragma GCC optimize("Ofast") // Ofast 等效于 -O3 -ffast-math
  6. #include <bits/stdc++.h>
  7.  
  8. using namespace std;
  9.  
  10. #define FAST_READ ios::sync_with_stdio(false); cin.tie(nullptr);
  11. #ifdef LOCAL
  12. #define see(x) cout <<"<DBG> " << #x << ": " << (x) << endl
  13. #endif
  14. #ifndef LOCAL
  15. #define see(x)
  16. #endif
  17.  
  18. template<typename T> void dprintln(const T &t) { cout << t << endl; } // for debug use
  19. template<typename T, typename ...Args> void dprintln(const T &t, const Args &...rest) { cout << t << ' '; dprintln(rest...); }
  20. template<typename T> void println(const T &t) { cout << t << '\n'; }
  21. template<typename T, typename ...Args> void println(const T &t, const Args &...rest) { cout << t << ' '; println(rest...); }
  22.  
  23. template<typename T> void print(const T &t) { cout << t << ' '; }
  24.  
  25. template<typename T, typename ...Args> void print(const T &t, const Args &...rest) { cout << t << ' '; print(rest...); }
  26.  
  27. // this overload is chosen when there's only one argument
  28. template<class T> void scan(T &t) { cin >> t; }
  29. template<class T, class ...Args> void scan(T &a, Args &...rest) { cin >> a; scan(rest...); }
  30. template <typename T> void read(vector<T> &vec) {for(auto &x: vec) cin >> x;}
  31. using ull = unsigned long long;
  32. using ll = long long;
  33. using ul = unsigned long;
  34. using vl = vector<ll>;
  35. using vi = vector<int>;
  36. using pii = pair<int,int>;
  37. using vb = vector<bool>;
  38. using vpii = vector<pii>;
  39. template <typename int_t>
  40. inline int_t lowbit(int_t x) {return x & -x;}
  41. #define rng(i, a, b) for(int i = (a); i < (b); ++i)
  42. #define down(i, b, a) for (int i = (b); i >= (a); i--)
  43. #define rep(n) for(int _ = 0, __ = (int)n; _ < __; _++)
  44. #define stp(i, a, b, c) for (int i = (a); i < (b); i += (c))
  45. #define FOR(x, cont) for (const auto &x: cont)
  46. #define For(x, cont) for (auto &x: cont)
  47. #define all(x) begin(x), end(x)
  48. #define pb push_back
  49. #define mp make_pair
  50. #define eb emplace_back
  51. #define SZ(x) (int)(x).size()
  52. template <typename T, typename Comp = less<T>>
  53. using pq = priority_queue<T, vector<T>, Comp>;
  54. #define popcnt(x) __builtin_popcountll((x))
  55. #define SET(arr, v) memset(arr, (v), sizeof (arr))
  56. #define UNIQ(vec) (vec).erase(unique(all(vec)), end(vec))
  57. #define LB(cont, x) int(lower_bound(all(cont), x) - begin(cont))
  58. #define AI(name, n, m) vv<int> name(n, vi(m));
  59. #define AL(name, n, m) vv<ll> name(size_t(n), vl(size_t(m)));
  60. auto bet = [](const ll x, const ll y, const ll i) {
  61. return x <= i && i <= y;
  62. };
  63.  
  64. template<typename T1, typename T2>
  65. T1 ceil(T1 x, T2 y) { // y >= 1,是整数。需要注意 x + y - 1 是否会溢出
  66. return (x + y - 1) / y;
  67. }
  68.  
  69.  
  70. inline int h_bit(unsigned long long x) {
  71. return sizeof(unsigned long long) * 8 - 1 - __builtin_clzll(x);
  72. }
  73. int pow2(int x){ // power of 2
  74. return 1 << (h_bit((ull)x) + (x != lowbit(x)));
  75. }
  76.  
  77. template <typename T>
  78. struct bit {
  79. vector<T> a;
  80. function<T(T,T)> bin_op;
  81. explicit bit(int n, function<T(T,T)> bin_op = [](const T &a, const T &b){return max(a, b);}, int v = 0):bin_op(bin_op){
  82. a.resize((ul)n + 1);
  83. for (int i = 1; i <= n; ++i) a[i] = v;
  84. }
  85.  
  86. T prefix(T x) {
  87. T res = 0;
  88. while (x) {
  89. res = bin_op(res, a[x]);
  90. // res += a[x];
  91. x -= x & -x;
  92. }
  93. return res;
  94. }
  95.  
  96. T sum(int l, int r) {
  97. if (l > r) return 0;
  98. return sum(r) - sum(l - 1);
  99. }
  100.  
  101. void modify(int x, T v) {
  102. while (x < a.size()) {
  103. a[x] = bin_op(a[x], v);
  104. x += x & -x;
  105. }
  106. }
  107. void clear(){
  108. fill(a.begin(), a.end(), 0);
  109. }
  110. };
  111.  
  112. template <typename T>
  113. struct r_bit {
  114. vector<T> a;
  115. function<T(T,T)> bin_op;
  116. explicit r_bit(int n, function<T(T,T)> bin_op = [](const T &a, const T &b){return max(a, b);}, int v = 0):bin_op(bin_op){
  117. a.resize((ul)n + 1);
  118. for (int i = 1; i <= n; ++i) a[i] = v;
  119. }
  120.  
  121. T suffix(T x) {
  122. T res = 0;
  123. while (x < a.size()) {
  124. res = bin_op(res, a[x]);
  125. // res += a[x];
  126. x += x & -x;
  127. }
  128. return res;
  129. }
  130.  
  131. T sum(int l, int r) {
  132. if (l > r) return 0;
  133. return sum(r) - sum(l - 1);
  134. }
  135.  
  136. void modify(int x, T v) {
  137. while (x > 0) {
  138. a[x] = bin_op(a[x], v);
  139. x -= x & -x;
  140. }
  141. }
  142. void clear(){
  143. fill(a.begin(), a.end(), 0);
  144. }
  145. };
  146.  
  147.  
  148.  
  149. vi get_prime(int n) {
  150. vi minp((ul)n + 1), p;
  151. for (int i = 2; i <= n; i++) {
  152. if (!minp[i]) {
  153. minp[i] = i;
  154. p.pb(i);
  155. }
  156. FOR(x, p) {
  157. if (x <= minp[i] && x * i <= n)
  158. minp[x * i] = x;
  159. else break;
  160. }
  161. }
  162. return p;
  163. }
  164.  
  165. const int mod = 998244353; //记住 mod 的行号!
  166.  
  167. inline void mul(int &x, const int y) {
  168. x = int(ll(x) * y % mod);
  169. }
  170. inline int Mul(int x, const int y) {
  171. return (int)(ll(x) * y % mod);
  172. }
  173. inline void add(int &x, const int y) {
  174. x += y;
  175. if (x >= mod) x -= mod;
  176. }
  177. inline int Add(int x, int y) {
  178. x += y;
  179. if (x >= mod) x -= mod;
  180. return x;
  181. }
  182. inline int Sub(int x, int y) {
  183. x -= y; if (x < 0) x += mod;
  184. return x;
  185. }
  186. inline ll submod(ll x, ll y) {
  187. return x >= y ? x - y : x - y + mod;
  188. }
  189. inline ll addmod(ll x, ll y) {
  190. x += y;
  191. return x >= mod ? x - mod : x;
  192. }
  193.  
  194. inline void sub_mod(ll &x, const ll y) {
  195. x -= y;
  196. if (x < 0) x += mod;
  197. }
  198. // alias templates
  199. template<typename T> using vv = vector<vector<T>>;
  200. template <typename T1, typename T2 = T1> using vp = vector<pair<T1,T2>>;
  201. template<typename T, int n> using va = vector<array<T,n>>;
  202.  
  203. #ifdef __GNUC__
  204. #include <ext/pb_ds/assoc_container.hpp>
  205. #include <ext/pb_ds/tree_policy.hpp>
  206. template <typename T>
  207. using rank_tree = __gnu_pbds::tree<
  208. T,
  209. __gnu_pbds::null_type,
  210. less<T>,
  211. __gnu_pbds::rb_tree_tag,
  212. __gnu_pbds::tree_order_statistics_node_update>;
  213. #endif
  214.  
  215. //union-find 并查集
  216. struct UF{
  217. vi par, sz;
  218. int tot;
  219. explicit UF(int n) { //0-indexed
  220. par.assign(n, 0);
  221. sz.assign(n, 1);
  222. rng (i, 0, n) par[i] = i;
  223. tot = n;
  224. }
  225.  
  226. int find(int x) {
  227. return x == par[x] ? x : par[x] = find(par[x]);
  228. }
  229. bool unite(int x, int y) {
  230. int rx = find(x), ry = find(y);
  231. if(rx != ry) {
  232. par[rx] = ry;
  233. --tot;
  234. return true;
  235. }
  236. return false;
  237. }
  238. };
  239.  
  240. int fp (int x, ll n) { //fast power: hat off to quailty
  241. if(n < 0) {
  242. x = fp(x, mod - 2);
  243. n = -n;
  244. }
  245. int ans = 1;
  246. while(n){
  247. if(n&1) mul(ans, x);
  248. n /= 2;
  249. mul(x, x);
  250. }
  251. return ans;
  252. }
  253.  
  254.  
  255. vi get_popcnt(int n) {
  256. vi res((ul)n + 1);
  257. rng (i, 0, n) {
  258. if (i * 2 <= n) res[i * 2] = res[i];
  259. if ((i & 1) == 0) res[i + 1] = res[i] + 1;
  260. }
  261. return res;
  262. }
  263.  
  264.  
  265. using idx_t = long;
  266. template <typename T, typename Compare = std::less<T>>
  267. struct ST{
  268. size_t n; // 0-indexed
  269. vv<T> a;
  270. template <typename ptr_t>
  271. ST(ptr_t beg, ptr_t end):n(end-beg){
  272. a.resize((size_t)h_bit(n)+1); // 注意:不能写成 h_bit(n)
  273. a[0].assign(beg, end);
  274. rng (i, 0, SZ(a)-1) {
  275. a[i+1].resize(n);
  276. rng(j, 0, n - (1 << i)) {
  277. a[i+1][j] = min(a[i][j], a[i][j+(1<<i)], Compare());
  278. }
  279. rng(j, n-(1 << i), n) {
  280. a[i+1][j] = a[i][j];
  281. }
  282. }
  283. }
  284.  
  285. T query(idx_t l, idx_t r) { // l <= r
  286. int i = h_bit(r - l + 1ul);
  287. return min(a[i][l], a[i][r+1-(1<<i)], Compare());
  288. }
  289. };
  290.  
  291. using poly = vector<int>;
  292.  
  293. void bit_reverse_swap(vector<int> &a) {
  294. int n = SZ(a);
  295. for (int i = 1, j = n >> 1, k; i < n - 1; i++) {
  296. if (i < j) swap(a[i], a[j]);
  297. //tricky
  298. for (k = n >> 1; j >= k; j -= k, k >>= 1); //inspect the highest "1"
  299. j += k;
  300. }
  301. }
  302.  
  303. void FFT(vector<int> &a, int type) {
  304. bit_reverse_swap(a);
  305. int n = SZ(a);
  306. for (int i = 2; i <= n; i *= 2) {
  307. const auto wi = fp(3, type * (mod - 1) / i); // i次单位根
  308. for (int j = 0; j < n; j += i) {
  309. auto w(1);
  310. for (int k = j, h = i >> 1; k < j + h; k++) {
  311. auto t = Mul(w, a[k + h]), u = a[k];
  312. a[k] = Add(u, t);
  313. a[k + h] = Sub(u, t);
  314. mul(w, wi);
  315. }
  316. }
  317. }
  318. const int inv = fp(n, -1);
  319. if (type == -1) for (auto &x : a) mul(x, inv);
  320. }
  321.  
  322. void fp(poly &a, const int n) {
  323. a.resize((ul)pow2((SZ(a)-1) * n + 1));
  324. FFT(a, 1);
  325. for(auto &x: a) x = fp(x, n);
  326. FFT(a, -1);
  327. }
  328. // DEBUG code by tourist
  329. string to_string(string s) {
  330. return '"' + s + '"';
  331. }
  332.  
  333. string to_string(const char* s) {
  334. return to_string((string) s);
  335. }
  336.  
  337. string to_string(bool b) {
  338. return (b ? "true" : "false");
  339. }
  340.  
  341. template <typename A, typename B>
  342. string to_string(pair<A, B> p) {
  343. return "(" + to_string(p.first) + ", " + to_string(p.second) + ")";
  344. }
  345.  
  346. template <typename A>
  347. string to_string(A v) {
  348. bool first = true;
  349. string res = "{";
  350. for (const auto &x : v) {
  351. if (!first) {
  352. res += ", ";
  353. }
  354. first = false;
  355. res += to_string(x);
  356. }
  357. res += "}";
  358. return res;
  359. }
  360.  
  361. void debug_out() { cerr << endl; }
  362.  
  363. template <typename Head, typename... Tail>
  364. void debug_out(Head H, Tail... T) {
  365. cerr << " " << to_string(H);
  366. debug_out(T...);
  367. }
  368.  
  369. #ifdef LOCAL
  370. #define debug(...) cerr << "[" << #__VA_ARGS__ << "]:", debug_out(__VA_ARGS__)
  371. #else
  372. #define debug(...) 42
  373. #endif
  374. // end DEBUG
  375. //////////////////*^"^*/////////////////////////////////////
  376. const int N = 1e5+5, M = 1005;
  377.  
  378. struct arc {
  379. int to, next;
  380. int cap, cost;
  381. };
  382.  
  383. const int nax = 2505;
  384. vi head(nax, -1);
  385. vector<arc> e;
  386.  
  387. void add_edge(int u, int v, int cap, int cost) {
  388. e.pb({v, head[u], cap, cost});
  389. head[u] = SZ(e) - 1;
  390. e.pb({u, head[v], 0, -cost});
  391. head[v] = SZ(e) - 1;
  392. }
  393.  
  394.  
  395.  
  396. int d[nax];
  397. const int inf = 100;
  398.  
  399. bool dijkstra(int s, int t) {
  400. pq<pii, greater<pii>> que;
  401. fill(all(d), inf);
  402. que.push({d[s] = 0, s});
  403. while (!que.empty()) {
  404. auto top = que.top();
  405. que.pop();
  406. int u = top.second;
  407. if (d[u] != top.first) continue;
  408. for (int i = head[u]; i != -1; i = e[i].next) {
  409. int v = e[i].to;
  410. if (e[i].cap > 0 && d[v] > d[u] + e[i].cost) {
  411. d[v] = d[u] + e[i].cost;
  412. que.push({d[v], v});
  413. }
  414. }
  415. }
  416. return d[t] < inf;
  417. }
  418.  
  419. vi cur(nax);
  420. vb vis(nax);
  421. int dfs(int u, int flow, int t) {
  422. if (u == t) return flow;
  423. vis[u] = true; // "vis[u] = true;" 不能写在 "if (u == t) return flow;" 之前!
  424. int pushed = 0;
  425. for (int &i = cur[u]; i != -1; i = e[i].next) {
  426. int v = e[i].to;
  427. if (!vis[v] && e[i].cap && d[v] == d[u] + e[i].cost) {
  428. int tmp = dfs(v, min(flow - pushed, e[i].cap), t);
  429. if (tmp) {
  430. e[i].cap -= tmp;
  431. e[i ^ 1].cap += tmp;
  432. pushed += tmp;
  433. if (flow == pushed) break;
  434. }
  435. }
  436. }
  437. vis[u] = false;
  438. return pushed;
  439. }
  440.  
  441. ll mcmf(int s, int t) {
  442. ll ans = 0;
  443. while (dijkstra(s, t)) {
  444. cur = head;
  445. while (true) {
  446. int f = dfs(s, INT_MAX, t);
  447. if (f) ans += f * d[t];
  448. else break;
  449. }
  450. }
  451. return ans;
  452. }
  453.  
  454.  
  455. int main() {
  456. FAST_READ
  457. #ifdef LOCAL
  458. ifstream in("main.in");
  459. cin.rdbuf(in.rdbuf());
  460. #endif
  461.  
  462. int n, m; scan(n, m);
  463. vv<int> h(n, vi(m));
  464. int tot = 0;
  465. For (row, h) For(x, row) scan(x), tot += x;
  466. int tar = tot / (n * m);
  467.  
  468. const int dx[]={0, 0, 1, -1}, dy[] = {1, -1, 0, 0};
  469. int s = n * m, t = s + 1;
  470.  
  471.  
  472. rng (i, 0, n) rng (j, 0, m) {
  473. int u = i * m + j;
  474. if (h[i][j] > tar) {
  475. add_edge(s, u, h[i][j] - tar, 0);
  476. }
  477. if (h[i][j] < tar) {
  478. add_edge(u, t, tar - h[i][j], 0);
  479. }
  480. rng (k, 0, 4) {
  481. int ni = i + dx[k], nj = j + dy[k];
  482. if(bet(0, n-1, ni) && bet(0, m-1, nj)) {
  483. add_edge(u, ni * m + nj, INT_MAX, 1);
  484. }
  485. }
  486. }
  487. dijkstra(s, t);
  488. println(mcmf(s, t));
  489.  
  490.  
  491.  
  492.  
  493.  
  494. #ifdef LOCAL
  495. cout << "Time elapsed: " << 1.0 * clock() / CLOCKS_PER_SEC << " s.\n";
  496. #endif
  497. return 0;
  498. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement