CSU1516矩阵快速幂+费马小定理

　　通过这个题真的学到了不少东西，最起码矩阵快速幂算是入门了，普通快速幂也彻底明白了(以前都是打模板)，了解了费马小定理

　　关键点 求(a^fib[b])%p 的值其中p是素数，0<a<p, b在int范围内

　　先假设fib[b]>p-1 那么上式

　　(a^fib[b])%p = (a^(p-1)*a^(p-1)*....*a^(p-1)*a^m)%p

　　(这里 m = fib[b]%(p-1))

　　由于p是素数且a<p那么gcd(a,p)=1，所以由费小可得(a^(p-1))%p=1

　　那上式就 = (a^(fib[b]%(p-1)))%p

　　然后fib[b]%(p-1)可以用矩阵快速幂求得

　　最后在用普通的快速幂就可以求得(a^m)%p

　　得到了这些剩下的就只是一个区间dp了

　　一时着急写出来，代码写的很丑，LL和int到处乱定义。。。

　　下面是运行结果，时限2s，有点险。还有1s过的。。。

 1 #include <cstdio>
 2 #include <cstring>
 3 #include <iostream>
 4 #include <algorithm>
 5 using namespace std;
 6 const int maxn = 105;
 7 const int INF = 0x3f3f3f3f;
 8 typedef long long LL;
 9 int a[maxn],b[maxn];
10 int n;
11 LL p;
12 LL sum[maxn],dp[maxn][maxn];
13 struct node
14 {
15     LL a[2][2];
16     node ()
17     {
18         a[0][0] = a[1][0] = a[0][1] = 1;
19         a[1][1] = 0;
20     }
21 };
22 node mul(node A,node B)
23 {
24     node ans;
25     for(int i = 0;i<2;++i)for(int j = 0;j<2;++j)
26     {
27         ans.a[i][j] = 0;
28         for(int k = 0;k<2;++k)ans.a[i][j] = (ans.a[i][j]+A.a[i][k]*B.a[k][j])%(p-1);
29     }
30     return ans;
31 }
32 LL mypow(LL x,LL m,LL mod)
33 {
34     LL ret = 1;
35     while(m)
36     {
37         if(m&1)ret = (ret*x)%mod;
38         x = (x*x)%mod;
39         m/=2;
40     }
41     return ret;
42 }
43 int fib(int x)
44 {
45     node ans,s;
46     while(x)
47     {
48         if(x&1)ans = mul(ans,s);
49         s = mul(s,s);
50         x/=2;
51     }
52     return ans.a[1][0]%(p-1);
53 }
54 int solve(int x)
55 {
56     int m = fib(b[x]-1);
57     return mypow(a[x],m,p);
58 }
59 LL gcd(LL x,LL y)
60 {
61     return y?gcd(y,x%y):x;
62 }
63 int main()
64 {
65 //    freopen("in.txt","r",stdin);
66     int kase = 1;
67     while(cin>>n>>p)
68     {
69         for(int i = 1;i<=n;++i)for(int j = 1;j<=n;++j)dp[i][j] = i==j?0:INF;
70         for(int i = 1;i<=n;++i)scanf("%d",a+i);
71         for(int i = 1;i<=n;++i)scanf("%d",b+i);
72         for(int i = 1;i<=n;++i)a[i] = solve(i);
73
74         sum[0] = 0;sum[1] = a[1];
75         for(int i = 2;i<=n;++i)sum[i] = sum[i-1]+(LL)a[i];
76         printf("case %d: ",kase++);
77         if(n==1){printf("0\n");continue;}
78         for(int r = 2;r<=n;++r)for(int i = 1;i<=n-r+1;++i)
79         {
80             int j = i+r-1;
81             for(int k = i;k<j;++k)
82                 dp[i][j] = min(dp[i][j],dp[i][k]+dp[k+1][j]+gcd(sum[k]-sum[i-1],sum[j]-sum[k-1]));
83         }
84         cout<<dp[1][n]<<endl;
85     }
86     return 0;
87 }