HDU 5674 Function

Problem Description

There is a function  f(n)=(a+b√)n+(a−b√)n . a and b are integers  (1≤a,b≤1,000,000) .
Maybe the function looks complex but it is actually an integer.The question is to calculate  f(xy) .The answer can be
very large,so just output the answer mod  1,000,000,007 .

 

Input

There are multiple test cases. The first line of input contains an integer  T(1≤T≤200)  indicating the number of test cases. For each test case:

One line contains four integers  a,b    (1≤a,b≤1,000,000) ,   x(1≤x≤50),y(1≤y≤1018) .

 

Output

For each test case, output one integer.

 

Sample Input

  
  

   
   3
3 5 1 1
3 5 2 1
1 1 49 99999
  
  

 

Sample Output

  
  

   
   6
28

   
   
160106184
   
   
循环节的证明真的是厉害啊，这种东西不知道的话实在是没什么办法。。
   
   

   
   
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<iostream>
#include<algorithm>
using namespace std;
typedef long long LL;
const int low(int x) { return x&-x; }
const int INF = 0x7FFFFFFF;
const int mod = 1e9 + 7;
const int maxn = 1e6 + 10;
int T;
LL a, b, x, y, z;

LL mul(LL a, LL b, LL mod) 
{
	LL ret;
	for (ret = 0; b; b >>= 1, (a <<= 1) %= mod)
	{
		if (b & 1) ret = (ret + a) % mod;
	}
	return ret;
}

LL get(LL x, LL y, LL z)
{
	LL ans = 1;
	for (; y; y >>= 1)
	{
		if (y & 1) ans = mul(ans, x, z);
		x = mul(x, x, z);
	}
	return ans;
}

struct martix
{
	LL a[2][2];
}A, B, e;

martix operator*(const martix&a, const martix&b)
{
	martix c;
	for (int i = 0; i < 2; i++)
	{
		for (int j = 0; j < 2; j++)
		{
			c.a[i][j] = 0;
			for (int k = 0; k < 2; k++)
			{
				(c.a[i][j] += a.a[i][k] * b.a[k][j] % mod) %= mod;
			}
		}
	}
	return c;
}

int main()
{
	scanf("%d", &T);
	while (T--)
	{
		cin >> a >> b >> x >> y;
		z = get(x, y, (LL)mod*mod - 1);

		A.a[0][0] = 2;	A.a[0][1] = 2 * a % mod;
		B.a[0][0] = 0;	B.a[0][1] = ((b - a*a) % mod + mod) % mod;
		B.a[1][0] = 1;	B.a[1][1] = 2 * a % mod;
		for (int i = 0; i < 3; i++)
		{
			for (int j = 0; j < 3; j++)
			{
				e.a[i][j] = i == j;
			}
		}
		for (--z; z; z >>= 1)
		{
			if (z & 1) e = e * B;
			B = B * B;
		}
		A = A*e;
		cout << A.a[0][1] << endl;
	}
	return 0;
}