HDU 6211 Pythagoras（本原勾股数）

Description

给出一长度为2n2n的序列a0,...,a2n−1a0,...,a2n−1，求∑ay mod 2n∑ay mod 2n，其中yy满足存在x<y<z≤109,(x,y,z)=1x<y<z≤109,(x,y,z)=1使得x2+y2=z2x2+y2=z2

Input

第一行输入一整数TT表示用例组数，每组用例首先输入一整数nn，之后输入a0,...,a2n−1a0,...,a2n−1

(1≤T≤3,1≤n≤17,1≤ai≤255)(1≤T≤3,1≤n≤17,1≤ai≤255)

Output

求∑ay mod 2n∑ay mod 2n

Sample Input

3

2

0 0 0 1

2

1 0 0 0

2

1 1 1 1

Sample Output

39788763

79577506

159154994

Solution

只要找出所有不超过109109的本原勾股数即可，对于一组本原勾股数x2+y2=z2x2+y2=z2，先证明几个结论

1.x​x​和y​y​互素

若d=(x,y)>1d=(x,y)>1，由z2=x2+y2z2=x2+y2知d2|z2d2|z2，即d|zd|z，进而d|(x,y,z)=1d|(x,y,z)=1，矛盾，故(x,y)=1(x,y)=1

2.x,yx,y一奇一偶

由结论11知，如果x,yx,y不是一奇一偶则必然为两个奇数，进而zz为偶数，假设x=2a+1,y=2b+1,z=2cx=2a+1,y=2b+1,z=2c，由x2+y2=z2x2+y2=z2知4a2+4a+1+4b2+4b+1=4c24a2+4a+1+4b2+4b+1=4c2，即a2+b2−c2+a+b=12a2+b2−c2+a+b=12，矛盾，故x,yx,y一奇一偶

3.存在奇偶性不同的两个数s,ts,t满足s>t,(s,t)=1s>t,(s,t)=1且x=s2−t2,y=2st,z=s2+t2x=s2−t2,y=2st,z=s2+t2

假设xx为奇数，yy为偶数，则zz必然为奇数，x2=z2−y2=(z−y)(z+y)x2=z2−y2=(z−y)(z+y)，首先证明d=(z+y,z−y)=1d=(z+y,z−y)=1

由于z+y,z−yz+y,z−y均为奇数，故dd为奇数，由于d|(z+y+z−y)=2z,d|((z+y)−(z−y))=2yd|(z+y+z−y)=2z,d|((z+y)−(z−y))=2y，故d|y,d|zd|y,d|z，且注意到d2|(z+y)(z−y)=x2d2|(z+y)(z−y)=x2，故d|xd|x，进而d|(x,y,z)=1d|(x,y,z)=1

由于z+yz+y和z−yz−y互素，故z+y,z−yz+y,z−y为两个互素的完全平方数，设z+y=m2,z−y=n2z+y=m2,z−y=n2，则x=mn,y=m2−n22,z=m2+n22x=mn,y=m2−n22,z=m2+n22，其中m,nm,n均为奇数且m>n,(m,n)=1m>n,(m,n)=1

令s=m+n2,t=m−n2s=m+n2,t=m−n2，则s>ts>t均为正整数，且由于s+t=m,s−t=ns+t=m,s−t=n，故(s+t,s−t)=1(s+t,s−t)=1，进而(s,t)=1(s,t)=1，而s+t=ms+t=m为奇数，显然s,ts,t一奇一偶，而代入即得x=s2−t2,y=2st,z=s2+t2x=s2−t2,y=2st,z=s2+t2

由结论33，只要枚举所有的(s,t)(s,t)即可得到所有的本原勾股数，而由于z=s2+t2≤109z=s2+t2≤109，故只需要枚举小范围的ss，然后找小于ss、与ss奇偶性不同且与ss互素的tt即可

Code

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<vector>
#include<queue>
#include<map>
#include<set>
#include<ctime>
using namespace std;
typedef long long ll;
const int mod=(1<<17)-1;
int T,n,a[mod+5],f[mod+5];
void init(int N=1e9)
{
int p[33];
for(int i=1;i*i<=N;i++)
{
int res=0,ii=i;
for(int j=2;j*j<=ii;j++)
if(ii%j==0)
{
p[res++]=j;
while(ii%j==0)ii/=j;
}
if(ii>1)p[res++]=ii;
int c=i*i;
for(int j=1;j<i;j++)
{
c+=2*j-1;
if(c>N)break;
if(!((i-j)&1))continue;
int k=0;
while(k<res&&j%p[k]!=0)k++;
if(k==res)f[max(i*i-j*j,2*i*j)&mod]++;
}
}
}
int main()
{
init();
scanf("%d",&T);
while(T--)
{
scanf("%d",&n);
n=(1<<n)-1;
for(int i=0;i<=n;i++)scanf("%d",&a[i]);
ll ans=0;
for(int i=0;i<=mod;i++)ans+=(ll)a[i&n]*f[i];
printf("%I64d\n",ans);
}
return 0;
}