鬼故事 题解

· · 个人记录

T3鬼故事 题解

好像有人没有用矩阵快速幂AC了……这里只是参考解答,非最优解。

题意:定义一个新数列bb_i=a_ia_{i+1}...a_{i+K-1},也就是从第i项开始连续K个斐波拉契数列项之积,求第M个到第Nb_i的和。

对于5\%数据(data\space1):

暴力枚举,不用多说

对于10\%数据(data\space1,2):

注意到这一部分数据中,N-M<=10^6,所以我们可以先用矩阵快速幂得到斐波拉契数列的第M项,然后再从MN枚举得到答案。

对于20\%数据(data\space1,2,5,6):

## 对于另外一种$20\%$数据($data\space1,2,3,4$): 易知此题需要用矩阵快速幂做,问题就在于如何构造矩阵。 题目要求的即$s_N-s_{M-1}$,其中$s_n=\sum_{i=1}^{n}b_i$。问题转化为如何求$s_n$。 因为$K=4$,所以$s_n=s_{n-1}+a_na_{n+1}a_{n+2}a_{n+3}$。由于出题人最开始的设定是$S_n=S_{n-1}+a_na_{n-1}a_{n-2}a_{n-3}$,为了方便出题人写题解,我们在这里求$S_n$($s_n=S_{n+3}$,很容易转换)。 有$S_n=S_{n-1}+a_na_{n-1}a_{n-2}a_{n-3}$,而学过矩阵快速幂的同学肯定知道$S_{n-1}$很容易处理,关键问题就在于如何处理$a_na_{n-1}a_{n-2}a_{n-3}$。 因为矩阵快速幂中的单位矩阵只包含常数,所以我们要想办法把$a_na_{n-1}a_{n-2}a_{n-3}$表示成若干个**多项式$\times$常数**之和的形式。 将$a_na_{n-1}a_{n-2}a_{n-3}$拆分。$a_na_{n-1}a_{n-2}a_{n-3}=(a_{n-1}+a_{n-2})a_{n-1}a_{n-2}a_{n-3}=a_{n-1}^2a_{n-2}a_{n-3}+a_{n-1}a_{n-2}^2a_{n-3}$。同理,将$a_{n-1}$拆成$a_{n-2}+a_{n-3}$:$a_{n-1}^2a_{n-2}a_{n-3}+a_{n-1}a_{n-2}^2a_{n-3}=2a_{n-2}a_{n-3}^3+3a_{n-2}^2a_{n-3}^2+a_{n-2}a_{n-3}^3$。 为了方便,我们设$b_n=a_n^3a_{n-1},c_n=a_n^2a_{n-1}^2,d_n=a_na_{n-1}^3$,则$a_na_{n-1}a_{n-2}a_{n-3}=2b_{n-2}+3c_{n-2}+d_{n-2}$。 既然有了上面这个式子,我们就要想办法把$b_{n-2},c_{n-2},d_{n-2}$表示出来。首先是$b_{n-2}$: $$b_{n-2}=a_{n-2}^3a_{n-3}=(a_{n-3}+a_{n-4})^3a_{n-3}$$ $$=a_{n-3}^4+3a_{n-3}^3a_{n-4}+3a_{n-3}^2a_{n-4}^2+a_{n-3}a_{n-4}^3$$ 为了方便,我们设$e_n=a_n^4$,则$b_{n-2}=e_{n-3}+3b_{n-3}+3c_{n-3}+d_{n-3}$。我们发现这个式子是可以用矩阵乘法实现的! 同理有: $$c_{n-2}=e_{n-3}+2b_{n-3}+c_{n-3}$$ $$d_{n-2}=e_{n-2}+b_{n-3}$$ $$e_{n-2}=e_{n-3}+4b_{n-3}+6c_{n-3}+4d_{n-3}+e_{n-4}$$ 至此,$b_{n-2},c_{n-2},d_{n-2},e_{n-2}$就全部被我们表示出来了。所以,我们发现我们可以构造这样的一个矩阵: $$\begin{pmatrix}1&?&?&?&?&0\\0&3&3&1&1&0\\0&2&1&0&1&0\\0&1&0&0&1&0\\0&4&6&4&1&1\\0&0&0&0&1&0\end{pmatrix}\times \begin{pmatrix}S_{n-1}\\b_{n-3}\\c_{n-3}\\d_{n-3}\\e_{n-3}\\e_{n-4}\end{pmatrix}=\begin{pmatrix}S_{n}\\b_{n-2}\\c_{n-2}\\d_{n-2}\\e_{n-2}\\e_{n-3}\end{pmatrix}$$ 这样,我们就可以达到矩阵快速幂可以进行的要求了。 还有一个问题,就是那几个$?$填什么。因为$S_n=S_{n-1}+a_na_{n-1}a_{n-2}a_{n-3}=S_{n-1}+2b_{n-2}+3c_{n-2}+d_{n-2}$,但我们没有直接的$b_{n-2},c_{n-2},d_{n-2}$,只有$b_{n-3},c_{n-3},d_{n-3},e_{n-3}$,但是我们有之前推出的式子可以表示$b_{n-2},c_{n-2},d_{n-2}$与$b_{n-3},c_{n-3},d_{n-3},e_{n-3}$的关系,所以我们将这几个式子代入回$S_n=S_{n-1}+2b_{n-2}+3c_{n-2}+d_{n-2}$,得到$S_n=S_{n-1}+13b_{n-3}+9c_{n-3}+2d_{n-3}+6e_{n-3}$。所以我们就得到了一个完整的式子: $$\begin{pmatrix}1&13&9&2&6&0\\0&3&3&1&1&0\\0&2&1&0&1&0\\0&1&0&0&1&0\\0&4&6&4&1&1\\0&0&0&0&1&0\end{pmatrix}\times \begin{pmatrix}S_{n-1}\\b_{n-3}\\c_{n-3}\\d_{n-3}\\e_{n-3}\\e_{n-4}\end{pmatrix}=\begin{pmatrix}S_{n}\\b_{n-2}\\c_{n-2}\\d_{n-2}\\e_{n-2}\\e_{n-3}\end{pmatrix}$$ 所以 $$\begin{pmatrix}1&13&9&2&6&0\\0&3&3&1&1&0\\0&2&1&0&1&0\\0&1&0&0&1&0\\0&4&6&4&1&1\\0&0&0&0&1&0\end{pmatrix}^{n-4}\times \begin{pmatrix}S_{4}\\b_{2}\\c_{2}\\d_{2}\\e_{2}\\e_{1}\end{pmatrix}=\begin{pmatrix}S_{n}\\...\\...\\...\\...\\...\end{pmatrix}$$ 回到一开始,$s_n=S_{n+3}$,所以: $$\begin{pmatrix}1&13&9&2&6&0\\0&3&3&1&1&0\\0&2&1&0&1&0\\0&1&0&0&1&0\\0&4&6&4&1&1\\0&0&0&0&1&0\end{pmatrix}^{n-1}\times \begin{pmatrix}s_{1}\\b_{2}\\c_{2}\\d_{2}\\e_{2}\\e_{1}\end{pmatrix}=\begin{pmatrix}s_{n}\\...\\...\\...\\...\\...\end{pmatrix}$$ 即: $$\begin{pmatrix}1&13&9&2&6&0\\0&3&3&1&1&0\\0&2&1&0&1&0\\0&1&0&0&1&0\\0&4&6&4&1&1\\0&0&0&0&1&0\end{pmatrix}^{n-1}\times \begin{pmatrix}6\\1\\1\\1\\1\\1\end{pmatrix}=\begin{pmatrix}s_{n}\\...\\...\\...\\...\\...\end{pmatrix}$$ ## 对于$30\%$数据($data\space1,2,3,4,5,6$): 加一个高精度就可以了。 ## 对于$50\%$数据($data\space1,2,3,4,5,6,7,8,9,10$): 压位高精。压$17$位可以通过此部分数据。 ## 对于$40\%$数据($data\space1,2,3,4,5,6,13,14$): 分$9$种情况讨论…… 对于每一种情况,求出对应的矩阵。 **从出题人心理角度入手:** $K$肯定越大越好,所以只用计算$K=8,9,10$的情况即可。 用普通高精度即可(其实$data\space11,12$)用常数小一点的不压位高精度也可以通过。 ## 对于$70\%$数据($data\space1\sim14$): 用压位高精分$9$种情况讨论。 # 对于$100\%$数据: ~~分49类讨论~~ 因为分类讨论过于复杂,所以我们要找到一个一般的求单位矩阵算法。 我们发现,只有$a_n^m,a_n^{m-1}a_{n-1}^1,a_n^{m-2}a_{n-1}^2,...,a_n^1a_{n-1}^{m-1}$与它们的上一项有互推关系,所以我们要用这些项来表示$s_n$,就如之前我们求$K=4$时问号填什么一样。 ## _如果这一部分理解有困难,可以结合$K=4$的部分一起看_ **第一步**:求$a_n^m,a_n^{m-1}a_{n-1}^1,a_n^{m-2}a_{n-1}^2,...,a_n^1a_{n-1}^{m-1}$与它们的上一项的递推式 例如,我们要求$a_n^{k}a_{n-1}^{m-k}$的递推式。按照套路,有$a_n^{k}a_{n-1}^{m-k}=(a_{n-1}+a_{n-2})^ka_{n-1}^{m-k}$。由二项式定理,上式等于$C_k^0a_{n-1}^m+C_k^1a^{m-1}_{n-1}a_{n-2}+...+c_k^ka_{n-1}^{m-k}a_{n-2}^k$。 于是,我们维护一个杨辉三角,对于每一个$a_n^{k}a_{n-1}^{m-k}$单独计算即可。我们称这些求出来的递推式为**那一堆式子**……(如$K=4$中的$b_{n-2}=e_{n-3}+3b_{n-3}+3c_{n-3}+d_{n-3},c_{n-2}=e_{n-3}+2b_{n-3}+c_{n-3},d_{n-2}=e_{n-2}+b_{n-3},e_{n-2}=e_{n-3}+4b_{n-3}+6c_{n-3}+4d_{n-3}+e_{n-4}$) **第二步**:求$a_na_{n-1}...a_{n-k+1}$与$a_n^m,a_n^{m-1}a_{n-1}^1,a_n^{m-2}a_{n-1}^2,...,a_n^1a_{n-1}^{m-1}$的关系式。 同样我们拆分$a_na_{n-1}...a_{n-k+1}$。在此之前,我们先证明一个东西。 求证:$a_n=a_{k+1}a_{n-k}+a_{k}a_{n-(k+1)}

证明:用数学归纳法证明

  1. k=1时,即证a_n=1a_{n-1}+1a_{n-2},显然成立

  2. 假设k=s(s>=1)时,a_n=a_ka_{n-k}+a_{k+1}a_{n-(k+1)}成立,则a_{s+1}a_{n-s}+a_{s}a_{n-(s+1)}=a_{s+1}(a_{n-(s+1)}+a_{n-(s+2)})+a_{s}a_{n-(s+1)}=a_{n-(s+1)}(a_{s+1}+a_s)+a_{n-(s+2)}a_{s+1}=a_{s+2}a_{n-(s+1)}+a_{s+1}a_{n-(s+2)},即k=s+1时命题成立

  3. 综上:命题成立

既然有这个结论,我们就可以得到a_na_{n-1}...a_{n-k+1}=(a_{k-1}a_{n-k+2}+a_{k-2}a_{n-k+1})(a_{k-2}a_{n-k+2}+a_{k-3}a_{n-k+1})...(a_2a_{n-k+2}+a_1a_{n-k+1})a_{n-k+2}a_{n-k+1})。观察这个式子,发现它是由很多个类似于(a_{k-1}a_{n-k+2}+a_{k-2}a_{n-k+1})的多项式组成的,我们称这样的括号为一个单位多项式,a_{k-1},a_{k-2}为一个单位多项式的两个系数。

现在,我们的目标就是求出这个式子展开后每一项的系数是多少。我们考虑用类似于求杨辉三角的方法求每一项的系数。

K=3,则原式=(a_{n-1}+a_{n-2})a_{n-1}a_{n-2},即a_{n-1}^2a_{n-2}+a_{n-1}a_{n-2}^2,即前面的系数为1,1K=4的我们之前已经推导过了,系数为2,3,1K=5的为6,13,9,2。我们把这个三角列出来:

0,1,1(1,1) 0,2,3,1(2,1) 0,6,13,9,2(3,2)

我们在每行前面补了一个0,后面括号中的数表示这一层所对应的那一个单位多项式的两个系数,我们设第i层的两个系数分别为p1_i,p2_i。我们发现k_{i,j},即第i层第j个数等于k_{i-1,j}p1_i+k_{i-1,j-1}p2_i。所以,我们就可以用类似于求杨辉三角的方式来求出每一项的系数了。

具体来说,我们可以从a_{n-1}...a_{n-k+1}=(a_{k-2}a_{n-k+2}+a_{k-3}a_{n-k+1})...(a_2a_{n-k+2}+a_1a_{n-k+1})a_{n-k+2}a_{n-k+1})的系数推到a_na_{n-1}...a_{n-k+1}=(a_{k-1}a_{n-k+2}+a_{k-2}a_{n-k+1})(a_{k-2}a_{n-k+2}+a_{k-3}a_{n-k+1})...(a_2a_{n-k+2}+a_1a_{n-k+1})a_{n-k+2}a_{n-k+1})的系数。后者比前者多了一个(a_{k-1}a_{n-k+2}+a_{k-2}a_{n-k+1}),而a_{k-1},a_{k-2}恰好就是三角里的p1_i,p2_i。于是,类似于杨辉三角的思路,每一个系数只与其上一层同样位置的数与上一层前一位的数有关系;而由于新乘上了p1_i,p2_i,每一个系数就等于上一层同样位置的数的p2_i倍与上一层前一位的数的p1_i倍之和(结合上面列出的三角形与K=4的情况来理解)。

强调一下,现在求出来的系数都是带在a_n^m,a_n^{m-1}a_{n-1}^1,a_n^{m-2}a_{n-1}^2,...,a_n^1a_{n-1}^{m-1}之前的。

系数求出来了,但还没有完。是否记得K=4时,将2,3,1求出来之后,还要分别把那一堆式子代入回原式?没错,求出来的系数还要分别乘上那一堆式子的系数,而那一堆式子我们已经求出来了。所以a_na_{n-1}...a_{n-k+1}a_n^m,a_n^{m-1}a_{n-1}^1,a_n^{m-2}a_{n-1}^2,...,a_n^1a_{n-1}^{m-1}的关系我们就求出来了。于是s_na_n^m,a_n^{m-1}a_{n-1}^1,a_n^{m-2}a_{n-1}^2,...,a_n^1a_{n-1}^{m-1}的关系我们也求出来了。注意,这里写的a_n^m,a_n^{m-1}a_{n-1}^1,a_n^{m-2}a_{n-1}^2,...,a_n^1a_{n-1}^{m-1}都有可能表示的是所有单项式下角标整体减12对应的一串另外的单项式,就比如把a_{n-2}a_{n-3}写成a_{n}a_{n-1},因为出题人太懒啦,只会复制黏贴!

然后矩阵快速幂,加一个压17位的高精度,这一题就做完了。是不是很简单啊!

如果仍有不理解的地方,可以结合代码理解。

当然还可以用一些毒瘤的方法来优化高精度矩阵快速幂,比如什么10进制的矩阵快速幂之类的(但是常数像\pi的位数一样大的出题人打的10进制的矩阵快速幂TLE了哦~)。

请大家注意常数优化!

标程(这是一份加了一些常数优化的代码,但仍然跑的非常慢):

#include<cstdio>
#include<algorithm>
#include<cstring>
#include<string>
#include<cmath>
#include<iostream>
using namespace std;
const long long mod=1000000007,gjm=1e18;
int k;
struct GJD
{
    long long shu[585];
    int len;
    GJD(){len=0;memset(shu,0,sizeof(shu));}
}m,n,lcgjd;
struct JZ
{
    long long shu[53][53];
    JZ(){memset(shu,0,sizeof(shu));}
}dw,lcjz,s,lcjzjz;
void mem1()
{
    memset(s.shu,0,sizeof(s.shu));
    for(register int i=1;i<=110;i++)s.shu[i][i]=1;
}
void times(JZ *a,JZ *b)
{
    for(register int i=1;i<=k+2;i++)
    {
        for(register int j=1;j<=k+2;j++)
        {
            lcjzjz.shu[i][j]=0;
            for(register int s=1;s<=k+2;s++)
            {
                lcjzjz.shu[i][j]=(lcjzjz.shu[i][j]+(*a).shu[i][s]*(*b).shu[s][j])%mod;
            }
        }
    }
    for(register int i=1;i<=k+2;i++)for(register int j=1;j<=k+2;j++)(*a).shu[i][j]=lcjzjz.shu[i][j];
}
long long x;
GJD d2(GJD &a)
{
    memset(lcgjd.shu,0,sizeof(lcgjd));
    lcgjd.len=a.len;x=0;
    for(register int i=a.len-1;i>=0;i--)
    {
        x=x*gjm+a.shu[i];
        lcgjd.shu[i]=x/2;
        x&=1;
    }
    while(lcgjd.len&&lcgjd.shu[lcgjd.len-1]==0)lcgjd.len--;
    return lcgjd;
}
void j1(GJD &a)
{
    for(register int i=0;i<a.len;i++)
    {
        a.shu[i]--;
        if(a.shu[i]>=0)break;
        a.shu[i]+=gjm;
    }
    while(a.len&&a.shu[a.len-1]==0)a.len--;
}
char lcm[10000],lcn[10000];
long long lcshu,lcshushu;
long long cmul[200]={0},fb[200]={0},exs[200][200]={0};
long long cal(GJD &a)
{
    if(!a.len)return 0;
    j1(a);
    lcjz=dw;mem1();
    for(;a.len;a=d2(a))
    {
        if(a.shu[0]&1)times(&s,&lcjz);
        times(&lcjz,&lcjz);
    }
    long long ans=0,b=1;
    for(register int i=1;i<=k;i++)b=b*fb[i]%mod;
    for(register int i=1;i<=k+1;i++)ans=(ans+s.shu[k+2][i])%mod;
    ans=(ans+s.shu[k+2][k+2]*b)%mod;
    return ans;
}
int main()
{
    scanf("%d%s%s",&k,lcm,lcn);
    int lclen=strlen(lcm);lcshu=0;lcshushu=1;
    for(register int i=lclen-1;i>=0;i--)
    {
        lcshu=lcshu+((long long)(lcm[i]-'0'))*lcshushu;
        lcshushu*=10;
        if((lclen-i)%18==0)
        {
            m.shu[m.len]=lcshu;
            lcshu=0;lcshushu=1;
            m.len++;
        }
    }
    if(lcshu){m.shu[m.len]=lcshu;m.len++;}
    lclen=strlen(lcn);lcshu=0;lcshushu=1;
    for(register int i=lclen-1;i>=0;i--)
    {
        lcshu=lcshu+((long long)(lcn[i]-'0'))*lcshushu;
        lcshushu*=10;
        if((lclen-i)%18==0)
        {
            n.shu[n.len]=lcshu;
            lcshu=0;lcshushu=1;
            n.len++;
        }
    }
    if(lcshu){n.shu[n.len]=lcshu;n.len++;}
    j1(m);
    cmul[1]=1;fb[1]=1;
    for(register int i=2;i<=k-1;i++)
    {
        fb[i]=(fb[i-1]+fb[i-2])%mod;
        for(register int j=i;j>=1;j--)
        {
            cmul[j]=(cmul[j]*fb[i]+cmul[j-1]*fb[i-1])%mod;
        }
    }
    fb[k]=(fb[k-1]+fb[k-2])%mod;
    exs[1][1]=1;
    for(register int i=2;i<=53;i++)
    {
        for(register int j=1;j<=i;j++)
        {
            exs[i][j]=(exs[i-1][j]+exs[i-1][j-1])%mod;
        }
    }
    for(register int i=1;i<=k;i++)
    {
        for(register int j=1;j<=k+1;j++)
        {
            dw.shu[i][j]=exs[k-i+2][j];
        }
        if(i==1)continue;
        for(register int j=1;j<=k+1;j++)
        {
            dw.shu[k+2][j]=(dw.shu[k+2][j]+dw.shu[i][j]*cmul[i-1])%mod;
        }
    }
    dw.shu[k+1][1]=1;dw.shu[k+2][k+2]=1;
    printf("%lld\n",((cal(n)-cal(m))%mod+mod)%mod);
    return 0;
}

(据说原来的标程被卡常了……)

PS:

应该也可以用斐波拉契数列模1000000007的循环节(长度为2000000016)来做。怎么算出来的?慢慢枚举就可以了,最多20s。记得不是直接取模就可以了,还要乘上中间被忽略的几大块。这样可以直接忽略高精时间复杂度的影响,应该会跑的飞快,大家可以试一试。