hdu4565So Easy!+矩阵快速幂+广义斐波拉契

Problem Description
　　A sequence Sn 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 Sn.
　　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 < 215, (a-1)2< b < a2, 0 < b, n < 231.The input will finish with the end of file.

Output
　　For each the case, output an integer Sn.

Sample Input

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

Sample Output

4
14
4

Source
2013 ACM-ICPC长沙赛区全国邀请赛——题目重现

直接递推构造矩阵，然后快速幂。参考上一篇博客：http://blog.csdn.net/xtulollipop/article/details/52382791
$(a-1)^2< b < a^2$所以依然可以$(a+\sqrt b)^n$与$(a-\sqrt b)^n$相加。
向上取整就是向下取整+1，所以有答案2*An-1+1=2*An

#include<cstdio>
#include<cmath>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cstdlib>
#include<queue>
#include<vector>
#include<map>
#include<stack>
#include<set>
using namespace std;
#define pi acos(-1.0)
#define EPS 1e-6    //log(x)
#define e exp(1.0); //2.718281828
//#define mod 1000000007
#define INF 0x7fffffff
#define inf 0x3f3f3f3f
#pragma comment(linker,"/STACK:102400000,102400000")
typedef long long LL;

#define debug(x) cout<<x<<endl;
#define debug2(x) cout<<x<<" ";

//#define MOD 10000007
LL MOD;
struct Mat{
    int n,m;
    LL mat[9][9];
};
Mat operator *(Mat a,Mat b){
    Mat c;
    memset(c.mat,0,sizeof(c.mat));
    c.n = a.n,c.m = b.m;

    for(int i=1;i<=a.n;i++){
        for(int j=1;j<=b.m;j++){
            for(int k=1;k<=a.m;k++){
                c.mat[i][j] += (a.mat[i][k]*b.mat[k][j])%MOD;
                c.mat[i][j] %= MOD;
            }
        }
    }
    return c;
}
Mat operator +(Mat a,Mat b){
    Mat c;
    memset(c.mat,0,sizeof(c.mat));
    c.n = a.n,c.m = a.m;

    for(int i=1;i<=a.n;i++){
        for(int j=1;j<=a.m;j++){
            c.mat[i][j] = (a.mat[i][j]+b.mat[i][j])%MOD;
        }
    }
    return c;
}
Mat operator ^(Mat a,LL k){
    Mat c;
    memset(c.mat,0,sizeof(c.mat));
    c.n = a.n,c.m = a.n;
    for(int i=1;i<=a.n;i++)c.mat[i][i] = 1;

    while(k){
        if(k&1){
            c = c*a;
        }
        a = a*a;
        k>>=1;
    }
    return c;
}
void out(Mat a){
    for(int i=1;i<=a.n;i++){
        for(int j=1;j<=a.m;j++){
            printf(j==a.m? "%I64d\n":"%I64d ",a.mat[i][j]);
        }
    }
}
LL quickPow(LL x, LL n, LL mm)
{
    LL a = 1;
    while (n)
    {
        a *= n&1 ? x : 1;
        a %= mm;
        n >>= 1 ;
        x *= x;
        x %= mm;
    }
    return a;
}
int main()
{
    LL a,b,n;
    while(scanf("%I64d %I64d %I64d %I64d",&a,&b,&n,&MOD)!=EOF){
        if(n==1){
            printf("%I64d\n",((LL)(a+ceil(sqrt(b*1.0)))%MOD));
            continue;
        }
        Mat pp;
        pp.n=pp.m=2;
        pp.mat[1][1]=a%MOD;
        pp.mat[1][2]=b%MOD;
        pp.mat[2][1]=1;
        pp.mat[2][2]=a%MOD;

        Mat A0;
        A0.n=2,A0.m=1;
        A0.mat[1][1]=a%MOD;
        A0.mat[2][1]=1;

        Mat ans=pp^(n-1);
        ans=ans*A0;
        printf("%I64d\n",(2*ans.mat[1][1])%MOD);
    }
    return 0;
}


/*
                   _ooOoo_
                  o8888888o
                  88" . "88
                  (| -_- |)
                  O\  =  /O
               ____/`---'\____
             .'  \\|     |//  `.
            /  \\|||  :  |||//  \
           /  _||||| -:- |||||-  \
           |   | \\\  -  /// |   |
           | \_|  ''\---/''  |   |
           \  .-\__  `-`  ___/-. /
         ___`. .'  /--.--\  `. . __
      ."" '<  `.___\_<|>_/___.'  >'"".
     | | :  `- \`.;`\ _ /`;.`/ - ` : | |
     \  \ `-.   \_ __\ /__ _/   .-` /  /
======`-.____`-.___\_____/___.-`____.-'======
                   `=---='
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
         I have a dream!A AC deram!!
 orz orz orz orz orz orz orz orz orz orz orz
 orz orz orz orz orz orz orz orz orz orz orz
 orz orz orz orz orz orz orz orz orz orz orz

*/