286DIV1E. Mr. Kitayuta's Gift

Posted 2020-08-30 ichneumon

题目大意

给定一个由小写字母构成的字符串$s$，要求添加$n(n\le 10^9)$个小写字母，求构成回文串的数目。

简要题解

$n$辣么大，显然要矩阵快速幂嘛。

考虑从两端开始构造以s ss为子串的回文串，该回文串长度为$N=n+s$，每次添加相同的字符，则需要$(N+1)/2$次，则用dp来计算并使用矩阵乘法来优化转移会得到一个$O(|s|^6\log N)$的算法，显然是不可接受的。

考虑这个dp做法，设$f[i][j][k]$表示从两端添加了$k$次字符，原来的$s$的子串$s_{ij}$任然需要被构造出来时的方案数，我们发现会存在两种转移，一种是添加了一个字符，$s_{ij}$转移到$s_{i+1,j-1}$，这是因为$s_i=s_j$，添加一个字符可以消去$s$中两个字符，记为n25转移；另一种是$s_i\not = s_j$只消去一个字符，记为n25转移。同时要注意到，每种转移都有到自身的转移（即在回文串两端添加字符后$s_{ij}$不变，n25转移有25种添加方式转移到自身，n24转移有24种添加方式转移到自身）

我们发现其实dp可以拆成两部分来计算，一部分是计算存在$k$个n24转移的转移路径方案数，这样一来，n25的转移数为$(|s|-k+1)/2$，那么还需要添加$(N+1)/2 - k - (|s|-k+1)/2$个到自身的转移，这么一来，n24转移和n25转移的顺序就没有意义了，不妨把n24转移放到前面，n25转移放到后面，然后建图，从start到end的转移方案数就是这个图的临接矩阵的$(N+1)/2$次方。这么一来，复杂度是$O(|s|^4 \log N)$，依旧不可接受。但注意到这$k$种建图是可以合并成一个有$O(|s|)$个起点，$O(|s|)$个终点的图的，所以只需要在这一个图上做矩阵快速幂就好，复杂度为$O(|s|^3\log N)$，问题解决。还有一个优化是，该矩阵是上三角矩阵，可以利用这点来加速矩阵乘法。以上讨论基于$N$是偶数，当$N$是奇数时，注意到只有$s_{ii}$出发的n25转移是不合法的，减去即可。

  1 #include <bits/stdc++.h>
  2 using namespace std;
  3 namespace my_header {
  4 #define pb push_back
  5 #define mp make_pair
  6 #define pir pair<int, int>
  7 #define vec vector<int>
  8 #define pc putchar
  9 #define clr(t) memset(t, 0, sizeof t)
 10 #define pse(t, v) memset(t, v, sizeof t)
 11 #define bl puts("")
 12 #define wn(x) wr(x), bl
 13 #define ws(x) wr(x), pc(‘ ‘)
 14     const int INF = 0x3f3f3f3f;
 15     typedef long long LL;
 16     typedef double DB;
 17     inline char gchar() {
 18         char ret = getchar();
 19         for(; (ret == ‘\n‘ || ret == ‘\r‘ || ret == ‘ ‘) && ret != EOF; ret = getchar());
 20         return ret; }
 21     template<class T> inline void fr(T &ret, char c = ‘ ‘, int flg = 1) {
 22         for(c = getchar(); (c < ‘0‘ || ‘9‘ < c) && c != ‘-‘; c = getchar());
 23         if (c == ‘-‘) { flg = -1; c = getchar(); }
 24         for(ret = 0; ‘0‘ <= c && c <= ‘9‘; c = getchar())
 25             ret = ret * 10 + c - ‘0‘;
 26         ret = ret * flg; }
 27     inline int fr() { int t; fr(t); return t; }
 28     template<class T> inline void fr(T&a, T&b) { fr(a), fr(b); }
 29     template<class T> inline void fr(T&a, T&b, T&c) { fr(a), fr(b), fr(c); }
 30     template<class T> inline char wr(T a, int b = 10, bool p = 1) {
 31         return a < 0 ? pc(‘-‘), wr(-a, b, 0) : (a == 0 ? (p ? pc(‘0‘) : p) : 
 32             (wr(a/b, b, 0), pc(‘0‘ + a % b)));
 33     }
 34     template<class T> inline void wt(T a) { wn(a); }
 35     template<class T> inline void wt(T a, T b) { ws(a), wn(b); }
 36     template<class T> inline void wt(T a, T b, T c) { ws(a), ws(b), wn(c); }
 37     template<class T> inline void wt(T a, T b, T c, T d) { ws(a), ws(b), ws(c), wn(d); }
 38     template<class T> inline T gcd(T a, T b) {
 39         return b == 0 ? a : gcd(b, a % b); }
 40     
 41 };
 42 using namespace my_header;
 43 const int maxS = 200 + 20;
 44 const int maxM = maxS * 2;
 45 const int mod = 10007;
 46 char str[maxS];
 47 int n, N, s, ans, f[maxS][maxS][maxS];
 48 inline void add(int &a, int b) {
 49     (a += b) %= mod;
 50 }
 51 inline void sub(int &a, int b) {
 52     (a -= b) %= mod;
 53     if (a < 0) a += mod;
 54 }
 55 int dp(int l, int r, int k) {
 56     int &ret = f[l][r][k];
 57     if (ret != -1)
 58         return ret;
 59     ret = 0;
 60     if (l > r)
 61         return 0;
 62     if (l == r)
 63         return ret = k == 0;
 64     if (str[l] == str[r]) {
 65         if (l == r - 1)
 66             return ret = k == 0;
 67         add(ret, dp(l + 1, r - 1, k));
 68     } else {
 69         add(ret, dp(l + 1, r, k - 1));
 70         add(ret, dp(l, r - 1, k - 1));
 71     }
 72     return ret;
 73 }
 74 struct Matrix {
 75     int d[maxM][maxM];
 76     int n;
 77     int* operator [] (int a) {
 78         return d[a];
 79     }
 80     void clear() {
 81         memset(d, 0, sizeof d);
 82     }
 83     void init() {
 84         memset(d, 0, sizeof d);
 85         for (int i = 1; i <= n; ++i)
 86             d[i][i] = 1;
 87     }
 88 };
 89 Matrix operator * (Matrix a, Matrix b) {
 90     Matrix c;
 91     c.n = a.n;
 92     c.clear();
 93     for (int i = 1; i <= a.n; ++i)
 94         for (int j = i; j <= a.n; ++j)
 95             for (int k = i; k <= j; ++k)
 96                 add(c[i][j], a[i][k] * b[k][j]);
 97     return c;
 98 }
 99 template<class T> inline T fpm(T b, int i) {
100     T r;
101     r.n = b.n;
102     r.init();
103     for (; i; i >>= 1, b = b * b)
104         if (i & 1)
105             r = r * b;
106     return r;
107 }
108 int main() {
109 #ifdef lol
110     freopen("E.in", "r", stdin);
111     freopen("E.out", "w", stdout);
112 #endif
113     cin >> (str + 1) >> n;
114     s = strlen(str + 1);
115     memset(f, -1, sizeof f);
116     for (int i = 0; i <= s; ++i)
117         f[1][s][i] = dp(1, s, i);
118     N = (s + n + 1) >> 1;
119     int n25 = (s + 1) / 2;
120     int num = s + n25;
121     Matrix x, res;
122     x.clear();
123     x.n = num + n25;
124     x[1][1] = 24;
125     for (int i = 2; i <= s; ++i) {
126         x[i][i] = 24;
127         x[i - 1][i] = 1;
128     }
129     
130     x[s][s + 1] = 1;
131     x[s + 1][s + 1] = 25;
132     x[s + 1][num + 1] = 1;
133     x[num + 1][num + 1] = 26;
134     for (int i = s + 2; i <= s + n25; ++i) {
135         x[i][i] = 25;
136         x[i - 1][i] = 1;
137         x[i][num + i - s] = 1;
138         x[num + i - s][num + i - s] = 26;
139     }
140     res = fpm(x, N);
141     for (int k = 0; k <= s - 1; ++k) {
142         int a = s - k + 1;
143         int b = s + n25 + (s - k + 1) / 2;
144         add(ans, dp(1, s, k) * res[a][b]);
145     }
146     if ((s + n) % 2 == 1) {
147         res = fpm(x, N - 1);
148         for (int k = 0; k <= s - 1; ++k)
149             if ((s - k) % 2 == 0) {
150                 int a = s - k + 1;
151                 int b = s + (s - k + 1) / 2;
152                 sub(ans, dp(1, s, k) * res[a][b]);
153             }
154     }
155     cout << ans << endl;
156     return 0;
157 }