首页 > 代码库 > [hdu 4959]Poor Akagi 数论(卢卡斯数,二次域运算,等比数列求和)
[hdu 4959]Poor Akagi 数论(卢卡斯数,二次域运算,等比数列求和)
Poor Akagi
Time Limit: 30000/15000 MS (Java/Others) Memory Limit: 32768/32768 K (Java/Others)Total Submission(s): 131 Accepted Submission(s): 29
Problem Description
Akagi is not only good at basketball but also good at math. Recently, he got a sequence Ln from his teacher. Ln is defined as follow:
$$\Large L(n)=\begin{cases}
2 & \text{ if } n=0 \\
1 & \text{ if } n=1 \\
L(n-1)+L(n-2) & \text{ if } n>1
\end{cases}$$
And Akagi’s teacher cherishes Agaki’s talent in mathematic. So he wants Agaki to spend more time studying math rather than playing basketball. So he decided to ask Agaki to solve a problem about Ln and promised that as soon as he solves this problem, he can go to play basketball. And this problem is:
Given N and K, you need to find \(\Large\sum\limits_{0}^{N}L_i^K\)
And Agaki needs your help.
$$\Large L(n)=\begin{cases}
2 & \text{ if } n=0 \\
1 & \text{ if } n=1 \\
L(n-1)+L(n-2) & \text{ if } n>1
\end{cases}$$
And Akagi’s teacher cherishes Agaki’s talent in mathematic. So he wants Agaki to spend more time studying math rather than playing basketball. So he decided to ask Agaki to solve a problem about Ln and promised that as soon as he solves this problem, he can go to play basketball. And this problem is:
Given N and K, you need to find \(\Large\sum\limits_{0}^{N}L_i^K\)
And Agaki needs your help.
Input
This problem contains multiple tests.
In the first line there’s one number T (1 ≤ T ≤ 20) which tells the total number of test cases. For each test case, there an integer N (0 ≤ N ≤ 10^18) and an integer K (1 ≤ K ≤ 100000) in a line.
In the first line there’s one number T (1 ≤ T ≤ 20) which tells the total number of test cases. For each test case, there an integer N (0 ≤ N ≤ 10^18) and an integer K (1 ≤ K ≤ 100000) in a line.
Output
For each test case, you need to output the answer mod 1000000007 in a line.
Sample Input
3 3 1 2 2 4 3
Sample Output
10 14 443
Source
BestCoder Round #5
题目大意
求卢卡斯数的k次方的前n项和
卢卡斯数为L[0]=2,L[1]=1,L[n]=L[n-2]+L[n-1](n>=2)
题目思路
当时看到题还以为直接根据 zoj 3774 找出二次剩余…… 结果发现1e9+7不存在二次剩余
最后发现了一种很巧妙的做法
直接根据卢卡斯数的通项公式
设
则求和公式为定义二次域
此时直接对二次域进行加、乘操作即可(最后的结果为整数,故根号五不会存在在结果之中)
重载二次域的加号和乘号,定义二次域的快速幂运算,全部带入公式即可。
=.=好像这一题的杭电的数据还没有修正
公比为一时直接返回n+1(可能带来溢出)竟然AC了
然后正解依然WA……
这里只放正解代码
/** **author : ahm001 ** **source : hdu 4959** **time : 08/21/14 ** **type : math ** **/ #include <cstdio> #include <iostream> #include <algorithm> #include <ctime> #include <cctype> #include <cmath> #include <string> #include <cstring> #include <stack> #include <queue> #include <list> #include <vector> #include <map> #include <set> #define sqr(x) ((x)*(x)) #define LL long long #define INF 0x3f3f3f3f #define PI acos(-1.0) #define eps 1e-10 #define mod 1000000007 using namespace std; int cnt=0; typedef pair <LL,LL> qf; qf operator + (qf a,qf b) { return make_pair((a.first+b.first)%mod,(a.second+b.second)%mod); } qf operator * (qf a,qf b) { // if ((((LL)a.first*(LL)b.first)%mod+((LL)a.second*(LL)b.second)%mod*5ll)%mod<0) // printf("%d %d %d %d\n",a.first,a.second,b.first,b.second); if (a.first<0) a.first+=mod; if (b.first<0) b.first+=mod; if (a.second<0) a.second+=mod; if (b.second<0) b.second+=mod; return make_pair((((LL)a.first*(LL)b.first)%mod+((LL)a.second*(LL)b.second)%mod*5ll)%mod, (((LL)a.first*(LL)b.second)%mod+((LL)a.second*(LL)b.first)%mod)%mod); } qf pow(qf a, LL x) { qf res(1,0); qf multi=a; while (x) { if (x&1) { res=res*multi; } multi=multi*multi; x/=2; } return res; } LL pow(LL a,LL b) { LL res=1; LL multi=a; while (b) { if (b&1) { res=res*multi%mod; } multi=multi*multi%mod; b/=2; } return res; } qf acce(qf a,LL b) { qf ans=make_pair(1,0); // if (a==ans) return make_pair(b+1,0);//这条语句去掉后AC了,但是n+1不取模将会造成后面的结果爆掉 qf powe=a; qf sum=a; qf multi=make_pair(1,0); while (b) { if (b&1) { ans=ans+(multi*sum); multi=multi*powe; } sum=sum*(powe+make_pair(1,0)); powe=powe*powe; b/=2; } return ans; } LL inv[100005]; qf r1[100005],r2[100005]; void egcd (LL a,LL b,LL &x,LL &y) { if (b==0) { x=1,y=0; return ; } egcd(b,a%b,x,y); LL t=x; x=y;y=t-a/b*y; } int main() { LL x,y; for (LL i=1;i<=100000;i++) { egcd(i,mod,x,y); x=(x+mod)%mod; inv[i]=x; } r1[0]=make_pair(1,0); r2[0]=make_pair(1,0); for (int i=1;i<=100000;i++) { r1[i]=r1[i-1]*make_pair(1,1); r2[i]=r2[i-1]*make_pair(1,-1); } int T; scanf("%d",&T); while (T--) { cnt=0; LL n,m; scanf("%I64d%I64d",&n,&m); // n=1e18; // m=1e5; qf ans=make_pair(0,0); LL Ca=1; LL v=pow(inv[2],m); for (LL i=0;i<=m;i++) { // printf("%lld\n",Ca); qf p(Ca,0); qf tmp=r1[i]*r2[m-i]*make_pair(v,0); tmp=acce(tmp,n); tmp=tmp*p; ans=ans+tmp; Ca=Ca*(m-i)%mod; Ca=Ca*inv[i+1]%mod; } LL aa=(LL)ans.first; printf("%I64d\n",aa); // printf("%d %d \n",ans.first,ans.second); // printf("%d\n",cnt); } return 0; }
声明:以上内容来自用户投稿及互联网公开渠道收集整理发布,本网站不拥有所有权,未作人工编辑处理,也不承担相关法律责任,若内容有误或涉及侵权可进行投诉: 投诉/举报 工作人员会在5个工作日内联系你,一经查实,本站将立刻删除涉嫌侵权内容。