[LOJ6257][CodePlus 2017 12 月赛]可做题2-矩阵快速幂-乘法逆元

可做题2

题目描述

“codeplus比赛的时候在做什么？有没有空？能来解决丢番图方程问题吗？”sublinekelzrip这样问qmqmqm。

当然，qmqmqm并不会丢番图方程问题，所以sublinekelzrip改为提出了另一个题目，现在请你帮助qmqmqm解决这个题目。

这个问题是这样的：

若一个数列aaa满足条件$a_n=a_{n-1}+a_{n-2},n \geq 3$,而$a_1,a_2$为任意实数，则我们称这个数列为广义斐波那契数列。

现在请你求出满足条件$a_1=i$,$a_2$为区间[l,r][l,r][l,r]中的整数，且$a_k\ mod\ p=m$的广义斐波那契数列有多少个。

输入格式

从标准输入读入数据。

本题包含多组数据，输入第一行包含一个正整数TTT，表示数据组数。对于每组数据：

一行六个用空格隔开的整数$i,l,r,k,p,m$，意义如题目描述所示。

输出格式

输出到标准输出。

输出共$T$行，每行一个数表示该组数据的答案。

样例输入

6
2 17 68 3 23 1
1 17 68 3 57 1
5 17 68 10 11 9
5 17 68 10 71 9
10 17 68 11 12 3
10 17 68 8 6 4

样例输出

3
1
4
1
5
9

数据范围与提示

$k ,r \leq 10^{18}$

对于所有数据，$0 \leq l \leq r$$1 \leq p \leq 10^9$,$0 \leq m < p$,$T=10$,$0 \leq i \leq 10^{18}$,$k \geq 3$。

---

看题花了1min想到，1h打完初稿，怎么拍都拍不错却一直WA。
于是一直WA了整整两个小时一直到下考，直到考试结束了还是没有过。
结果晚上看数据发现

i!要!取!模!!!

(╯‵□′)╯︵┻━┻

T2明明一眼就看出来了却没有写……

---

思路:
很坑（没错咱说的就是那个$i$），但确实不是一道难题。
首先找规律，设$a_1=i,a_2=j$
那么有
$a_3=i+j$
$a_4=i+2j$
$a_5=2i+3j$
$\cdot \cdot\cdot$
观察$i$和$j$的系数:

k= 1 2 3 4 5 6 7 8
i  1 0 1 1 2 3 5 7
j  0 1 1 2 3 5 7 9

可以发现这些系数就是斐波那契数列。
那么用矩阵快速幂加速求斐波那契数，令系数分别为$a$和$b$，那么问题变成:
$a*i+b*j==m\ (mod\ p)$
$a*i$是个定值，那么变成
$b*j == m-a*i \ (mod\ p)$
然后就是求这个方程在$[l,r]$之间的解的个数了~
先考虑令场上的所有数除掉$gcd(b,p)$，若$m-a*i$不整除$gcd(b,p)$，则无解。
然后求逆元得:
$j==\frac{m-a*i}{b} \ (mod\ p)$
统计一下，然后就没有了~

注意一定要首先让$i$对$p$取模!!!

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<algorithm>

using namespace std;

typedef long long ll;

ll xx,l,r,k,p,m;

namespace fibo
{
    struct mat
    {
        ll a[4][4],lena,lenb;

        inline void init(){memset(a,0,sizeof(a));}
        mat(ll a=0,ll b=0){lena=a;lenb=b;init();}
        inline void e()
        {
            init();for(int i=0;i<lena;i++)a[i][i]=1;
        }
        inline mat operator * (mat o)
        {
            mat ret(lena,o.lenb);
            for(int i=0;i<lena;i++)
                for(int j=0;j<o.lenb;j++)
                {
                    for(int k=0;k<lenb;k++)
                        ret.a[i][j]+=a[i][k]*o.a[k][j]%p;
                    ret.a[i][j]%=p;
                }
            return ret;
        }
    };

    inline mat mqpow(mat a,ll b)
    {
        mat ret(a.lena,a.lenb);ret.e();
        while(b)
        {
            if(b&1)
                ret=ret*a;
            a=a*a;
            b>>=1;
        }
        return ret;
    }

    inline ll fib(ll x,ll a1=1,ll a2=2)
    {
        if(x==0)return 1;
        if(x==-2)return 1;
        if(x<0)return 0;
        mat base(2,2),mul(2,2);
        base.a[0][0]=a1;
        base.a[0][1]=a2;
        mul.a[0][1]=mul.a[1][0]=mul.a[1][1]=1;
        base=base*mqpow(mul,x-1);
        return base.a[0][0];
    }
}

using namespace fibo;

namespace invs
{
    inline ll gcd(ll a,ll b){if(!b)return a;return gcd(b,a%b);}
    inline ll exgcd(ll a,ll b,ll &x,ll &y)
    {
        if(!b){x=1;y=0;return a;}
        ll ret=exgcd(b,a%b,y,x);
        y=y-a/b*x;
        return ret;
    }
    inline ll inv(ll a,ll b)
    {
        ll x,y;
        exgcd(a,b,x,y);
        return (x%b+b)%b;
    }
}

using namespace invs;

inline ll calc(ll pos,ll v)
{
    return pos/p+(pos%p>=v%p);
}

inline void mian()
{
    scanf("%lld%lld%lld%lld%lld%lld",&xx,&l,&r,&k,&p,&m);
    m=((m-fib(k-3)*(xx%p)%p)%p+p)%p;
    ll b=fib(k-2);
    ll gcds=gcd(b,p);
    if(m%gcds){puts("0");return;}
    else b/=gcds,m/=gcds,p/=gcds;

    m=m*inv(b,p)%p;
    printf("%lld\n",calc(r,m)-calc(l-1,m));
}

int main()
{
    int t;
    scanf("%d",&t);
    while(t--)
        mian();

    return 0;
}