@雅礼集训01/06 - T3@ math


@description@

给出 n, m, x,你需要求出下列式子的值:

\[\sum_{(\sum_{i=1}^mk_i)=n}(\prod_{i=1}^{m}\sin(k_i*x)) \]

其中 ki 为正整数。由于答案非常大,你只需要输出答案(保证不为 0)的正负(如果是负数输出负号,否则输出正号)和从左往右第一个非 0 数位上的数字即可。

input
第一行一个整数 T 表示数据组数。
对于每组数据,每行有两个整数 m,n 和一个两位小数 x。
对于 100%的数据,T ≤ 10,n ≤ 10^9,m ≤ 30。
output
输出共 T 行,每行两个字符表示答案。

sample input
2
3 5 0.01
3 6 0.02
sample output
+2
+4

@solution@

其实题目要你求解的说白了就是卷积。即对于多项式 \(f(p) = \sum_{i=1}^{+\infty}\sin(i*x)p^i\),求 \(f^m(p)\) 中第 n 项的系数为多少。然而 n 最大可以为 10^9,FFT 再怎么优化也过不了。

我们发现对于这个卷积,我们想要得到的其实只有第 n 项的系数,其他项的系数并不重要,而 FFT 必然要求出所有的系数,所以时间复杂度肯定降不下来。我们不得不换一种思路。
如果是一般多项式,FFT 就是时间复杂度的下限。但是我们的系数是三角函数,也就是说,我们要利用三角函数的一些性质。

数据范围 n <= 10^9,且对象是三角函数。如果熟悉的话,很容易想到三角函数的递推式:

\[\sin kx = 2\cos x\sin(k-1)x -\sin(k-2)x \]

边界为已知量 \(\sin 0, \sin x\)\(\cos x\) 是一个我们可以预先知道的常数,这样子就可以和矩阵乘法扯上关系了。

这个式子是怎么来的呢?可以参考下面的推导过程(但是如果你很熟悉了可以直接跳过)。
首先是三角函数的和差公式:

\[\sin kx = \sin x\cos (k-1)x+\cos x\sin (k-1)x \]

这个式子中的 \(\cos (k-1)x\) 是我们想要消去的,因此再对它使用一次和差公式。

\[\sin kx =\sin x\cos(k-2)x\cos x-\sin^2 x\sin(k-2)x + \cos x\sin (k-1)x \]

又出现了 \(\cos (k-2)x\)。但是,注意到其实可以把 \(\sin x\cos(k-2)x\) 当作一个整体,而这个整体出现在 \(\sin (k-1)x\) 的和差公式中。
因此,我们把 \(\sin (k-1)x = \sin x\cos (k-2)x+\cos x\sin (k-2)x\) 变形代入上式。

\[\sin kx =\cos x\sin(k-1)x -\cos^2 x\sin(k-2)x-\sin^2 x\sin(k-2)x + \cos x\sin (k-1)x \]

恒等变形,就可以得到我们的结果:

\[\sin kx =2\cos x\sin(k-1)x -\sin(k-2)x \]

好的我们推完式子再回到我们刚刚的思路。

考虑几种特殊情况吧。
当 m = 1 时,就是求 \(\sin nx\), 直接矩阵加速即可(当然最直接的还是暴力算)。
当 m = 2 时,令 \(g[i] = \sum_{j=1}^{i-1}(\sin(j*x)*\sin((i-j)*x))\),我们相当于是求 \(g[n]\)
我们尝试建立递推关系。代入三角函数的递推式得到(注意边界情况的存在):

\[g[i] = \sin x*\sin((i-1)*x)+\sum_{j=2}^{i-1}((2\cos x\sin(j-1)x -\sin(j-2)x)\sin((i-j)*x))\\ =\sin x*\sin((i-1)*x)+2*\cos x*g[i-1]-g[i-2]\]

其中 \(\sin((i-1)*x)\) 虽然不是常数,但是也可以通过矩阵乘法得到。

更一般地,令 \(dp[i][j]\) 表示 j 个多项式卷积第 i 项的系数,我们可以得到如下的递推关系:

\[\begin{cases}dp[i][j]=\sin x*dp[i-1][j-1]+2*\cos x*dp[i-1][j]-dp[i-2][j] & i \geq 2 \\ dp[i][j] = \sin x & i = 1 且 j = 1\\ dp[i][j] = 0 & otherwise \end{cases}\]

然后就可以矩阵加速了。

@accepted code@

不知道为什么,求解非零位时必须按照标程写才能过。
如果 p 一直作除法改成 ans 一直作乘法就过不了。
是因为后者浮点误差太大了吗?如果有知道的拜托请评论在下面告诉我好吗 QAQ。

#include<cstdio>
#include<cmath>
const int MAXN = 30;
struct matrix{
    int r, c;
    double m[MAXN*2 + 5][MAXN*2 + 5];
}M, B;
matrix operator * (matrix A, matrix B) {
    matrix C; C.r = A.r, C.c = B.c;
    for(int i=0;i<C.r;i++)
        for(int j=0;j<C.c;j++) {
            C.m[i][j] = 0;
            for(int k=0;k<A.c;k++)
                C.m[i][j] += A.m[i][k]*B.m[k][j];
        }
    return C;
}
/*
void Print(matrix M) {
    puts("");
    for(int i=0;i<M.r;i++) {
        for(int j=0;j<M.c;j++)
            printf("%lf ", M.m[i][j]);
        puts("");
    }
    puts("");
}
*/
matrix qpow(matrix A, int p) {
    matrix ret; ret.r = A.r, ret.c = A.c;
    for(int i=0;i<ret.r;i++)
        for(int j=0;j<ret.c;j++)
            ret.m[i][j] = (i == j);
    while( p ) {
        if( p & 1 ) ret = ret * A;
        A = A * A;
        p >>= 1;
    }
    return ret;
}
void solve() {
    int m, n; double x;
    scanf("%d%d%lf", &m, &n, &x);
    B.r = 2*m, B.c = 1, B.m[0][0] = 1;
    for(int i=0;i<2*m;i++)
        B.m[i][0] = 0;
    B.m[m-1][0] = sin(x);
    M.r = M.c = 2*m;
    for(int i=0;i<2*m;i++)
        for(int j=0;j<2*m;j++)
            M.m[i][j] = 0;
    for(int i=0;i<m;i++) {
        M.m[i][i] = 2*cos(x);
        M.m[i][m + i] = -1;
        M.m[m + i][i] = 1;
    }
    for(int i=2;i<=m;i++)
        M.m[i-2][i-1] = sin(x);
    //Print(M); Print(B);
    //Print(qpow(M, n-1)*B);
    double ans = (qpow(M, n-1)*B).m[0][0];
    putchar(ans > 0 ? '+' : '-');
    ans = fabs(ans);
    if( ans > 1 ) {
        while( ans >= 10 ) ans /= 10;
        printf("%d\n", int(ans));
    }
    else {
        double p = 0.1;
        while( p > ans ) p /= 10;
        printf("%d\n", int(ans/p));
    }
}
int main() {
    int T; scanf("%d", &T);
    for(int i=1;i<=T;i++) solve();
}

@details@

不得不说这道题还真的挺容易让人走错方向来着。
首先这是一个卷积形式,一开始肯定是思考能不能用 FFT(特别是像我一样最近才入多项式这个坑什么都想先来 FFT 一下)。
然后 m 个正整数的和为 n,又是一个组合数学的经典问题。这又是一个大坑。
然后 sin(x),一样是因为最近学了多项式,搞得我都想泰勒展开了……

当然如果你是神犇肯定不会像我一样去想上面那些错误思路而是一眼就秒出了这道题的正解。

posted @ 2019-01-09 13:50  Tiw_Air_OAO  阅读(219)  评论(0编辑  收藏  举报