tangjz的博客

判断回文可以用 hash 判定原串与逆串是否相等，找某个中心扩展的回文串可以用原串和逆串最长公共前缀。
允许字符串重排时，回文与每个字符出现的奇偶性有关。
一个串的本质不同回文子串个数是关于串长线性的，长度为 $n$ 的串在末尾加一个字符 $c$ 后如果想多一个本质不同的回文子串，那么必然是当前串的最长回文后缀。
没有长度为 2 的（或偶数长度）回文子串意味着任意连续 2 个字符不同，没有长度大于 1 的回文子串意味着任意连续 3 个字符不同。
最长回文子序列是原串和逆串的最长公共子序列。
回文只限制了等于和不等于的关系，可以用最小表示法表示回文。规模较大，限制较少时，也可采用下述方法表示 border。
从回文中心入手比较方便解决问题。例如从一个串里选择两个不相交的子区间拼成最长回文串，如果回文中心所在的那个区间没有包含极大回文子串，那么可以调整为包含极大回文子串。
回文 border 是回文串的回文前缀，也是它的回文后缀。由于周期的性质，border 可以按照长度之差分为 $\mathcal{O}(\log |S|)$ 组，每一组称为一个 border series，其中同组 border 的长度之差相等，具体来说是将 border 与更小 border 的长度之差相等的 border 分为一组（虽然倒过来分也没关系）。
在字符串末尾增加一个字符后，回文后缀只会有三种事件：后缀变长、后缀消失、后缀出现。其中后缀消失和后缀出现的总次数关于长度线性，而同一个 border series 的后缀除最长的 border 外其他 border 要么同时变长、要么同时消失。
回文树可以维护字符串的本质不同回文子串，并支持在两端添加字符。
AC 自动机计算 fail 的复杂度基于所有节点的深度之和不超过总串长，在后缀树和回文树上则没有这样的性质。如果对 kmp、AC 自动机、后缀树、回文树需要支持末尾增删操作，有两种做法可以选择，一种是记录 border series 强行计算，一种是额外记录每个节点增加某个字符能到达的节点信息，前者是 $\mathcal{O}(\log |S|)$ 的，后者是 $\mathcal{O}(\sum)$ 或 $\mathcal{O}(\log \sum)$ 的。
考虑任意一组 border（不只是回文 border），不妨设 border 长度从小到大依次为 $l_1, l_2, \cdots, l_{m - 1}$ ，并且我们约定 $l_0 = 0, l_m = |S|$，那么对于任意的 $i = 1, 2, \cdots, m - 1$，要么有 $l_{i + 1} > 2 l_i$，要么存在 $j = 0, 1, \cdots, i - 1$ 使得 $l_{i + 1} = 2 l_i - l_j$，而且当 $l_{i + 1} - l_i$ 是 $l_i - l_{i - 1}$ 的倍数时，必有 $l_{i + 1} = 2 l_i - l_{i - 1}$。此外，这些性质限制的序列一定能表示成某个字符串 $S$ 的 border。

BZOJ 4044 [Cerc2014] Virus synthesis & BZOJ 4163 字符串合成

2018-03-01 20:05:34 By tangjz

网上题解都在说什么啊……

题意

BZOJ 4044 & BZOJ 4163：

初始有一个空串，利用下面的操作构造给定串 $S$。

串开头或末尾加一个字符
串开头或末尾加一个该串的逆串

最小化操作数，$|S| \leq 10^5$。

分析

操作 2 会生成长度为偶数的回文串，而一个字符串的本质不同回文子串个数是 $\mathcal{O}(|S|)$ 的，尝试将操作步骤按照操作 2 进行划分。

在最后一次使用操作 2 后，只会在偶回文子串的两端添加字符，设 $f(T)$ 表示构造偶回文子串 $T$ 的最小操作数，则对于初始的空串 $\epsilon$ 有 $f(\epsilon) = 0$，而答案 $\mathrm{ans}$ 为 $\mathop{\min}\limits_{T \in \mathrm{substring}(S)}\{f(T) + |S| - |T|\}$。

如果从未使用操作 2，则类似地有 $f(T) \leq |T|$ 和 $\mathrm{ans} \leq |S|$。

在第一次使用操作 2 前，只会通过添加字符获得偶回文子串的前一半，即 $f(T) \leq \dfrac{|T|}{2} + 1$。

在相邻两次操作 2 之间，有可能在偶回文子串 $P$ 的两端添加字符，再通过操作 2 变成 $T$，其中 $P$ 是 $T$ 的前半部分的子串（也是后半部分的）。

如果 $P$ 不是 $T$ 的前缀，那么一定是在 $P$ 的前端增加了字符。去掉生成 $T$ 的最后一次操作 1 而保留之后的操作 2，则可以用 $P$ 生成字符串 $T'$，其中 $T'$ 是 $T$ 去掉首尾各一个字符组成的串。这意味着，如果生成 $T'$ 经过至少一次操作 2，即 $T' \neq \epsilon$，那么有 $f(T) \leq f(T') + 1$。

如果 $P$ 是 $T$ 的前缀，那么 $P$ 是 $T$ 的 border，在 $P$ 的后端增加字符后使用一次操作 2 可以得到 $T$。事实上，$P$ 是 $T$ 的前半部分中最长的偶 border，因为更短的 border 只会使得操作 1 的操作次数增多。因此可以找到 $T$ 的长度不超过 $\dfrac{|T|}{2}$ 的最长 border $P'$，找到 $P'$ 的最长偶 border $P$，那么有 $f(T) \leq f(P) + \dfrac{|T|}{2} - |P| + 1$。

事实上，不需要找到 $P'$ 的最长偶 border $P$。如果 $P'$ 不是 $T$ 的偶 border，那么构造 $T'$ 的最优操作一定会是使用 $P$ 构造的，直接忽略不是 $P'$ 的长度不是偶数的转移即可。

附：一个串的 border 既是它的前缀又是它的后缀，一个串的 border 的 border 也是这个串的 border。

综上所述，我们可以对 $S$ 构建回文树，对每个状态维护长度不超过一半的最长 border，并进行 DP，时间复杂度 $\mathcal{O}(|S|)$。

代码

贴一个 BZOJ 4044 的代码仅供交流学习。

#include <bits/stdc++.h>
using namespace std;
const int maxn = (int)1e5 + 3, maxc = 4;
int t, n, tot, cur, len[maxn], fail[maxn], nxt[maxn][maxc + 1];
char buf[maxn], str[maxn];
int newNode(int length, int prev = 0) {
    len[tot] = length;
    fail[tot] = prev;
    memset(nxt[tot], 0, sizeof nxt[0]);
    return tot++;
}
void init() {
    n = tot = 0;
    cur = newNode(0);
    fail[cur] = newNode(-1, cur);
    buf[0] = '\0';
}
int slen, trans[257], half[maxn], f[maxn], ans;
void append(char ch) {
    int o = trans[buf[++n] = ch];
    for( ; buf[n - len[cur] - 1] != buf[n]; cur = fail[cur]);
    if(!nxt[cur][o]) {
        int p = newNode(len[cur] + 2), q;
        for(q = fail[cur]; buf[n - len[q] - 1] != buf[n]; q = fail[q]);
        fail[p] = nxt[q][o];
        fa[p] = cur;
        for(q = half[cur]; buf[n - len[q] - 1] != buf[n]; q = fail[q]);
        for(q = nxt[q][o]; len[q] > len[p] / 2; q = fail[q]);
        half[p] = q;
        nxt[cur][o] = p;
        if(!(len[p] & 1)) {
            f[p] = (cur ? f[cur] : len[p] / 2) + 1;
            if(!(len[q] & 1))
                f[p] = min(f[p], f[q] + len[p] / 2 - len[q] + 1);
            ans = min(ans, f[p] + slen - len[p]);
        }
    }
    cur = nxt[cur][o];
}
int main() {
    trans['A'] = 0;
    trans['G'] = 1;
    trans['C'] = 2;
    trans['T'] = 3;
    scanf("%d", &t);
    while(t--) {
        scanf("%s", str);
        ans = slen = strlen(str);
        init();
        for(int i = 0; i < slen; ++i)
            append(str[i]);
        printf("%d\n", ans);
    }
    return 0;
}

BZOJ 3452 Tyvj1955 Lunatic & ProjectEuler 368 A Kempner-like series

2018-02-27 21:51:39 By tangjz

这是一道数学分析背景的数位 DP 经典题……

为什么 BZOJ 上没人做呢……

题意

BZOJ 3452：给定 $n$ 个不含前导零的数字串 $s_1, s_2, \cdots, s_n$，对于所有满足十进制表示中不含任何一个给定数字串的正整数 $x$，求它们的倒数和，四舍五入精确到小数点后第四位。保证 $n, |s_i| \leq 3~(i = 1, 2, \cdots, n)$。

ProjectEuler 368：对于所有满足十进制表示中不含连续三个相同数位的正整数 $x$，求它们的倒数和，四舍五入精确到小数点后第十位。

~~注意：Tyvj 上的数据是错的~~

分析

对于任意一个正整数 $x$，我们可以利用字符串匹配检查它是否对答案产生贡献，而且可以发现这个正整数的任意前缀数字都是可以产生贡献的，因此可以考虑数位 DP 套字符串匹配，按照数位划分状态的同时将字符串匹配的信息存储在状态中。

由于题目涉及到多个字符串间的匹配，因此对于 BZOJ 3452 我们采用 AC 自动机进行字符串匹配，而对于 ProjectEuler 368 我们则只需要记录最后一个数位以及它是否与倒数第二个数位相同。这里假设读者已经具备相关知识储备。

在数位 DP 时，我们需要考虑在一些不含前导 $0$ 的数字 $x~(x \in S)$ 后添加一个数位 $d$ 并维护倒数和，注意到有 $$\frac{1}{10 x + d} = \sum_{k \geq 1}{\frac{(-d)^{k - 1}}{(10 x)^k}},$$ 于是有 $$\sum_{x \in S}{\frac{1}{(10 x + d)^n}} = \sum_{k \geq n}{{k - 1 \choose n - 1} \frac{(-d)^{k - n}}{10^k} \sum_{x \in S}{\frac{1}{x^k}}},$$ 那么当 $x$ 和 $k$ 足够大时， $\sum\limits_{x \in S}{\dfrac{1}{x^k}}$ 是极小的，因此可以为 $x$ 确定长度下界、为 $k$ 确定上界，对答案进行近似计算。

具体来说，注意到 $x = 1$ 时 $k$ 无法取到足够好的上界，因此我们可以设定阈值 $L$，对于 $x < 10^L$ 的情况采用枚举的方法，对于 $x \geq 10^L$ 的情况采用数位 DP 近似。设定 $L$ 后，可以直观地设定 $k$ 的上限 $K$，只需让 $\sum\limits_{x \in S}{\dfrac{1}{x^K}}$ 的最大值也远远小于答案允许的误差即可。

现在继续考虑如何用数位 DP 计算 $x \geq 10^L$ 的贡献，我们定义 $f(i, j, k)$ 表示长度为 $i$、匹配状态为 $j$ 的所有能够产生贡献的数字的 $k$ 次方倒数和，若可行的匹配状态集合为 $P$，则对答案的贡献为 $$\sum_{i > L} \sum_{j \in P}{f(i, j, 1)},$$ 可以再次设定长度的上界进行近似计算，不过这里为了保证精度更好可以采用更好的方法。

设 $M = |P|$，数位 DP 的初始状态 $f(L, j, k)$ 可以写成一个长度为 $M K$ 的行向量 $A$，而 $f(i + 1, j, k)$ 与 $f(i, j, k)$ 之间的转移也可以写成一个大小为 $M K \times M K$ 的矩阵 $B$，则 $x \geq 10^L$ 的全部 DP 信息为 $$\sum_{k \geq 1}{A B^k} = A B (I - B)^{-1},$$ 其中 $I$ 是大小为 $M K \times M K$ 的单位矩阵，所求即结果矩阵的 $M$ 个特定元素之和。

于是本题就在 $\mathcal{O}(10^L + M^3 K^3)$ 的时间复杂度中解决了。

代码

贴一个 ProjectEuler 368 的代码仅供交流学习。

#include <bits/stdc++.h>
typedef double DB;
const int maxd = 10, maxs = 101, maxl = 6;
int m, e, sz;
DB lpw[maxd + 1][maxs], c[maxs][maxs];
DB f[maxs], mat[maxs][maxs << 1 | 1], ans;
void dfs(int dep, int val, int msk) {
    if(dep == maxl) {
        DB prd = 1;
        for(int i = 0; i < e; ++i)
            f[i * m + msk] += 1 / (prd *= val);
        return;
    }
    ans += 1.0 / val;
    for(int i = 0; i < maxd; ++i) {
        if(msk & 1 && (msk >> 1) == i)
            continue;
        dfs(dep + 1, val * maxd + i, i << 1 | ((msk >> 1) == i));
    }
}
int main() {
    m = maxd << 1;
    e = (maxs - 1) / m;
    sz = e * m;
    lpw[0][0] = 0;
    for(int i = 1; i <= maxd; ++i) {
        lpw[i][0] = 0;
        lpw[i][1] = log(i);
        for(int j = 2; j < e; ++j)
            lpw[i][j] = lpw[i][j - 1] + lpw[i][1];
    }
    lpw[maxd][e] = lpw[maxd][e - 1] + lpw[maxd][1];
    for(int i = 0; i < e; ++i) {
        c[i][0] = c[i][i] = 1;
        for(int j = 1; j < i; ++j)
            c[i][j] = c[i - 1][j - 1] + c[i - 1][j];
    }
    for(int i = 0; i < sz; ++i)
        mat[i][i] = mat[i][sz + i] = 1;
    for(int i = 0; i < m; ++i)
        for(int d = 0; d < maxd; ++d) {
            if(i & 1 && (i >> 1) == d)
                continue;
            int ii = d << 1 | ((i >> 1) == d);
            for(int j = 0; j < e; ++j)
                for(int jj = 0; jj <= j; ++jj) {
                    if(!d && jj < j)
                        continue;
                    mat[j * m + i][jj * m + ii] -= ((j - jj) & 1 ? -1 : 1) * exp(log(c[j][jj]) + lpw[d][j - jj] - lpw[maxd][j + 1]);
                }
        }
    for(int i = 0, j, k; i < sz; ++i) {
        for(j = i, k = j + 1; k < sz; ++k)
            if(fabs(mat[k][i]) > fabs(mat[j][i]))
                j = k;
//        assert(fabs(mat[j][i]) > 1e-20);
        DB prd = mat[j][i];
        for(k = i; k < sz << 1; ++k) {
            std::swap(mat[j][k], mat[i][k]);
            mat[i][k] /= prd;
        }
        for(j = 0; j < sz; ++j) {
            if(i == j)
                continue;
            DB prd = mat[j][i];
            for(int k = i; k < sz << 1; ++k)
                mat[j][k] -= prd * mat[i][k];
        }
    }
    for(int i = 1; i < maxd; ++i)
        dfs(1, i, i << 1);
    for(int i = 0; i < sz; ++i)
        for(int j = 0; j < m; ++j)
            ans += f[i] * mat[i][sz + j];
    printf("%.10f\n", (double)ans);
    return 0;
}

Written with StackEdit.

字符串数位DP 高斯消元
阅读全文

共 3 篇博客