HDU 6442 GuGu Convolution（快速幂）

Description

给出两个整数$A,B$，定义序列$f(x)=\sum _{i=0}^{\infty }\frac{A^i}{i!}{x}^i$，$g(x)=\sum _{i=0}^{\infty }\frac{(\sqrt{B}{)}^{2i+1}}{(2i+1)!}{x}^{2i+1}$，求$n!\cdot [x^n]f\ast g$

Input

首先输入一整数$T$表示用例组数，每组用例输入四个整数$A,B,n,p$

$(1\le T\le 10^5,1\le A,B\le 10^6,1\le n\le 10^{18},1\le p\le 10^9)$

Output

求$n!\cdot [x^n]f\ast g$，把结果写成$\sum _{i=1}^q{a}_i\sqrt{b_i}$的形式，其中$b_i$为无平方因子数且$b_i̸=b_j$，输出$q$，$a_i\%p$和$b_i$

Sample Input

3
1 1 1 7
523 12 2 2100
1 1 1000000000000000000 998244353

Sample Output

1 1 1
1 2092 3
1 121099884 1

Solution

答案即为$\sum _{i=0,i\&1=1}^n{C}_n^i{A}^{n-i}(\sqrt{B}{)}^i$，也即$\frac{(A+\sqrt{B}{)}^n-(A-\sqrt{B}{)}^n}{2}$，分子这两部分都可以写成$X+Y\sqrt{B}$的形式，且两者的$X$相同，$Y$相反，答案即为$Y\sqrt{B}$，故对于该二元组$(X,Y)$定义加法和乘法之后用快速幂即可得到$Y$值，注意要先把$B$拆分为一个无平方因子数和另一个完全平方数的乘积

Code

#include<cstdio>
#include<algorithm>
#include<vector>
using namespace std;
typedef long long ll;
typedef pair<int,int>P;
int mod,w;
int add(int x,int y)
{
	x+=y;
	if(x>=mod)x-=mod;
	return x;
}
int mul(int x,int y)
{
	ll z=1ll*x*y;
	return z-z/mod*mod;
}
#define x1 first
#define x2 second
P Add(P A,P B)
{
	return P(add(A.x1,B.x1),add(A.x2,B.x2));
}
P Mul(P A,P B)
{
	return P(add(mul(A.x1,B.x1),mul(mul(A.x2,B.x2),w)),add(mul(A.x1,B.x2),mul(A.x2,B.x1)));
}
int Pow(P A,ll n)
{
	P ans=P(1,0);
	while(n)
	{
		if(n&1)ans=Mul(ans,A);
		A=Mul(A,A);
		n>>=1;
	}
	return ans.x2;
}
int main()
{
	int T,A,B,C;
	ll n;
	scanf("%d",&T);
	while(T--)
	{
		scanf("%d%d%lld%d",&A,&B,&n,&mod);
		w=1,C=1;
		for(int i=2;i*i<=B;i++)
			if(B%i==0)
			{
				int num=0;
				while(B%i==0)B/=i,num++;
				for(int j=0;j<num/2;j++)C*=i;
				if(num&1)w*=i;
			}
		if(B>1)w*=B;
		printf("1 %d %d\n",Pow(P(A,C),n),w);
	}
	return 0;
}