【BZOJ1876】[SDOI2009]SuperGCD（数论，高精度）

【BZOJ1876】[SDOI2009]SuperGCD（数论，高精度）

题面

BZOJ
洛谷

题解

那些说数论只会\(gcd\)的人呢？我现在连\(gcd\)都不会，谁来教教我啊？
显然\(gcd\)除了辗转相除之外还可以辗转相减，然而辗转相减对于这题而言显然还不够优秀。
我们这样子来做。
如果当前\(a,b\)都是\(2\)的倍数，那么我们就把\(2\)直接同时除掉，直接在\(gcd\)中乘上一个\(2\)。否则如果只有一个数是\(2\)的倍数，显然可以直接把这个\(2\)给除掉。
这样子可以大大减少复杂度，这个似乎叫做\(Stein\)算法。
给个小提醒，判断一个数是不是\(2\)的倍数的时候，用\(\&1\)判断比用\(\%2\)判断快了\(20\)倍。

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
using namespace std;
#define ll long long
inline int read()
{
	int x=0;bool t=false;char ch=getchar();
	while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
	if(ch=='-')t=true,ch=getchar();
	while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
	return t?-x:x;
}
const ll yw=1000000000000000000;
char ch[11000];
struct BigNum
{
	ll s[800];int ws;
	void output()
		{
			printf("%lld",s[ws]);
			for(int i=ws-1;i;--i)
				printf("%018lld",s[i]);
			puts("");
		}
	void clear(){memset(s,0,sizeof(s));ws=0;}
	void init(char *ch)
		{
			int l=strlen(ch+1);reverse(&ch[1],&ch[l+1]);
			for(int i=1;i<=l;i+=18)
			{
				++ws;ll ss=1;
				for(int j=0;j<18;++j)
					if(i+j<=l)s[ws]+=(ch[i+j]-48)*ss,ss*=10;
					else break;
			}
		}
	void Div2()
		{
			for(int i=ws;i;--i)s[i-1]+=(s[i]&1)*yw,s[i]>>=1;s[0]=0;
			while(!s[ws])--ws;
		}
	void Multi2()
		{
			for(int i=1;i<=ws;++i)s[i]=s[i]<<1;
			for(int i=1;i<=ws;++i)s[i+1]+=s[i]/yw,s[i]%=yw;
			while(s[ws+1])++ws,s[ws+1]+=s[ws]/yw,s[ws]%=yw;
		}
}A,B,One,tmp;
BigNum operator-(BigNum a,BigNum b)
{
	int ws=max(a.ws,1);
	for(int i=1;i<=ws;++i)a.s[i]-=b.s[i];
	for(int i=ws-1;i;--i)if(a.s[i]<0)a.s[i]+=yw,a.s[i+1]-=1;
	while(!a.s[ws])--ws;
	a.ws=ws;return a;
}
bool operator<(BigNum a,BigNum b)
{
	if(a.ws!=b.ws)return a.ws<b.ws;
	for(int i=a.ws;i;--i)
		if(a.s[i]!=b.s[i])return a.s[i]<b.s[i];
	return false;
}int main()
{
	scanf("%s",ch+1);A.init(ch);
	scanf("%s",ch+1);B.init(ch);
	One.s[1]=One.ws=1;int two=0;
	if(B<A)swap(A,B);
	while(233)
	{
		if(A<One)break;
		if(!(A.s[1]&1)&&!(B.s[1]&1))A.Div2(),B.Div2(),++two;
		else if(!(A.s[1]&1))A.Div2();
		else if(!(B.s[1]&1))B.Div2();
		else B=B-A;if(B<A)swap(A,B);
	}
	while(two)B.Multi2(),--two;B.output();
	return 0;
}