Problem
在 A 国,一年有 (50!=1 imes2 imesdots imes50) 天,每天的编号从 (1) 到 (50!) 。
A 国的土地上经常发生自然灾害。据 A 国科学家的妍究, A 国有 (n) 种自然灾害,每种自然灾害的发生都有周期性。
具体地,第 (i) 种自然灾害可以用一个只包含 Y
和 N
组成的字符串组成。
设这个字符串的长度为 (l) 。如果这个字符串的第 (j) 个字符为 Y
,那么这就表示对于所有的 (kin ext{N}) ,在第 (k imes l+j) 天会发生第 (i) 种自然灾害,否则不会发生。
一天同时发生的自然灾害种数决定了该天整个 A 国境内的危险程度。
所以,你需要对于所有的 (0le ile n) ,求出一年内有多少天同时发生的自然灾害种数恰好为 (i) 。
(n leq 30),字符串的长度不超过 (50)。
有几个部分分:
- 所有字符串长度的乘积不超过 (1111111);
- 所有字符串长度两两互质或相等;
- 所有字符串长度不超过 (48)。
Solution
校内模拟赛的题,被锤爆了,深有体会。
部分分 1
将第 (i) 个字符串的长度称为 (l_i)。
将所有串的长度的 (mathrm{lcm}) 记为 (L)。
那么我们可以考虑将所有串倍长到 (L) 的长度,然后就可以直接计算 (1 sim L) 内的答案。
所有答案只要乘一个 (frac{50!}{L}) 的系数即可。
部分分 2
我们考虑一下,对于一个 (x),决定最终这一天发生的所有自然灾害数量实际上是所有的 (x mod p(1 leq p leq 50)) 的值。
如果我们反过来,去考虑每个 (x mod l_i) 的值,如果我们钦定了 (x equiv a_i pmod{l_i}),那么就相当于列出一个同余方程组,并且这个同余方程组在模 (L) 意义下最多只有一个解。有解的充要条件是:
因为部分分里所有字符串长度两两互质或相等,互质就忽略上面的条件,相等就直接枚举。
如果忽略了上面的条件,那么只要设计一个 DP,(f(i,j)) 表示在 (1 sim mathrm{lcm}(l_1,l_2,dots ,l_i)) 天中,前 (i) 个灾害,一共发生了 (j) 次的天数。
接下来的思路
像上面那样直接 DP 相当于没有限制下面这个方程一定有解。
考虑如何限制有解的条件。有一个思路是,取一个参数 (K),令 (t_i=frac{l_i}{mathrm{gcd}(l_i,K)}),使得不同的 (t_i) 两两互质。相当于加一个同余方程,使得方程组变成:
这样的目的是枚举第一个方程之后,后面的方程之间互相就不会矛盾了。这样我们就能通过将第一个方程与其他每个方程联立的方式,知道具体在原来的字符串是哪个位置,从而计算 DP 时需要的贡献:
于是现在 (K) 的选取就是关键的问题。再次注意 (K) 需要满足的条件:令 (t_i=frac{l_i}{mathrm{gcd}(l_i,K)}),不同的 (t_i) 需要两两互质。
算法一:根号分类的思想
遇到这种小范围的数论题,可以考虑一下按根号将质因子分类的经典思路(套路):[NOI2015]寿司晚宴
大于 (sqrt {50}) 的质因子在一个 (l_i) 最多只会出现一次。不妨取 (K=2^5 imes 3^3 imes 5^2 imes 7^2) 先试试。
求出每个 (l_i) 唯一大于根号的质因子(没有就看成 (1)),即 (t_i=frac{l_i}{mathrm{gcd}(l_i,K)}),然后根据 (t_i) 将所有 (l_i) 分类,这样不同的 (t_i) 显然是互质的。用和部分分 2 一样的 DP,(t_i) 相同的一起转移即可。最后同样要乘上一个系数
复杂度大概是 (mathcal Oleft(K(mcdot n^2+sum t_i) ight)),其中 (m) 表示不同的 (t_i) 的个数。发现这样不太过得去。
我们考虑把 (K) 取成 (K=2^5 imes 3^3 imes 5^2),然后强制把 (7|t_i) 的 (t_i) 都设成 (49)(可以认为是将满足 (7|t_i, t_i eq 49) 的那些 (7) 倍长),这样就把时间复杂度降低了。
实现的时候有一些技巧,可以避免使用中国剩余定理,可以看代码。
#include <bits/stdc++.h>
template <class T>
inline void read(T &x)
{
static char ch;
while (!isdigit(ch = getchar()));
x = ch - '0';
while (isdigit(ch = getchar()))
x = x * 10 + ch - '0';
}
template <class T>
inline void putint(T x)
{
static char buf[25], *tail = buf;
if (!x)
putchar('0');
else
{
for (; x; x /= 10) *++tail = x % 10 + '0';
for (; tail != buf; --tail) putchar(*tail);
}
}
const int MaxL = 55;
const int MaxN = 32;
const int mod = 1e4 + 7; //prime
const int M = 32 * 27 * 25; //Solution 中说的参数 K
int n, L = 50;
bool occ[MaxL];
int cnt[MaxL][MaxL];
std::vector<int> vec[MaxL];
int ans[MaxN];
int f[MaxN], g[MaxN];
bool vis[MaxL];
inline void add(int &x, const int &y)
{
x += y;
if (x >= mod)
x -= mod;
}
int main()
{
freopen("disaster.in", "r", stdin);
freopen("disaster.out", "w", stdout);
scanf("%d", &n);
for (int i = 1; i <= n; ++i)
{
static char s[MaxL];
scanf("%s", s);
int len = strlen(s);
occ[len] = true;
for (int j = 0; j < len; ++j)
if (s[j] == 'Y')
++cnt[len][j];
}
for (int i = 1; i <= L; ++i)
if (occ[i])
{
int t = i % 7 == 0 ? 49 : i / std::__gcd(i, M);
vec[t].push_back(i);
vis[t] = true;
}
vis[32] = vis[27] = vis[25] = true;
for (int tM = 0; tM < M; ++tM)
{
memset(f, 0, sizeof(f));
f[0] = 1;
for (int P = 1; P <= L; ++P)
if (!vec[P].empty())
{
static int tf[MaxN], cur[MaxN];
memset(tf, 0, sizeof(tf));
memset(g, 0, sizeof(g));
int im = vec[P].size();
for (int i = 0; i < im; ++i)
cur[i] = tM % vec[P][i];
for (int tim = 0; tim < P; ++tim)
{
int count = 0;
for (int i = 0; i < im; ++i)
{
int len = vec[P][i];
count += cnt[len][cur[i]];
cur[i] = (cur[i] + M) % len; //注意这里的处理
//备注:这里由于 K 取得不好,实现得不够精细,实际的复杂度不是上面所说的那样,比较好的实现可以看算法二
}
++g[count];
}
for (int i = 0; i <= n; ++i)
if (f[i])
{
for (int j = 0; i + j <= n; ++j)
add(tf[i + j], 1LL * f[i] * g[j] % mod);
}
for (int i = 0; i <= n; ++i)
f[i] = tf[i];
}
for (int i = 0; i <= n; ++i)
add(ans[i], f[i]);
}
for (int i = 0; i <= n; ++i)
{
for (int j = 1; j <= L; ++j)
if (!vis[j])
ans[i] = 1LL * ans[i] * j % mod;
putint(ans[i]);
putchar('
');
}
return 0;
}
算法二:省去冗余枚举
我们发现算法一的 (K) 取的并不太好,因为我们发现,并没有必要让 (t_i) 一定是大于根号的质因子,如果是不超过根号的也是可以的,但是要满足不同的 (t_i) 互质的条件。
形式化地,我们考虑每个质因子 (p_i) 在 (1sim 50) 范围内的最大可能的次数 (q_i)。取 (K=prod p_i^{q_i-1}=2^4 imes 3^2 imes 5 imes 7),不难发现这个 (K) 取得恰到好处。
学习这种做法的时候顺便学到了一种简化的实现方法,我们将第 (i) 个串的长度倍长到 (t_i imes K),即 (mathrm{lcm}(K,l_i)),并且记一个数组 (cnt_{i,j}) 表示倍长后所有 (t_k=i) 的串的第 (j) 位的总和。那么对相同的 (t_i) 转移的时候,直接利用对应位置的 (cnt) 转移即可。
时间复杂度是 (mathcal Oleft(K(mcdot n^2+sum t_i) ight))。
//orz pen
//What a magical solution !
#include <bits/stdc++.h>
const int K = 16 * 9 * 5 * 7;
const int MaxN = K * 50;
const int MaxM = 20;
const int mod = 1e4 + 7;
int n, m;
char s[MaxN];
bool vis[MaxN];
int l[MaxN], pos[MaxN];
int ans[MaxN];
int f[MaxM][32];
int cnt[MaxM][MaxN], g[32];
inline void add(int &x, const int &y)
{
x += y;
if (x >= mod)
x -= mod;
}
inline int get_pos(int x)
{
if (!vis[x])
vis[l[pos[x] = ++m] = x] = true;
return pos[x];
}
inline int qpow(int x, int y)
{
int res = 1;
for (; y; y >>= 1, x = 1LL * x * x % mod)
if (y & 1)
res = 1LL * res * x % mod;
return res;
}
int main()
{
freopen("disaster.in", "r", stdin);
freopen("disaster.out", "w", stdout);
scanf("%d", &n);
for (int i = 1; i <= n; ++i)
{
scanf("%s", s);
int len = strlen(s);
int x = K * len / std::__gcd(K, len);
int pos = get_pos(x);
for (int j = 0; j < x; ++j)
cnt[pos][j] += s[j % len] == 'Y';
}
for (int k = 0; k < K; ++k)
{
f[0][0] = 1;
for (int i = 1; i <= m; ++i)
{
for (int j = 0; j <= n; ++j)
g[j] = 0;
for (int j = k; j < l[i]; j += K)
++g[cnt[i][j]];
for (int j = 0; j <= n; ++j)
{
f[i][j] = 0;
for (int p = 0; p <= j; ++p)
add(f[i][j], 1LL * f[i - 1][j - p] * g[p] % mod);
}
}
for (int i = 0; i <= n; ++i)
add(ans[i], f[m][i]);
}
int cur = 1;
for (int i = 1; i <= 50; ++i)
cur = 1LL * cur * i % mod;
for (int i = 1; i <= m; ++i)
cur = 1LL * cur * qpow(l[i] / K, mod - 2) % mod;
cur = 1LL * cur * qpow(K, mod - 2) % mod;
for (int i = 0; i <= n; ++i)
{
ans[i] = 1LL * ans[i] * cur % mod;
printf("%d
", ans[i]);
}
return 0;
}