1181 字
6 分钟
BZOJ3771: Triple

Description#

我们讲一个悲伤的故事。 从前有一个贫穷的樵夫在河边砍柴。 这时候河里出现了一个水神,夺过了他的斧头,说: “这把斧头,是不是你的?” 樵夫一看:“是啊是啊!” 水神把斧头扔在一边,又拿起一个东西问: “这把斧头,是不是你的?” 樵夫看不清楚,但又怕真的是自己的斧头,只好又答:“是啊是啊!” 水神又把手上的东西扔在一边,拿起第三个东西问: “这把斧头,是不是你的?” 樵夫还是看不清楚,但是他觉得再这样下去他就没法砍柴了。 于是他又一次答:“是啊是啊!真的是!” 水神看着他,哈哈大笑道: “你看看你现在的样子,真是丑陋!” 之后就消失了。 樵夫觉得很坑爹,他今天不仅没有砍到柴,还丢了一把斧头给那个水神。 于是他准备回家换一把斧头。 回家之后他才发现真正坑爹的事情才刚开始。 水神拿着的的确是他的斧头。 但是不一定是他拿出去的那把,还有可能是水神不知道怎么偷偷从他家里拿走的。 换句话说,水神可能拿走了他的一把,两把或者三把斧头。 樵夫觉得今天真是倒霉透了,但不管怎么样日子还得过。 他想统计他的损失。 樵夫的每一把斧头都有一个价值,不同斧头的价值不同。总损失就是丢掉的斧头价值和。 他想对于每个可能的总损失,计算有几种可能的方案。 注意:如果水神拿走了两把斧头a和b,(a,b)和(b,a)视为一种方案。拿走三把斧头时,(a,b,c),(b,c,a),(c,a,b),(c,b,a),(b,a,c),(a,c,b)视为一种方案。

Input#

第一行是整数N,表示有N把斧头。 接下来n行升序输入N个数字Ai,表示每把斧头的价值。

Output#

若干行,按升序对于所有可能的总损失输出一行x y,x为损失值,y为方案数。

Sample Input#

4
4
5
6
7

Sample Output#

4 1
5 1
6 1
7 1
9 1
10 1
11 2
12 1
13 1
15 1
16 1
17 1
18 1

样例解释#

11 有两种方案是4+7和5+6,其他损失值都有唯一方案,例如4=4,5=5,10=4+6,18=5+6+7.

HINT#

所有数据满足:Ai<=40000

题解#

一把,两把或者三把

首先写出生成函数F(x)=xa1+xa2+xa3+...+xanF(x) = x^{a_1}+x^{a_2}+x^{a_3}+...+x^{a_n}
两把就是F2(x)F^2(x)减去B(x)=x2a1+x2a2+x2a3+...+x2anB(x)=x^{2a_1}+x^{2a_2}+x^{2a_3}+...+x^{2a_n}
三把就是F3(x)F^3(x)减去3F(x)B(x)3*F(x)*B(x)加上C(x)=x3a1+x3a2+x3a3+...+x3anC(x)=x^{3a_1}+x^{3a_2}+x^{3a_3}+...+x^{3a_n} FFT优化乘法

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
inline int read()
{
    int x=0,f=1;char ch=getchar();
    while (ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
    while (ch>='0'&&ch<='9'){x=x*10+ch-'0';ch=getchar();}
    return x*f;
}
const double pi = acos(-1.);
struct Complex
{
    double x, y;
    Complex() { x = y = 0; }
    Complex(double a, double b) : x(a), y(b) {}
    Complex operator+(const Complex &a) { return Complex(x + a.x, y + a.y); }
    Complex operator-(const Complex &a) { return Complex(x - a.x, y - a.y); }
    Complex operator*(const Complex &a) { return Complex(x * a.x - y * a.y, x * a.y + y * a.x); }
    Complex operator*(const double a) { return Complex(x * a, y * a); }
    Complex Get() { return Complex(x, -y); }
};
double Inv;
const int MAXN = 5e5 + 5;
int N;
int rev[MAXN];
void FFt(Complex *a, int op)
{
    Complex wn, w, t;
    for (int i = 0; i < N; i++)
        if (i < rev[i])
            swap(a[i], a[rev[i]]);
    for (int k = 2; k <= N; k <<= 1)
    {
        wn = Complex(cos(pi / (k >> 1)), op * sin(pi / (k >> 1)));
        for (int j = 0; j < N; j += k)
        {
            w = Complex(1, 0);
            for (int i = 0; i < (k >> 1); i++, w = w * wn)
            {
                t = a[i + j + (k >> 1)] * w;
                a[i + j + (k >> 1)] = a[i + j] - t;
                a[i + j] = a[i + j] + t;
            }
        }
    }
    if (op == -1)
        for (int i = 0; i < N; i++)
            a[i] = a[i] * Inv;
}
Complex A[MAXN], B[MAXN], C[MAXN];
Complex tmp[MAXN];
long long ans[MAXN];
int m = 0;
int main()
{
    int n = read();
    for (int i = 1; i <= n; i++)
    {
        int x = read();
        ans[x]++;
        A[x].x = 1;
        B[2 * x].x = 1;
        C[3 * x].x = 1;
        m = max(m, 3 * x);
    }
    m = m + m + 1;
    N = 1;
    while (N < m) 
        N <<= 1;
    for (int i = 1; i < N; i++) 
        if (i & 1)
            rev[i] = (rev[i >> 1] >> 1) | (N >> 1);
        else
            rev[i] = (rev[i >> 1] >> 1);
    Inv = 1. / N;
    FFt(A, 1);
    for (int i = 0; i < N; i++) tmp[i] = A[i] * A[i];
    FFt(tmp, -1);
    for (int i = 0; i < N; i++) ans[i] += round((tmp[i].x - B[i].x) / 2);
    // for (int i = 0; i < N; i++)
    //     if (ans[i])
    //         printf ("%d %lld\n", i, ans[i]);
    // printf ("======================================\n");
    FFt(B, 1);
    for (int i = 0; i < N; i++) B[i] = B[i] * A[i];
    FFt(B, -1);
    for (int i = 0; i < N; i++) tmp[i] = A[i] * A[i] * A[i];
    FFt(tmp, -1);
    for (int i = 0; i < N; i++) ans[i] += round((tmp[i].x - B[i].x * 3 + C[i].x * 2) / 6);
    for (int i = 0; i <= m / 2; i++)
        if (ans[i])
            printf ("%d %lld\n", i, ans[i]);
    // while (1);
}
BZOJ3771: Triple
https://www.nekomio.com/posts/128/
作者
NekoMio
发布于
2018-02-26
许可协议
CC BY-NC-SA 4.0