hnu13150（hdu 4565）SO EASY!（矩阵快速幂）

本文主要是介绍hnu13150（hdu 4565）SO EASY!（矩阵快速幂），希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！

So Easy!

Time Limit: 1000ms, Special Time Limit:2500ms, Memory Limit:65536KB

Total submit users: 18, Accepted users: 8

Problem 13150 : No special judgement

Problem description

A sequence S_n is defined as:

Where a, b, n, m are positive integers.┌x┐is the ceil of x. For example, ┌3.14┐=4. You are to calculate S_n. You, a top coder, say: So easy!

Input

There are several test cases, each test case in one line contains four positive integers: a, b, n, m. Where 0< a, m < 2¹⁵, (a-1)²< b < a², 0 < b, n < 2³¹.
The input will finish with the end of file.

Output

For each the case, output an integer S_n.

Sample Input

2 3 1 2013
2 3 2 2013
2 2 1 2013

Sample Output

4
14
4

Problem Source

HNU Contest

原来博士出的SO EASY就是这个题。。。之前听过了但是没做过

之前做过一个类似的所以一看到就知道是用矩阵快速幂做

但是这个多了一个取上界第一直觉就是按之前的想法求值之后加一

按照这个想法写了之后跑了样例一样。。。就交了发现WA了

就意识到应该是那个上界的处理上出了问题

因为 ceil((a+sqrt(b))^n)=(a+sqrt(b))^n+(a-sqrt(b)^n)

令 a+sqrt(b)=x a-sqrt(b)=y

则(a+sqrt(b))^n+(a-sqrt(b)^n)=x^n+y^n=(x+y)(x^(n-1)+y^(n-1))-x*y(x^(n-2)+y^(n-2))

设 g(n)=x^n+y^n x+y=2*a=p x*y=a*a-b=q

|g(n) g(n-1)|=|g(n-1) g(n-2)|* |p 1|

|-q 0|

然后就可以由g(1) g(0) 2阶矩阵的n-1次幂求值

g(1)=2*a=p g(0)=2

注意数据范围所有的都要是__int64

代码如下：

#include <cstdio>
#include <iostream>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <string.h>
#include <string>#define MEM(a,x) memset(a,x,sizeof a)
#define eps 1e-8
//#define MOD 10009
#define MAXN 10010
#define INF 0x7fffffff
#define ll __int64
#define bug cout<<"here"<<endl;
#define fread freopen("ceshi.txt","r",stdin)
#define fwrite freopen("out.txt","w",stdout)using namespace std;ll mod,c,d;struct matrix
{ll a[2][2];
};matrix ori,res;void init()
{res.a[0][1]=res.a[1][0]=0;res.a[0][0]=res.a[1][1]=1;
}matrix multiply(matrix x,matrix y)
{matrix temp;memset(temp.a,0,sizeof(temp.a));for(int i=0;i<2;i++)for(int j=0;j<2;j++)for(int k=0;k<2;k++){temp.a[i][j]+=x.a[i][k]*y.a[k][j];temp.a[i][j]%=mod;}return temp;
}void calc(ll n)
{while(n){if(n&1)   res=multiply(ori,res);n>>=1;ori=multiply(ori,ori);}
}int main()
{
//    fread;ll n;while(scanf("%I64d%I64d%I64d%I64d",&c,&d,&n,&mod)!=EOF){init();ll x=2*c;ll y=c*c-d;ori.a[0][0]=x;ori.a[1][1]=0;ori.a[0][1]=1;ori.a[1][0]=-y;
//        ori.a[0][0]=ori.a[1][1]=c;
//        ori.a[0][1]=d;
//        ori.a[1][0]=1;calc(n-1);ll ans;ans=((x*res.a[0][0]+2*res.a[1][0])%mod+mod)%mod;//        double ans=0.0;
//        ll an=res.a[0][0]*c%mod+res.a[0][1]%mod;
//        an%=mod;
//        ll bn=res.a[1][0]*c%mod+res.a[1][1]%mod;
//        bn%=mod;
//        ans+=(double)1.0*res.a[0][0]*c%mod;
//        ans+=(double)1.0*res.a[0][1]%mod;
//        ans+=(double)(res.a[1][0]*c*1.0+1.0*res.a[1][1])*sqrt((double)5)%mod;
//        ans=ceil((double)((double)an+(double)bn*1.0*sqrt((double)d)));
//        printf("%llf\n",ans);
//        ans=(double)an+(double)bn*1.0*sqrt((double)d);
//        ll pp=(int)ans;
//        pp++;
//        pp%=mod;printf("%I64d\n",ans);}return 0;
}

这篇关于hnu13150（hdu 4565）SO EASY!（矩阵快速幂）的文章就介绍到这儿，希望我们推荐的文章对编程师们有所帮助！