DoorKickers的博客

博客

[AH2017/HNOI2017] 礼物

2020-04-19 12:22:55 By DoorKickers

首先, 两个手环增加非负整数亮度, 等价于其中一个增加一个整数亮度.

答案应该是 $$ \sum_{i = 1}^{n}(x_i - y_i + c)^2 $$ 展开 $$ \sum_{i = 1}^{n}(x_i^2 + y_i^2 - 2x_iy_i + 2x_ic - 2y_ic + c^2) \\= \sum_{i = 1}^{n}(x_i^2 + y_i^2 + c^2) + 2c\sum_{i = 1}^{n}{(x_i - y_i)} - 2\sum_{i = 1}^{n}x_iy_i $$ 然后发现我们只需要算最后一项. $$ \sum_{i = 1}^{n}x_iy_i $$ 设$x_0 = y_0 = 0$

考虑将$x$翻转, 然后左移一位 原式等于 $$ \sum_{i = 0}^{n}{x_{n - i}y_i} $$ 是个卷积.

考虑将x 破环成链(主要是卷积的性质)

假设

1 2 3 1 2 3 0 :x

0 4 5 6 0 0 0 :y

设结果为$z(x)$ $$ z_i = \sum_{k = 0}^{i}x_{i - k}y_k $$ 发现卷起来以后n ~ 2n - 1取个max就是答案. 代码:

#include <bits/stdc++.h>

using namespace std;

#define int long long

const int g = 3;
const int N = 1e6 + 10;
const int mod = 1004535809;

int omega[N], omegaInv[N], a[N], b[N], rev[N];

inline int fp(int a, int n, int mod)
{
    int res = 1;
    while (n)
    {
        if (n & 1) res = (res * a) % mod;
        n >>= 1;
        a = (a * a) % mod;
    }
    return res;
}

inline int inv(int x)
{
    return fp(x, mod - 2, mod);
}

inline void init(const int n)
{
    int k = 0;
    while ((1 << k) < n) k++;
    int x = fp(g, (mod - 1) / n, mod);
    omega[0] = omegaInv[0] = 1;
    for (int i = 1; i < n; i++)
    {
        rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << k - 1);
        omega[i] = (omega[i - 1] * x) % mod;
        omegaInv[i] = inv(omega[i]);
    }
}

inline void ntt(int *a, const int n, int *omega)
{
    for (int i = 0; i < n; i++)
    {
        if (i > rev[i]) continue;
        swap(a[i], a[rev[i]]);
    }
    for (int l = 2; l <= n; l *= 2)
    {
        int m = l / 2;
        for (int *p = a; p != a + n; p += l)
        {
            for (int i = 0; i < m; i++)
            {
                int t = (omega[n / l * i] * p[m + i]) % mod;
                p[m + i] = (p[i] - t + mod) % mod;
                p[i] = (p[i] + t) % mod;
            }
        }
    }
}

inline void dft(int *a, const int n)
{
    ntt(a, n, omega);
}

inline void idft(int *a, const int n)
{
    ntt(a, n, omegaInv);
    int x = inv(n);
    for (int i = 0; i < n; i++)
        a[i] = (a[i] * x) % mod;
}

inline void multiply(int *a1, const int n1, int *a2, const int n2)
{
    int n = 1;
    while (n < n1 + n2 + 2) n *= 2;
    init(n);
    dft(a1, n);
    dft(a2, n);
    for (int i = 0; i < n; i++)
        a1[i] = (a1[i] * a2[i]) % mod;
    idft(a1, n);
    idft(a2, n);
}

inline void squared(int *a, int n1)
{
    int n = 1;
    while (n < n1 * 2 + 2) n *= 2;
    init(n);
    dft(a, n);
    for (int i = 0; i < n; i++)
        a[i] = (a[i] * a[i]) % mod;
    idft(a, n);
}

inline void poly_print(int *a, int n)
{
    for (int i = 0; i <= n; i++)
        cout << a[i] << ' ';
    cout << endl;
}

inline void test()
{
    double m;
    cin >> m;
    int t = m;
    cout << "fk " << t << endl;
    int n = floor(m);
    cout << n << endl;
    // n = ceil(m);
    cout << n << endl;
    exit(0);
}


signed main()
{
    // freopen("input.txt", "r", stdin);
    // test();
    int n, m;
    cin >> n >> m;
    int X = 0;
    int Y = 0;
    int XX = 0;
    int YY = 0;
    for (int i = 0; i < n; i++)
    {
        cin >> a[i];
        X += a[i];
        XX += fp(a[i], 2, mod);
        a[i + n] = a[i];
    }
    reverse(a, a + n * 2);
    // for (int i = 0; i < n * 2; i++)
    // {
    //     cout << a[i] << ' ';
    // }
    // cout << endl;
    for (int i = 1; i <= n; i++)
    {
        cin >> b[i];
        Y += b[i];
        YY += fp(b[i], 2, mod);
    }
    int ans1 = ((Y - X) / n) + 1;
    // cout << "fk " << (double)((Y - X) / n);
    int ans2 = (Y - X) / n;
    int ans3 = (Y - X) / n - 1;
    // // cout << endl;
    // cout << "ans1 " << ans1 << endl;
    // cout << "ans2 " << ans2 << endl;
    // cout << "ans3 " << ans3 << endl;
    // cout << "x " << (Y - X) << endl;
    // cout << "n " << n << endl;
    // cout << fixed << setprecision(50) << (double)((Y - X) / n) << endl;
    // double t1 = Y - X;
    // double t2 = n;
    // cout << t1 / t2 << endl;
    // cout << "double " << (double)((Y - X) / n) << endl;
    // cout << (int)(ceil((double)((Y - X) / n))) << endl;
    int res1 = XX + YY;
    int res2 = min(ans1 * 2 * (X - Y) + n * ans1 * ans1, ans2 * 2 * (X - Y) + n * ans2 * ans2);
    res2 = min(res2, ans3 * 2 * (X - Y) + n * ans3 * ans3);
    int len = n * 2 - 1;
    multiply(a, len, b, len);
    int res3 = 0;
    for (int i = n; i < 2 * n; i++)
        res3 = max(res3, a[i]);
    cout << res1 + res2 - 2 * res3 << endl;
    return 0;
}

[SDOI2017]序列计数

2020-04-19 08:31:39 By DoorKickers

而且, 这题用NTT + FFT貌似不可行.

NTT模数不对, FFT爆精度.

首先想一个暴力算法:

设$f[i][j]$ 表示选$i$个, 在模$p$ 意义下是j的方案数. $$ f[i][j] = f[i - 1][k] \times l \quad (k + l) \bmod p = j $$ 复杂度$O(nm)$

然后发现一个优化, 可以只记录一下每个数模$p$的结果, 记录一下$cnt$

然后可以优化到差不多$O(np)$

然后又发现一个性质 $$ f[k] * f[j] = f[j + k] $$ 然后类比 $$ x^a * x^b = x^{a + b} $$ 如何快速计算$x^{a + b}$呢

快速幂就好了.

同理, 计算$f$ 也可以快速幂.

然后数组可以用滚动数组优化.

最终复杂度是$O(\log n * p)$ 再算上一些奇奇怪怪的预处理.

#include <bits/stdc++.h>

using namespace std;

// #define int long long
#define ll long long

const int N = 2e2 + 10;
const int L = 2e7 + 10;
const int mod = 20170408;

ll F[N], G[N], f[N], g[N], n, m, p;

int v[L], q[L / 10], cnt[N];


inline void multiply(ll *a, ll *b)
{
    static ll res[N];
    for (int i = 0; i < p; i++)
        for (int j = 0; j < p; j++)
        {
            res[i + j] = (a[i] * b[j] % mod + res[i + j]) % mod;
        }
    for (int i = 0; i < p; i++)
    {
        a[i] = (res[i] + res[i + p]) % mod;
        res[i] = res[i + p] = 0;
    }
}

signed main()
{
    cin >> n >> m >> p;
    F[1] = G[1] = f[1] = g[1] = 1;
    for (int i = 2; i <= m; i++)
    {
        f[i % p]++;
        g[i % p]++;
        F[i % p]++;
        G[i % p]++;
        if (v[i] == 0)
        {
            v[i] = i;
            // cout << i << ' ';
            G[i % p]--;
            g[i % p]--;
            q[++q[0]] = i;
        }
        for (int j = 1; j <= q[0]; j++)
        {
            if (q[j] > v[i] || q[j] * i > m) break;
            v[i * q[j]] = q[j];
        }
    }
    n--;
    while (n)
    {
        if (n & 1)
        {
            multiply(F, f);
            multiply(G, g);
        }
        n >>= 1;
        multiply(f, f);
        multiply(g, g);
    }
    cout << (F[0] - G[0] + mod) % mod << endl;
    return 0;
}

多项式入门推式子好题

2020-04-18 16:05:01 By DoorKickers

题目大意

给定$n$个数$q_1, q_2 \dots q_n$ 定义 $$ F_j = \sum_{i = 1}^{j - 1}\frac{q_i \times q_j}{(i - j)^ 2} - \sum_{i = j + 1}^{n}\frac{q_i \times q_j}{(i - j) ^ 2} $$

同时 $$ E_i = \frac{F_i}{q_i} $$

你的任务是对于所有的$1 \leq i \leq n$, 求出 $E_i$


Solution

推下公式: $$ E_j = \frac{F_j}{q_j} = \sum_{i = 1}^{j - 1}{\frac{q_i}{(i - j)^2}} - \sum_{i = j + 1}^{n}{\frac{q_i}{(i - j) ^2}} $$ 然后如果要用FFT加速的话, 要把这个推成卷积式

卷积式就是$C_k = \sum_{0}^{k}{A[i] * B[k - i]}$

继续推 $$ E_j = \sum_{i = 1}^{j}{\frac{q_i}{(i - j)^2}} - \sum_{i = j}^{n}{\frac{q_i}{(i - j) ^ 2}} \\ $$ 设$f[i] = q_i$, $g[i] = \frac{1}{i ^ 2}$

原式变成 $$ E_j = \sum_{i = 1}^{j}{f[i] * g[j - i]} - \sum_{i = j}^{n}{f[i] * g[i - j]} $$ 然后设$f[0] = g[0] = 0$ $$ E_j = \sum_{i = 0}^{j}f[i] * g[j - i] - \sum_{i = j}^{n}f[i] * g[i - j] $$ 发现左边已经是一个卷积式子了, 继续推右边的.

将他展开, 发现 $$ \sum_{i = j}^{n}f[i] * g[i - j] = \\ f[j] * g[0] + f[j + 1] * g[1] + \dots + f[j + (n - j)] * g[n - j] \\ = \sum_{i = 0}^{n - j}f[j + i] * g[i] $$ 然后发现$f$还不太行, 翻转一下.

设$f'[i] = f[n - i]$ $$ \sum_{i = 0}^{n - j}f[j + i] * g[i] = \sum_{i = 0}^{n - j}f'[n - (j + i)] * g[i] \\ = \sum_{i = 0}^{n - j}f'[n - j - i] * g[i] \\ $$ 然后设$t = n - j$

原式 $$ \sum_{i = 0}^{t}f'[t - i] * g[i] $$ 答案就是 $$ E_j = \sum_{i = 0}^{j}{f[i] * g[j - i]} - \sum_{i = 0}^{t}f[t - i] * g[i]\quad (t = n - j) $$ 然后设$A(x) = \sum_{i = 0}^{n}f[i], B(x) = \sum_{i = 0}^{n} g[i] , C(x) = \sum_{i = 0}^{n}f'[i]$

设$D(x) = A(x) * B(x), E(x) = B(x) * C(x)$

然后 $E_j = D(j) - E(n - j)$

就做完了

代码:

#include <bits/stdc++.h>

using namespace std;

#define int long long

const long double PI = 3.14159265358978323846264338;
const int N = 1e6 + 10;

complex<double> omega[N], omegaInv[N];
int rev[N];
double f[N], g[N], ff[N], res1[N], res2[N];

inline void init(const int n)
{
    int k = 0;
    while ((1 << k) < n) k++;
    for (int i = 0; i < n; i++)
    {
        rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (k - 1));
        omega[i] = complex<double>(cos(2 * PI / n * i), sin(2 * PI / n * i));
        omegaInv[i] = conj(omega[i]);
    }
    rev[0] = 0;
}

void fft(complex<double> *a, const int n, const complex<double> *omega)
{
    for (int i = 0; i < n; i++)
    {
        if (i > rev[i]) continue;
        swap(a[i], a[rev[i]]);
    }
    for (int l = 2; l <= n; l *= 2)
    {
        int m = l / 2;
        for (complex<double> *p = a; p != a + n; p += l)
            for (int i = 0; i < m; i++)
            {
                complex<double> t = omega[n / l * i] * p[m + i];
                p[m + i] = p[i] - t;
                p[i] += t;
            }
    }
}

inline void dft(complex<double> *a, const int n)
{
    fft(a, n, omega);
}

inline void idft(complex<double> *a, const int n)
{
    fft(a, n, omegaInv);
    for (int i = 0; i < n; i++)
    {
        a[i] /= n;
    }
}

inline void squared(int *a1, const int n1)
{
    int n = 1;
    while (n < n1 + n1 + 2) n *= 2;
    static complex<double> c1[N];
    for (int i = 0; i <= n1; i++) c1[i].real(a1[i]);
    init(n);
    dft(c1, n);
    for (int i = 0; i < n; i++)
        c1[i] *= c1[i];
    idft(c1, n);
    for (int i = 0; i <= n1 * 2; i++)
        a1[i] = static_cast<int>(floor(c1[i].real() + 0.5));
}

inline void multiply(double *a1, const int n1, double *a2, const int n2, double *res)
{
    int n = 1;
    while (n < n1 + n2 + 2) n *= 2;
    static complex<double> c1[N], c2[N];
    memset(c1, 0, sizeof(c1));
    memset(c2, 0, sizeof(c2));
    for (int i = 0; i <= n1; i++) c1[i].real(a1[i]);
    for (int i = 0; i <= n2; i++) c2[i].real(a2[i]);
    init(n);
    dft(c1, n);
    dft(c2, n);
    for (int i = 0; i < n; i++)
        c1[i] *= c2[i];
    idft(c1, n);
    // idft(c2, n);
    for (int i = 0; i <= n1 + n2; i++)
        res[i] = c1[i].real();
}

signed main()
{
    int n;
    cin >> n;
    for (int i = 1; i <= n; i++)
    {
        cin >> f[i];
        ff[n - i] = f[i];
        g[i] = (double)(1.0 / i / i);
    }
    // for (int i = 1; i <= n; i++)
    //     cin >> g[i];
    multiply(f, n, g, n, res1);
    // for (int i = 0; i <= n * 2; i++) 
    //     cout << fixed << setprecision(10) << res1[i] << ' ';
    // cout << endl;
    // exit(0);
    multiply(ff, n, g, n, res2);
    for (int i = 1; i <= n; i++)
        cout << fixed << setprecision(3) << res1[i] - res2[n - i] << endl;
    return 0;
}

SAM Bilion 老师的SAM 入门小题题

2020-04-10 16:15:55 By DoorKickers

首先 发现$L$满足单调性, 可以二分

考虑一个DP​

设$f_i$ 表示将前$i$个元素分成若干段最大有多少被匹配的

转移

$$ f_i = \max(f_j + i - j, f_{i - 1}) \quad (i - len(i) \leq j \leq i - L) $$

$len(i)$ 表示以$i$ 结尾最长能被匹配到的子串的长度, $L$是当前二分的值

如何求出以$i$结尾最长能被匹配的子串长度呢?

一个朴素的方法是对每个位置暴力走, 然后能走多少就是多少, 但复杂度太高

发现走的路程都是重复的, 可以记录一些信息来节省时间.

我们可以从位置1开始暴力匹配, 如果走不下去, 就跳后缀链接再匹配

这样可以充分利用之前的信息, 减少复杂度.

然后观察一下转移方程, 发现状态不能优化, 只能优化转移, 发现一个性质.

$$ i - len(i) \leq i + 1 - len(i + 1) \\i - L < i + 1 - L $$

发现决策的位置是递增的, 于是我们可以维护$f_j - j$ 的单调递减队列.

总体思路就是对$M$个串建广义SAM, 然后对每个询问串, 先预处理$len(i)$

然后二分答案, 用$dp$ check.

总复杂度为$O(\sum{|S|} + \sum{|P|\log|P|})$

代码:

#include <bits/stdc++.h>

using namespace std;

#define int long long

const int N = 1e6 + 10;

char s[N];
int tot, fa[N * 2], ch[N * 2][26], len[N * 2], mx[N * 2], f[N * 2];
pair<int, int> qu[N * 2];
inline int insert(int x, int las)
{
    if (ch[las][x] && len[las] + 1 == len[ch[las][x]])
        return ch[las][x];
    int cur = ++tot;
    int p = las;
    int clone = 0;
    int fl = 0;
    len[cur] = len[las] + 1;
    while (p && !ch[p][x])
    {
        ch[p][x] = cur;
        p = fa[p];
    }
    if (p == 0)
        fa[cur] = 1;
    else 
    {
        int q = ch[p][x];
        if (len[p] + 1 == len[q])
            fa[cur] = q;
        else 
        {
            if (p == las) fl = 1;
            clone = ++tot;
            len[clone] = len[p] + 1;
            fa[clone] = fa[q];
            memcpy(ch[clone], ch[q], sizeof(ch[q]));
            while (p && ch[p][x] == q)
            {
                ch[p][x] = clone;
                p = fa[p];
            }
            fa[q] = fa[cur] = clone;
        }
    }
    return fl ? clone : cur;
}

inline bool check(int k, int len)
{
    int res = 0;
    int l = 1, r = 0;
    qu[++r] = make_pair(0, 0);
    for (int i = 0; i < k; i++)
        f[i] = 0;
    for (int i = k; i <= len; i++)
    {
        while (l <= r && qu[l].second < i - mx[i])
            l++;
        if (l > r)
            f[i] = f[i - 1];
        else 
            f[i] = max(f[i - 1], qu[l].first + i);
        int t = i - k + 1;
        while (l <= r && qu[l].first <= f[t] - t)
            r--;
        qu[++r] = make_pair(f[t] - t, t);
        res = max(f[i], res);
    }
    double temp = res;
    double t = (double)(len) * 0.9;
    if (temp >= t)
        return true;
    else 
        return false;
}

signed main()
{
    int n, m;
    cin >> n >> m;
    tot = 1;
    for (int i = 1; i <= m; i++)
    {
        cin >> s + 1;
        int len = strlen(s + 1);
        int las = 1;
        for (int j = 1; j <= len; j++)
            las = insert(s[j] - '0', las);
    }
    for (int i = 1; i <= n; i++)
    {
        cin >> s + 1;
        int leng = strlen(s + 1);
        int cur = 1;
        int res = 0;
        for (int j = 1; j <= leng; j++)
        {
            int x = s[j] - '0';
            while (cur != 1 && !ch[cur][x])
            {
                cur = fa[cur];
                res = len[cur];
            }
            if (!ch[cur][x])
                mx[j] = 0;
            else 
            {
                res += 1;
                mx[j] = res;
                cur = ch[cur][x];
            }
        }
        int l = 0;
        int r = leng;
        while (l < r)
        {
            int mid = (l + r + 1) >> 1;
            if (check(mid, leng))
                l = mid;
            else 
                r = mid - 1;
        }
        cout << l << endl;
    }
    return 0;
}

SBzlt's Math note

2020-03-25 22:12:05 By DoorKickers

Binomial Coefficients

Let $\begin{pmatrix}r \\ k \end{pmatrix}$ to denote that choose $k$ elements from a set containing r terms without considering order.

With combinational interpretation. $$ \begin{pmatrix} r \\ k \end{pmatrix} = \frac{n(n - 1) \dots (n - k + 1)}{k(k - 1)\dots 1} = \frac{n!}{(n - k)!\ k!} $$

There are some identities that hold when r is a positive real number and k is non-negative integer.

The first is symmetry identity $$ {r \choose k} = {r \choose r - k} $$ Next, we got *absorption identity* $$ {r \choose k} = \frac{r}{k}{r - 1 \choose k - 1} $$ If we multiply both sides by k, then we have $$ k{r \choose k} = r{r - 1 \choose k - 1} $$ There is a companion identity of the above one. $$ \begin{align}(r - k){r \choose k} &= (r - k) {r \choose r - k} \\&=r{r - 1 \choose r - k - 1} \\&=r{r - 1 \choose k}\end{align} $$ One of the other important formula is called *addition formula* $$ {r \choose k} = {r - 1\choose k - 1} + {r - 1 \choose k} $$ And it's easy to prove.

Beside these, we also have two identities of summation form. $$ \sum_{k = -r}^{n}{{r + k \choose k}} = {r + n + 1 \choose n} $$ and $$ \sum_{k = 0}^{n}{k \choose m} = {n + 1 \choose m + 1} $$ Notice that the summation $0 + 1 + 2 + \dots + n$ is a special case of former identity when $m = 1$ , and its result is equal to $\frac{n(n+ 1)}{2}$ which is ${n + 1 \choose 2}$ manipulated by that formula.

Here is one of the most important identity in binomial coefficients. $$ (x + y)^r =\sum_{k}{{r \choose k}x^ky^{r - k}} $$ And when $x = y = 1$ and $r = n$ is nonnegative integer, we got $$ {n \choose 0} + {n \choose 1} + {n \choose 2} + \dots + {n \choose n} = 2^n $$ If we replace x by -1, then $$ {n \choose 0} - {n \choose 1} + {n \choose 2} - {n \choose 3} + \dots + (-1)^n{n \choose n} = 0 $$

Assumed that $y = 1$ and writing z instead of x to emphasize that an arbitrary complex number can be involved here. $$ (1 + z)^r = \sum_{k}{{r \choose k}z^k} $$

DoorKickers Avatar