Tyvj4876：近似排列计数 （矩阵快速幂）

题目传送门：http://tyvj.cn/p/4876

---

题目分析：比赛的时候见到这题作为T3出现，想了5min就知道是个矩阵乘法，然而由于T1T2花了太久时间，只剩下半小时，就没有敲代码。比赛之后又想清楚了一些细节，过来把这题补了。
由于k最大只有2，当k=2时，符合条件的排列第i位一定是i+2,i+1,i,i-1,i-2中的一个，我们不妨用一个状态压缩将这几位有没有选记下来；又考虑到n高达$10^9$，不妨将转移矩阵记下来，对它做矩阵快速幂；要强制某一个位置取某一个值的话，就将n分成m段，每一段做矩阵快速幂，再暴力合并段与段之间即可。
经过仔细研究，我发现k=2时，只需要记i+2,i+1,i,i-1有没有选过即可。我们记f[i][s]表示选完前i个位置的数，i+2,i+1,i,i-1的状态为s的方案数。接下来要考虑第i+1位选什么数，如果i-1还没选过，第i+1位就一定要选i-1，否则它可以选i~i+3中没选过的任意一个。若第x位一定要选y，就先将f[x-1]构出来，暴力判断每一个f[x-1][s]是否合法，合法才能成功转移到f[x]的对应位置。时间复杂度$O(Tm\log(n)2^{6k})$。（其实不难想，就是写起来很烦）
话说这题我刚码完代码就过样例，交上去立马AC了，调都没调过，我是不是应该去买彩票？

---

CODE：

#include<iostream>
#include<string>
#include<cstring>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<stdio.h>
#include<algorithm>
using namespace std;

const int maxs=16;
const int maxm=110;
const long long M=1000000007;
typedef long long LL;

struct mat
{
    LL val[maxs][maxs];
} ;
mat e,a[2];

struct data
{
    int x,y;
} b[maxm];

LL st[maxs];
LL ed[maxs];
int t,n,m,k,ms;

bool Comp(data p,data q)
{
    return p.x<q.x;
}

int Abs(int p)
{
    if (p>=0) return p;
    return -p;
}

mat Times(mat p,mat q)
{
    mat r;
    for (int i=0; i<ms; i++)
        for (int j=0; j<ms; j++)
            r.val[i][j]=0;
    for (int i=0; i<ms; i++)
        for (int j=0; j<ms; j++)
            for (int k=0; k<ms; k++)
                r.val[i][j]=(r.val[i][j]+p.val[i][k]*q.val[k][j]%M)%M;
    return r;
}

mat Fast_power(int p)
{
    if (!p) return e;
    mat temp=Fast_power(p>>1);
    temp=Times(temp,temp);
    if (p&1) temp=Times(temp,a[k]);
    return temp;
}

int main()
{
    freopen("count.in","r",stdin);
    freopen("count.out","w",stdout);

    for (int i=0; i<maxs; i++) e.val[i][i]=1;
    for (int i=0; i<4; i++)
        if (i&1) a[0].val[i][(i>>1)|2]=1;
        else
        {
            int j=(i>>1)|2;
            if (j&1) a[0].val[i][j^1]=1;
            if (j&2) a[0].val[i][j^2]=1;
        }
    for (int i=0; i<maxs; i++)
        if (i&1) a[1].val[i][(i>>1)|8]=1;
        else
        {
            int j=(i>>1)|8;
            if (j&1) a[1].val[i][j^1]=1;
            if (j&2) a[1].val[i][j^2]=1;
            if (j&4) a[1].val[i][j^4]=1;
            if (j&8) a[1].val[i][j^8]=1;
        }

    scanf("%d",&t);
    while (t--)
    {
        scanf("%d%d%d",&n,&m,&k);
        for (int i=1; i<=m; i++) scanf("%d%d",&b[i].x,&b[i].y);
        sort(b+1,b+m+1,Comp);

        bool sol=true;
        for (int i=1; i<=m; i++)
            if ( b[i].x<1 || b[i].x>n || b[i].y<1 || b[i].y>n || Abs(b[i].x-b[i].y)>k )
            {
                printf("0\n");
                sol=false;
                break;
            }
        if (!sol) continue;

        for (int i=2; i<=m; i++)
            if (b[i-1].x==b[i].x)
            {
                printf("0\n");
                sol=false;
                break;
            }
        if (!sol) continue;

        if (!k)
        {
            printf("1\n");
            continue;
        }
        k--;

        int h=1;
        b[0].x=0;
        memset(st,0,sizeof(st));
        if (b[1].x==1)
        {
            h=2;
            if (k) st[14^(1<<b[1].y)]=1;
            else st[3^(1<<(b[1].y-1))]=1;
        }
        else
            if (k) st[12]=1;
            else st[2]=1;
        if (k) ms=maxs;
        else ms=4;

        for (int u=h; u<=m; u++)
        {
            int len=b[u].x-b[u-1].x-1;
            mat temp=Fast_power(len);
            memset(ed,0,sizeof(ed));
            for (int i=0; i<ms; i++)
                for (int j=0; j<ms; j++)
                    ed[j]=(ed[j]+st[i]*temp.val[i][j]%M)%M;
            for (int i=0; i<ms; i++) st[i]=ed[i];

            memset(ed,0,sizeof(ed));
            if (k)
            {
                if (b[u].y==b[u].x-k-1)
                {
                    for (int i=0; i<ms; i++)
                        if (i&1) ed[(i>>1)|8]=(ed[(i>>1)|8]+st[i])%M;
                }
                else
                {
                    int j=1<<(b[u].y-b[u].x+1);
                    for (int i=0; i<ms; i++)
                        if ( (!(i&1)) && (((i>>1)|8)&j) )
                            ed[((i>>1)|8)^j]=(ed[((i>>1)|8)^j]+st[i])%M;
                }
            }
            else
            {
                if (b[u].y==b[u].x-k-1)
                {
                    for (int i=0; i<ms; i++)
                        if (i&1) ed[(i>>1)|2]=(ed[(i>>1)|2]+st[i])%M;
                }
                else
                {
                    int j=1<<(b[u].y-b[u].x);
                    for (int i=0; i<ms; i++)
                        if ( (!(i&1)) && (((i>>1)|2)&j) )
                            ed[((i>>1)|2)^j]=(ed[((i>>1)|2)^j]+st[i])%M;
                }
            }
            for (int i=0; i<ms; i++) st[i]=ed[i];
        }

        mat temp=Fast_power(n-b[m].x);
        memset(ed,0,sizeof(ed));
        for (int i=0; i<ms; i++)
            for (int j=0; j<ms; j++)
                ed[j]=(ed[j]+st[i]*temp.val[i][j]%M)%M;
        for (int i=0; i<ms; i++) st[i]=ed[i];

        if (k) printf("%lld\n",st[12]);
        else printf("%lld\n",st[2]);
    }

    return 0;
}

---

一些后话：
做题的时候tututu跟我讲了另一道题：
有两个长度为n的数组a,b，定义其乘法运算为：

$$\sum_{j=0}^{n-1}\sum_{k=0}^{n-1}c[(j+k)\mod{n}]+=a[j]*b[k]$$
c为乘出来的结果。现在给定数组x和数字y，求 $x^y$。
这很明显可以用倍增来做，然而在转移的时候它只需要 $n^2$的时间，所以时间复杂度是 $O(n^2\log(p))$的。然而如果写成矩阵乘法的形式，时间复杂度就是 $O(n^3\log(p))$的，这样显然没有必要。
tututu：那什么时候快速幂可以做到 $n^2$转移，什么时候要 $n^3$转移呢？Tyvj4876这题又可不可以做到 $2^{4k}$转移呢？
关于这个问题我们讨论了很久，最后得出的答案是：如果转移的信息可以用一个数组+一些运算法则存下来，就没有必要存矩阵，转移的时间也就可以做到 $n^2$。而Tyvj4876这题并不能找到一个完全存储转移信息的数组和运算，也就不能做到 $2^{4k}$转移。
（如果路过的神犇有什么更好的想法欢迎指点）