寒假训练题,比赛全名 XVIII Open Cup named after E.V. Pankratiev. Eastern Grand Prix。

Description

在三维整点空间 \(E\) 中,每个整点存了个数字,一开始数字都是 \(0\),每次对这样一个区域 \(S(x,y,z,a)\) 进行两种操作

  • \(M\) 操作 所有整点上的数字 \(+1\)
  • \(Q\) 操作 求所有整点上的数字和

\(M\) 操作全部做完后才会有 \(Q\) 操作。

\[ E:=\{(x,y,z)|1\le z\le y\le x\le n\}, \]

\[ S(x,y,z,a):=\{(x+u,y+v,z+w)|0\le w\le v\le u< a\}. \]

\[ n\le 100,\quad |M|\le 10^5,\quad |Q|\le 10^5. \]

Analysis

比赛时想了蛮久,有点思路但最终没能做出来的题……

首先要观察这个 \(E\)\(S\) 到底是啥形状,\(S\) 显然是边长为 \(a\) 的形状平移顶点到 \((x,y,z)\),什么 \(<a\) 等价成 \(\le a-1\) 即可(下文都按 \(\le a\) 处理),所以其实本质都一样,只需要考虑那个不等式限制。

离散的很丑,我们当然是先按连续的来,并且把下界 \(1\) 换成 \(0\) 这种更合适的起始点,先考虑已知 \(z=z_0\),对它进行切片降维,限制变成 \(z_0\le y\le x\le a\),这是个腰为 \((a-z_0)\) 的等腰直角三角形。考虑 \(z_0\) 从左到右变化,于是它就是一个三角锥,来自一个立方体从两个面对角线垂直切割的结果。这个很好想,大概长这样

Mathematica 绘制的区域

然后保证的是改完才会查,这就很像前缀和差分那一套,但是我想了很久也不知道怎么猜到差分式子才会是正确的形状。后来对内讨论 wjh 说他的做法就是直接手构……具体来说:

我们先考虑二维的,也就是每次平面上加一个等腰直角三角形时我们要怎么做

一次修改是加上这样的一个数阵:

\[ \begin{array}{|ccccc|} y\\ \uparrow\\ \hline 0&0&0&0&1\\ 0&0&0&1&1\\ 0&0&1&1&1\\ 0&1&1&1&1\\ 1&1&1&1&1&\rightarrow x\\ \hline \end{array} \]

\[ (1\le y\le x\le 5) \]

怎么办?其实没那么复杂,我们可以胡乱设计差分方式,目标是把修改操作涉及的非 \(0\) 元素降低到 \(\mathcal O(1)\) 规模。比如这里我们先横着差分一次

\[ \begin{array}{|ccccc|c} \hline 0&0&0&0&1&-1\\ 0&0&0&1&0&-1\\ 0&0&1&0&0&-1\\ 0&1&0&0&0&-1\\ 1&0&0&0&0&-1\\ \hline \end{array} \]

为了缩减斜线的规模,再斜着沿着右上方向差分一次(左下减右上)

\[ \begin{array}{|ccccc|c} 0&0&0&0&0&-1&1\\ \hline 0&0&0&0&0&-1&1\\ 0&0&0&0&0&-1&1\\ 0&0&0&0&0&-1&1\\ 0&0&0&0&0&-1&1\\ 1&0&0&0&0&-1\\ \hline \end{array} \]

这个时候发现还有一竖数,小菜一碟……竖着差分一下就完事了

\[ \begin{array}{|ccccc|c} 0&0&0&0&0&1&-1\\ 0&0&0&0&0&0&0\\ \hline 0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0\\ -1&0&0&0&0&0&1\\ 1&0&0&0&0&-1\\ \hline \end{array} \]

最后只有六个位置要改,大体就是这么个过程,需要三次差分,而如果你要复原回去,也就是按倒着的顺序求前缀和罢了。实际上不管你按照什么顺序都可以降到这个规模,只要你横竖斜都有。我觉得最佳的差分顺序应该是先斜着差分,然后发现三角阵直接变成一横一竖,目标就很明确了。

好了,说了一大套,三维要怎么做呢?按照上面的思路,也是先用斜着差分思路会更清晰,按照 \(y=x\) 正方向,这个时候我们发现做了一次之后所有所有修改都集中到 \({\color{cyan} x=a+1},{\color{orange} z=y}\) 平面上了,

一次差分过后数字集中的两个面

前一个好办,已经是个方的了,关键是后一个。于是我们再沿着 \(z=y\) 正方向差分,如此云云。可是我们突然意识到不对,这个比起刚刚的二维情形怎么繁这么多,所选用的轴逐渐失去了其对称性。到了三维空间,“斜”的概念其实也应该发生改变,第一步是否有更好的选择呢?我们应该试试直接往 \((1,1,1)\) 向量方向差分,这才是最“斜”的方案。

\((1,1,1)\)-向差分后数字的大致排布

果不其然,采用这个方向后结果非常明晰,剩下的点都在面 \({\color{skyblue} xOy},{\color{pink} x=a+1}\) 上!这剩下的都是两个等腰三角形,上文已经讨论过了,具体怎么差分我就不再说啦~毕竟不做也能猜到,肯定是 \(y=x,z=y\) 两个方向各差分一次,然后 \(x,y,z\) 三个方向各差分一次啦,总共六次。

Generation

人力推导三维差分实在太容易出错,于是我写了个程序生成其最终关键点

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
// 静かな海
#include <bits/stdc++.h>

using namespace std;

const int N = 53, M = 998244353;

int n = 10, mp[N][N][N];

#define fa(var) for (int var=1; var<=n; var++)
#define fd(var) for (int var=n; var; var--)

inline void dbgmp() {
fa(z) {
printf("\nz = %d layer.\n", z);

fa(x) fa(y)
printf("%d%c", mp[x][y][z], "\t\n"[y==n]);

puts("=============\n");
}
}

inline void diff(int dx, int dy, int dz) {
fd(x) fd(y) fd(z)
mp[x][y][z] -= mp[x-dx][y-dy][z-dz];
}

inline string get_id(int pos) {
if (pos <= 3) return to_string(pos);
else {
int val = pos - (n - 4);
if (!val) return "n";
else if (val > 0) return "n + " + to_string(val);
else return "n - " + to_string(-val);
}
}

inline void output_key() {
fa(x) fa(y) fa(z) {
if (mp[x][y][z]) {
printf("mp[%s][%s][%s] = %d;\n",
get_id(x).c_str(), get_id(y).c_str(),
get_id(z).c_str(), mp[x][y][z]);
}
}
}

signed main() {

fa(z)
for (int y=z; y<=n; y++)
for (int x=y; x<=n; x++)
mp[x][y][z] = 1;

// for -1
n+=4;

diff(1,1,1);
diff(1,1,0);
diff(0,1,1);
diff(1,0,0);
diff(0,1,0);
diff(0,0,1);

// dbgmp();
output_key();


return 0;
}

diff 是差分函数,get_id 是输出的格式化,我们要根据该坐标离 \(1\)\(n\) 的相对远近猜出它是常数位置还是和 \(n\) 有关的位置。思路就是如此暴力~

然后拿生成的关键点进行了一个验证,保证它可以生成对应的三角锥并且不会多加:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
// 静かな海
#include <bits/stdc++.h>

using namespace std;

inline int rd() {
int x=0,c,f=1;while(!isdigit(c=getchar()))if(c=='-')f=-1;
for(;isdigit(c);c=getchar())x=x*10+c-'0';return x*f;
}
const int N = 103, M = 998244353;

#define fa(var) for (int var=1; var<=n; var++)
#define fd(var) for (int var=n; var; var--)

int n, mp[N][N][N];

inline void summ(int dx, int dy, int dz) {
fa(x) fa(y) fa(z)
mp[x][y][z] += mp[x-dx][y-dy][z-dz];
}

inline void dbgmp() {
fa(z) {
printf("\nz = %d layer.\n", z);

fa(x) fa(y)
printf("%d%c", mp[x][y][z], "\t\n"[y==n]);

puts("=============\n");
}
}

signed main() {

n = rd();

mp[1][1][1] = 1;
mp[1][1][2] = -1;
mp[1][2][1] = -1;
mp[1][2][3] = 1;
mp[1][3][2] = 1;
mp[1][3][3] = -1;
mp[n + 1][1][1] = -1;
mp[n + 1][1][2] = 1;
mp[n + 1][n + 2][1] = 1;
mp[n + 1][n + 2][n + 3] = -1;
mp[n + 1][n + 3][2] = -1;
mp[n + 1][n + 3][n + 3] = 1;
mp[n + 2][2][1] = 1;
mp[n + 2][2][3] = -1;
mp[n + 2][n + 2][1] = -1;
mp[n + 2][n + 2][n + 3] = 1;
mp[n + 2][n + 4][3] = 1;
mp[n + 2][n + 4][n + 3] = -1;
mp[n + 3][3][2] = -1;
mp[n + 3][3][3] = 1;
mp[n + 3][n + 3][2] = 1;
mp[n + 3][n + 3][n + 3] = -1;
mp[n + 3][n + 4][3] = -1;
mp[n + 3][n + 4][n + 3] = 1;

n += 5;

summ(0,0,1);
summ(0,1,0);
summ(1,0,0);
summ(0,1,1);
summ(1,1,0);
summ(1,1,1);


dbgmp();

return 0;
}

嘿嘿嘿,summ 函数的参数就是 diff 原先的顺序反过来,没有丝毫更改。

什么?你说后面询问怎么回答?说实话一开始我被这后面的前缀和难住了……我开始类比其它差分和前缀和的关系,希望得到一些规律,比如其相关的项应该是坐标相差 \(1\) 的,某些地方的符号还需要更改,但我在生成的差分修改位置上做了不少试验,都过不了样例。

后来我想,其实我已经得到原数组了,所以直接多暴力维护几个方向(几个直的几个斜的)的前缀和就能做了。三维太麻烦,可以直接只考虑求每一层的等腰三角形,而等腰三角形用两个斜向+纵/横的前缀和就能维护。

然后回头一看佳禾的代码,喂,好像修改时的差分也能用那个维护???那我不是白忙活一场嘛!!orzzzzzzzz

这样总复杂度是 \(\mathcal O(|M|+n^3+n|Q|)\)。至于那个 \(n\),就算是为了代码复杂度牺牲了……(

于是我们就“轻松”完成了这题。

Solution

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
// 静かな海
#include <bits/stdc++.h>

using namespace std;

inline int rd() {
int x=0,c,f=1;while(!isdigit(c=getchar()))if(c=='-')f=-1;
for(;isdigit(c);c=getchar())x=x*10+c-'0';return x*f;
}
const int N = 113, M = 998244353;

#define fa(var) for (int var=1; var<=n; var++)
#define fd(var) for (int var=n; var; var--)

int n, mp[N][N][N];

inline void summ(int dx, int dy, int dz) {
fa(x) fa(y) fa(z)
mp[x][y][z] += mp[x-dx][y-dy][z-dz];
}

inline void add(int x, int y, int z, int a) {
mp[x][y][z] += 1;
mp[x][y][z + 1] += -1;
mp[x][y + 1][z] += -1;
mp[x][y + 1][z + 2] += 1;
mp[x][y + 2][z + 1] += 1;
mp[x][y + 2][z + 2] += -1;
mp[x + a][y][z] += -1;
mp[x + a][y][z + 1] += 1;
mp[x + a][y + a + 1][z] += 1;
mp[x + a][y + a + 1][z + a + 2] += -1;
mp[x + a][y + a + 2][z + 1] += -1;
mp[x + a][y + a + 2][z + a + 2] += 1;
mp[x + a + 1][y + 1][z] += 1;
mp[x + a + 1][y + 1][z + 2] += -1;
mp[x + a + 1][y + a + 1][z] += -1;
mp[x + a + 1][y + a + 1][z + a + 2] += 1;
mp[x + a + 1][y + a + 3][z + 2] += 1;
mp[x + a + 1][y + a + 3][z + a + 2] += -1;
mp[x + a + 2][y + 2][z + 1] += -1;
mp[x + a + 2][y + 2][z + 2] += 1;
mp[x + a + 2][y + a + 2][z + 1] += 1;
mp[x + a + 2][y + a + 2][z + a + 2] += -1;
mp[x + a + 2][y + a + 3][z + 2] += -1;
mp[x + a + 2][y + a + 3][z + a + 2] += 1;
}

using ll = long long;
// failed...
inline ll query(int x, int y, int z, int a) {
ll res1 = 0, res2 = 0;
res2 += mp[x - 1][y - 1][z - 1] * 1;
res2 += mp[x - 1][y - 1][z] * -1;
res2 += mp[x - 1][y][z - 1] * -1;
res2 += mp[x - 1][y][z + 1] * 1;
res2 += mp[x - 1][y + 1][z] * 1;
res2 += mp[x - 1][y + 1][z + 1] * -1;
res1 += mp[x + a - 1][y - 1][z - 1] * -1;
res1 += mp[x + a - 1][y - 1][z] * 1;
res2 += mp[x + a - 1][y + a][z - 1] * 1;
res1 += mp[x + a - 1][y + a][z + a + 1] * -1;
res2 += mp[x + a - 1][y + a + 1][z] * -1;
res1 += mp[x + a - 1][y + a + 1][z + a + 1] * 1;
res1 += mp[x + a][y][z - 1] * 1;
res1 += mp[x + a][y][z + 1] * -1;
res2 += mp[x + a][y + a][z - 1] * -1;
res1 += mp[x + a][y + a][z + a + 1] * 1;
res2 += mp[x + a][y + a + 2][z + 1] * 1;
res1 += mp[x + a][y + a + 2][z + a + 1] * -1;
res1 += mp[x + a + 1][y + 1][z] * -1;
res1 += mp[x + a + 1][y + 1][z + 1] * 1;
res2 += mp[x + a + 1][y + a + 1][z] * 1;
res1 += mp[x + a + 1][y + a + 1][z + a + 1] * -1;
res2 += mp[x + a + 1][y + a + 2][z + 1] * -1;
res1 += mp[x + a + 1][y + a + 2][z + a + 1] * 1;
return res1 - res2;
}

inline void dbgmp() {
fa(z) {
printf("\nz = %d layer.\n", z);

fa(x) fa(y)
printf("%d%c", mp[x][y][z], "\t\n"[y==n]);

puts("=============\n");
}
}

inline void dbgsx(ll sx[N][N][N]) {
fa(z) {
printf("\nz = %d layer.\n", z);

fa(x) fa(y)
printf("%d%c", sx[z][x][y], "\t\n"[y==n]);

puts("=============\n");
}
}

// z x y
ll s[N][N][N], sh[N][N][N], sv[N][N][N];

signed main() {

// #ifndef ONLINE_JUDGE
// freopen("input.in", "r", stdin);
// freopen("output.out", "w", stdout);
// #endif

n = rd(); int m = rd(), q = rd();

while (m--) {
int x = rd(), y = rd(), z = rd(), a = rd();
add(x, y, z, a);
}

// recover
summ(0,0,1), summ(0,1,0), summ(1,0,0), summ(0,1,1), summ(1,1,0), summ(1,1,1);
// sum
// summ(0,0,1), summ(0,1,0), summ(1,0,0), summ(0,1,1), summ(1,1,0), summ(1,1,1);

// dbgmp();


fa(z) {
fa(x) fa(y) s[z][x][y] = s[z][x-1][y-1] + mp[x][y][z];
fa(y) fa(x) sh[z][x][y] = sh[z][x-1][y] + s[z][x][y];
fa(x) fa(y) sv[z][x][y] = sv[z][x][y-1] + s[z][x][y];
}

// dbgsx(sv);

while (q--) {
int x = rd(), y = rd(), z = rd(), a = rd();
ll res = 0;
for (int l=0; l<a; l++)
res += (sv[z+l][x+a-1][y+a-1] - sv[z+l][x+a-1][y+l-1])
- (sh[z+l][x+a-2][y+l-1] - (x+l-2<0 ? 0 :sh[z+l][x+l-2][y+l-1]));

printf("%lld\n", res);
}

return 0;
}