matrix fast pow

来源:互联网 发布:拍小成本电影知乎 编辑:程序博客网 时间:2024/06/07 20:07
ll mod=1e9+7;     //改modstruct matrix{    int n,m;    ll ma[15][15];    matrix(int x,int y):n(x),m(y)    {        memset(ma,0,sizeof(ma));    }    matrix operator*(matrix& x)    {        matrix res(n,x.m);        for(int i=1; i<=n; i++)            for(int j=1; j<=x.m; j++)                for(int k=1; k<=m; k++)                    res.ma[i][j]+=ma[i][k]*x.ma[k][j]%mod+mod,                                  res.ma[i][j]%=mod;        return res;    }    void pf(void)    {        for(int i=1; i<=n; i++)            for(int j=1; j<=m; j++)            {                printf("%d ",ma[i][j]);                if(j==m)printf("\n");            }    }};ll pv[105]= {0,1,1,0,1,1,2,1,3,2,5};  //改前几项, 前面用一个0占位const int num=10;                          //前几项个数double fn[num+1],wa[num+1];int bm(void){    int f=1,w=1,last=1;    fn[1]=0,wa[1]=1.0/pv[1];    for(int i=2; i<=num; i++)    {        double sum=-pv[i];        for(int j=1; j<=f; j++)            sum+=fn[j]*pv[i-j];        if(fabs(sum)<eps)continue;        double t[num+1]= {0,-1/sum};        for(int j=1; j<=f; j++)            t[j+1]=fn[j]/sum;        for(int j=1; j<=w; j++)            fn[i-last+j-1]+=-sum*wa[j];        memcpy(wa,t,sizeof(wa));        last=i,f++,w++;    }    while(fabs(fn[f])<eps)f--;    return f;}ll qpow(matrix x,matrix res,ll y){    for(; y; y>>=1,x=x*x)        if(y&1)res=res*x;    return res.ma[1][1];}bool flag=false;matrix res(0,0),x(0,0);ll cal(ll y){    if(y<=num)return pv[y]%mod;    if(!flag)    {        int siz=bm();        res=matrix(siz,siz);        for(int i=1; i<=siz; i++)        {            if(fn[i]<0)fn[i]-=eps;            else fn[i]+=eps;            res.ma[i][1]=(int)fn[i];            if(i<siz)res.ma[i][i+1]=1;        }        for(int i=1; i<=siz; i++)            res.ma[i][i+1]=1;        x=matrix(1,siz);        for(int i=num,j=1; i>num-siz; i--)            x.ma[1][j++]=pv[i];        flag=true;    }    return qpow(res,x,y-num);}int main(){    ll n;    while(~scanf("%lld",&n))               //改读入格式        printf("%lld\n",cal(n));    return 0;}
原创粉丝点击