题解：P8164 [JOI 2022 Final] 沙堡 2 (Sandcastle 2)

提供一个比较聪明的做法。

首先考虑一个子矩形合法的条件，这当且仅当他在值域上的前驱和后继都和他直接相邻，我们从最大的点出发，每次选小于自己的最大的点出发。根据这个条件暴力判断可以做到 O(H^3W^2)。

不妨动用我们的人类智慧，对每一个点的状态做一个哈希，具体地：

令 x 为一个极大的值，f_{i,j} 为和 (i,j) 这个点相邻的点中，小于 a_{i,j} 的最大的点的值，不存在值为 0。如果四周的点都不比这个点大，这个点的哈希值为 x-f_{i,j}，否则为 a_{i,j}-f_{i,j}。

这样哈希之后，我们发现：一个合法的矩形的权值和为 x。这其实就是存在一条哈密顿路径。

这个条件是充要的，如果矩形不能被一条哈密顿路径覆盖，他的权值和一定大于 x。否则，每个点都会向 f_{i,j} 所在的方向走，这样每个点的权值就可以被恰好抵消。

枚举想计算的子矩形的上下边界，对矩形的右端点坐标进行扫描，每个点有作为左边界上的点，右边界上的点，中间的点三种不同的权值，可以预处理。

在对右边进行扫描的过程中维护每个位置的权值，使用哈希表或者 std::map 即可实现 O(H^2W) 或者 O(H^2W\log W)。

注意到 \min(H,W)\le\sqrt{HW}，选择较短的那一边作为 H 可以做到 O(HW\sqrt{HW})。

这个做法不适用与矩形长/宽为 1，需要特殊处理。

#include<iostream>
#include<vector>
#include<algorithm>
#include<unordered_map>
using ll = long long; const int N = 1e5; const ll inf = 2e18;
std::vector<std::vector<int>> a; std::unordered_map<ll, int> mp;
int dx[4] = {0, 1, 0, -1}, dy[4] = {1, 0, -1, 0}, n, m, t[N], ct; ll sm, rs;
ll hs(int l, int r, int x, int y, int op){
    int mx = 0, mn = 0; for(int i = 0; i < 4; i++){
        int tx = x + dx[i], ty = y + dy[i];
        if(tx >= l && tx <= r && ty >= 1 && ty <= m &&
            (op == 1 ? (dy[i] != -1) : (op == 2 ? (dy[i] != 1) : 1)))
        a[tx][ty] > a[x][y] ? mx = a[tx][ty] : mn = std::max(mn, a[tx][ty]);
    } return mx ? a[x][y] - mn : inf - mn;
}
int main(){
    //freopen("sand.in", "r", stdin), freopen("sand.out", "w", stdout);
    std::ios::sync_with_stdio(false), std::cin.tie(0);
    std::cin >> n >> m; a.resize(n + 1, std::vector<int>(m + 1));
    for(int i = 1; i <= n; i++) for(int j = 1; j <= m; j++) std::cin >> a[i][j];
    if(n > m){
        std::vector<std::vector<int>> b(m + 1, std::vector<int>(n + 1));
        for(int i = 1; i <= m; i++) for(int j = 1; j <= n; j++)
            b[i][j] = a[j][i]; std::swap(n, m), std::swap(a, b);
    }
    for(int i = 1; i <= n; i++) for(int j = 1; j <= m; j++) t[++ct] = a[i][j];
    std::sort(t + 1, t + ct + 1); ct = std::unique(t + 1, t + ct + 1) - t;
    for(int i = 1; i <= n; i++) for(int j = 1; j <= m; j++)
        a[i][j] = std::lower_bound(t + 1, t + ct + 1, a[i][j]) - t + 1;
    std::vector<std::vector<ll>> f(n + 1, std::vector<ll>(m + 1)),
        g(n + 1, std::vector<ll>(m + 1)), h(n + 1, std::vector<ll>(m + 1));
    std::vector<ll> F(m + 1), G(m + 1), H(m + 1); for(int l = 1; l <= n; l++){ 
        for(int i = 1; i <= m; i++) F[i] = f[l][i] = hs(l, l, l, i, 0),
            G[i] = g[l][i] = hs(l, l, l, i, 1), H[i] = h[l][i] = hs(l, l, l, i, 2);
        for(int r = l + 1; sm = 0, mp = std::unordered_map<ll, int>(), r <= n; r++){
            for(int i = 1; i <= m; i++){
                F[i] += hs(l, r, r - 1, i, 0) - f[r - 1][i] + (f[r][i] = hs(l, r, r, i, 0));
                G[i] += hs(l, r, r - 1, i, 1) - g[r - 1][i] + (g[r][i] = hs(l, r, r, i, 1));
                H[i] += hs(l, r, r - 1, i, 2) - h[r - 1][i] + (h[r][i] = hs(l, r, r, i, 2));
                rs += mp[inf - sm - H[i]]; sm += F[i], ++mp[G[i] - sm];
            }
        }
    }
    for(int i = 1; i <= n; i++){
        int x = 0, y = 0; for(int j = 1; j < m; j++)
            a[i][j] < a[i][j+1] ? (rs += j-y, x = j) : (rs += j-x, y = j);
    }
    for (int i = 1; i <= m; ++i) {
        int x = 0, y = 0; for(int j = 1; j < n; j++)
            a[j][i] < a[j+1][i] ? (rs += j-y, x = j) : (rs += j-x, y = j);
    }
    std::cout << rs + 1ll * n * m << '\n';
}