hdu 3483 矩阵快速幂+二项式定理

Description

This is a very simple problem. Given three integers N, x, and M, your task is to calculate out the following value: 



Input

There are several test cases. For each case, there is a line with three integers N, x, and M, where 1 ≤ N, M ≤ 2*10 9, and 1 ≤ x ≤ 50. 

The input ends up with three negative numbers, which should not be processed as a case. 

Output

For each test case, print a line with an integer indicating the result.

Sample Input

100 1 10000
3 4 1000
-1 -1 -1

Sample Output

5050
444

例：

 

1. Sn=1^x * x^1 + 2^x * x^2 +...+ n^x * x^n; 
2. Sn+1=1^x * x^1 + 2^x * x^2 +...+ n^x * x^n+(n+1)^x * x^(n+1)=Sn+(n+1)^x * x^(n+1),将(n+1)^x二项式展开然后用矩阵快速幂，x^(n+1)则每次用乘x递推
 

构造矩阵: 

|1 xC(x,0) xC(x,1) xC(x,2) ... xC(x,x)|  |Sn       | |S(n+1)           | 

|0 xC(0,0) 0       0       ... 0      |  |x^n * n^0| |x^(n+1) * (n+1)^0| 

|0 xC(1,0) xC(1,1) 0       ... 0      | *|x^n * n^1|=|x^(n+1) * (n+1)^1| 

|0 xC(2,0) xC(2,1) xC(2,2) ... 0      |  |x^n * n^2| |x^(n+1) * (n+1)^2| 

|...                                  |  |...      | |...              | 

|0 xC(x,0) xC(x,1) xC(x,2) ... xC(x,x)|  |x^n * n^x| |x^(n+1) * (n+1)^x| 

#include<bits/stdc++.h>
#define inf 0x3f3f3f3f
#define ll long long
using namespace std;
struct node
{
ll mat[60][60];
};
ll n,m;
int x;
ll c[60][60];
void cc()
{
c[0][0]=c[1][0]=c[1][1]=1;
for(int i=2;i<=55;i++)
{
c[i][0]=c[i][i]=1;
for(int j=1;j<i;j++)
c[i][j]=c[i-1][j]+c[i-1][j-1];
}
}
node mult(node a,node b)
{
node c;
memset(c.mat,0,sizeof(c.mat));
for(int i=0;i<=x+1;i++)
for(int j=0;j<=x+1;j++)
for(int k=0;k<=x+1;k++)
c.mat[i][j]=(c.mat[i][j]+a.mat[i][k]*b.mat[k][j])%m;
return c;
}
int main()
{
cc();
while(cin>>n>>x>>m)
{
if(n==-1)
return 0;
node m1;
memset(m1.mat,0,sizeof(m1.mat));
m1.mat[0][0]=1;
for(int i=1;i<=x+1;i++)
m1.mat[0][i]=(x*c[x][i-1])%m;
for(int i=1;i<=x+1;i++)
{
for(int j=1;j<=i;j++)
m1.mat[i][j]=(x*c[i-1][j-1])%m;
}
node s;
memset(s.mat,0,sizeof(s.mat));
for(int i=0;i<=x+1;i++)
s.mat[i][i]=1;
int k=n-1;
while(k)
{
if(k%2==1)
s=mult(s,m1);
k=k/2;
m1=mult(m1,m1);
}
ll sum=0;
for(int i=0;i<=x+1;i++)
sum=(sum+x*s.mat[0][i])%m;
cout<<sum<<endl;
}
}