【总结】扩展欧几里得算法

前言

---

    最近的一连串考试中不幸考到了自己学习得最不扎实的部分——数学,还是想一点点把那些无论基础还是高深的知识点捡起来 ……

概念

    扩展欧几里德算法，在已知a, b求解一组x，y，使它们满足贝祖等式： ax+by = gcd(a, b) =d，常用于求解模线性方程及方程组。

　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　——百度百科
    贴上觉得简洁易懂的讲解链接http://blog.csdn.net/acmore_xiong/article/details/47694909
    代码

void Exgcd(int a,int b,int &x,int &y,int &g) {
    if(b == 0) {
        g = a; x = 1, y = 0; return;
    }
    Exgcd(b,a%b,y,x,g);
    y -= (a/b)*x;
}

    然后来看求解不定整数方程的应用吧
    对于 ax+by = c (a,b不全为0)
    g = gcd(a,b)
    若方程无解 c % g != 0
    有解时: c’ = c/g
    然后我们可以用上面的代码解 ax + by = g
    得到一组特殊解 $x_0,y_0$, 于是 $x_1 = x_0*c/g$, $y_1 = y_0*c/g$就是原方程的一组特殊解啦
    通解可以表示为 $x_i = x_1 + t*b/g$,$y_i = y_1 - t*a/g$ (t∈Z)

例题

---

    青蛙的约会这种经典题不多说
    上昨天考试题$Math$

题意

    T 组数据,每组问二项式(x^(a/b)+x^(c/d))^n 展开后有多少有理项,又有多少方案使得有理项不相邻,方案数%p.

数据范围

    $|a|,|b|,|c|,|d|,T< =10^5,n< =10^9, 0=< p< =10^7$

思路:

    60 分算法:
    暴力计算有理项数量 x,则 ans = x!*(N+1-x)!*C(N+2-x,x),费马小定理求逆元即可再搞组合数.注意快速幂要先将底数%Mod 否则易炸.
    对C(N+2-x,x)的理解:画图,有N+1-x个无理项,为满足要求,x个有理项只能放在两端或无理项之间的共N+2-x个位置上.
    满分算法:
    将(x^(a/b)+x^(c/d))^n 展开共 n+1 项,每一项可以表示为 x^((a/b)k+(c/d)(n-k))(k∈[0,n])= x^(((ad-bc)*k+bcn)/bd).
    于是问题按照如下转化
    $(ad-bc)*k+bcn \equiv 0 (mod(bd))->$
    $(bc-ad)*k+bdy \equiv bcn(mod(bd))->$
    $(bc-ad)*k+bdy \equiv bcn$
    可以使用扩展欧几里得求解 k,计算 k∈[0,n]的可行解个数.
    记可行解个数为 A,则 Ans = A! * (n+1-A)! * C(n+2-A,A). 组合数部分费马小定理+逆元计算,但是要注意 C 应单独写函数判断两种情况,对于 C(m,n): 1)m < n 2)n==0 ,否则在遇到模数P==n+2-A 而 A==0 的情况答案错误.

代码

60分暴力

#include<cstdio>
using namespace std;
typedef long long LL;
const int sm = 8e3 + 50;

int T,Mod;

LL Ans;
LL Fac[sm];

LL QMod(LL a,int b) {
    LL ret = 1, x = a % Mod;
    while(b) {
        if(b&1) ret = 1ll*ret*x%Mod;
        x = 1ll*x*x%Mod;
        b >>= 1;
    }
    return ret;
}

int main() {
    freopen("math.in","r",stdin);
    freopen("math.out","w",stdout);
    LL x,y; 
    LL N,A,B,C,D,Up,Down,Tmp;

    scanf("%d%d",&T,&Mod);

    Fac[0] = 1;
    for(int i = 1; i <= 8010; ++ i)
        Fac[i] = 1ll* Fac[i-1] * i % Mod;

    for(int i = 1; i <= T; ++ i) {
        scanf("%lld%lld%lld%lld%lld",&N,&A,&B,&C,&D);

        x = 0;

        for(int j = 0; j <= N; ++ j) {
            Up = 1ll*A*D*(N-j)+1ll*C*B*j;
            Down = 1ll*B*D;
            if(Up % Down == 0)++x;
        }
        y = N + 1 - x;

        if(N-2*x+2 <= 0) Ans = 0;
        else {
            Tmp = 1ll*Fac[N+2-x]*QMod(1ll*Fac[x]*Fac[N+2-2*x],Mod-2)%Mod;
            Ans = 1ll*Fac[x]*Fac[y]%Mod;
            Ans = Ans*Tmp%Mod;
        }
        printf("%lld %lld\n",x,Ans);
    }
    return 0;
}

AC代码

#include<cstdio>
#include<cmath>
#include<algorithm>
using namespace std;

typedef long long LL;

const int sm = 1e7 + 10;
const int B = 15, sn = sm/B + 10;

int T,P;

int A;
LL n, Ans;

LL gcd(LL u,LL v) { return (v==0)?u:gcd(v,u%v); }

LL ret;
LL Qmod(LL a,LL b) {
    ret = 1; a %= P;
    while(b) {
        if(b&1) ret = ret*a%P;
        a = a*a%P;
        b>>=1;
    }
    return ret;
}

LL x,y;
void Exgcd(LL a,LL b,LL &x,LL &y,LL &g) {
    if(b == 0) {
        g = a; x = 1, y = 0; return;
    }
    Exgcd(b,a%b,y,x,g);
    y -= (a/b)*x;
}

LL Abs,g;
LL a,b,c,d;
bool flag;
int Calc(LL aa,LL bb) { 

    Exgcd(aa,bb,x,y,g);

    if(b*c*n%g) return ret = 0; 

    x%=(b*d);

    x = x*b % abs(b*d) *c % abs(b*d) *n % abs(b*d)/g;

    x = (x % abs(b*d/g) + abs(b*d/g) ) % abs(b*d/g);

    ret = ((x<=n+1)?1:0) + (n-x)/(abs(b*d/g));

    return ret;
}

LL fac[sn];
void Init() {
    LL tmp = 1, Ct = 0, cnt = 0;
    fac[0] = 1;
    for(int i = 1; i <= 1e7; ++i) {
        if(cnt%B == 0) fac[++Ct] = fac[Ct-1] * i % P, cnt = 1;
        else fac[Ct] = fac[Ct] * i % P, ++cnt;
    } 
}

LL t,q,l;
LL Fac(LL u) {
    t = u/B; q = 1;
    if(u%B!=0) {
        l = u%B;
        for(int i = 1; i <= l; ++ i)    
            q = q * (t*B+i) % P; 
    }
    return fac[t] * q % P; 
}

LL C(LL m,LL n) {
    if(n==0)  return 1;
    if(n> m)  return 0;
    return Fac(m)*Qmod(Fac(n)*Fac(m-n)%P,P-2)%P;
}

int main() {

    freopen("math.in","r",stdin);
    freopen("math.out","w",stdout);

    scanf("%d%d",&T,&P);

    Init();

    for(int i = 1; i <= T; ++i){
        scanf("%lld%lld%lld%lld%lld",&n,&a,&b,&c,&d);

        A = Calc(b*c-a*d,b*d);

        if(A >= P || n + 1 - A >= P || A > n + 2 - A) Ans = 0;
        else Ans = Fac(A) * Fac(n+1-A) % P * C(n+2-A,A) % P;

        printf("%d %lld\n",A,Ans);
    }
    return 0;
}