Description
Input
第一行有四个整数 n, p, k, r，所有整数含义见问题描述。
1 ≤ n ≤ 10^9, 0 ≤ r < k ≤ 50, 2 ≤ p ≤ 2^30 ? 1
Output
一行一个整数代表答案。
Sample Input
2 10007 2 0
Sample Output
8
HINT
Source
黑吉辽沪冀晋六省联考
DP+矩阵乘法+思路~
可以发现要求的式子相当于是从n个物品中取出%k==r个物品的方案数。所以我们用f[i][j]表示从i个物品中取出%k==j个物品的方案数，那么f[i][j]=f[i-1][j]+f[i-1][j-1]（j==0时为f[i-1][k-1]）。
那么我们用矩阵乘法优化一下这个DP就好了。
注意矩阵乘法的顺序不能换；n要*k所以要开成long long。
#include<cstdio>
#include<cstring>
#include<iostream>
using namespace std;
#define ll long longll n,mod,k,r;struct node{ll a[50][50];node operator * (const node&u){node z;for(int i=0;i<k;i++)for(int j=0;j<k;j++){z.a[i][j]=0;for(int b=0;b<k;b++) z.a[i][j]=(z.a[i][j]+a[i][b]*u.a[b][j]%mod)%mod;}return z;}
}a,c,e;int read()
{int x=0,f=1;char ch=getchar();while(ch<'0' || ch>'9') {if(ch=='-') f=-1;ch=getchar();}while(ch>='0' && ch<='9') {x=x*10+ch-'0';ch=getchar();}return x*f;
}int main()
{n=read();mod=read();k=read();r=read();e.a[0][0]=1;for(int i=0;i<k;i++) a.a[(i-1+k)%k][i]++,a.a[i][i]++,c.a[i][i]=1;n*=k;for(;n;n>>=1,a=a*a) if(n&1) c=c*a;c=e*c;printf("%lld\n",c.a[0][r]);return 0;
}