HDU 3579 模餘方程

/*HDU 3579不一定互質情況，神牛公式解析連結：http://yzmduncan.iteye.com/blog/1323599把推導過程翻譯成代碼即可... */
#include<iostream>#include<cstdio>#include<cmath>using namespace std;typedef long long ll; 
ll gcd(ll a,ll b){    if(b==0) return a;    return gcd(b,a%b);} 
ll kzgcd(ll a,ll b,ll&x,ll&y){    if(b==0) {        x=1, y=0;        return a;    }    ll d = kzgcd(b,a%b,x,y);    ll k = y;    y = x - a/b*y;    x = k;    return d;}
ll inv(ll a,ll b){ /// 求 a 的逆元 整數最小x    ll x,y;    ll k = kzgcd(a,b,x,y);     if(k!=1) return -1;    return (x%b+b)%b;}
bool merge(ll a1,ll n1,ll&a2,ll&n2){    ll d = gcd(n1,n2);    ll c = a2 - a1;    if(c%d) return false;    ll k = inv(n1/d,n2/d)*c/d;    ll a,n;    n = n1/d*n2;    a = (n1*k+a1)%(n1*n2/d);    a2 = a;    n2 = n;    return true;}
ll china_reminder(ll k,ll*a,ll*b){ /// 餘數，mod數     ll n1=b[1],a1=a[1];    for(ll i=2;i<=k;i++){        ll nn,aa;        nn = b[i], aa = a[i];        if(!merge(a1,n1,aa,nn))            return -1;        n1 = nn;        a1 = aa;    }     return (a1%n1+n1)%n1;}
int main(){     ll t,a[1110],b[1100];     int f;     cin>>f;     for(int ca=1;ca<=f;ca++){         cin>>t;         for(ll i=1;i<=t;i++)             cin>>b[i];         for(ll i=1;i<=t;i++)             cin>>a[i];         printf("Case %d: ",ca);         ll sum = china_reminder(t,a,b);         if(!sum){   /// 特判...             sum=1;             for(int i=1;i<=t;i++)                 sum = sum*b[i]/gcd(sum,b[i]);         }         printf("%lld\n",sum);     }}