详解几个方程的解法

一：求解不定方程 ax+by=c

　　先给出两个介绍扩展欧几里德很详细的网站：

　　http://www.cnblogs.com/frog112111/archive/2012/08/19/2646012.html

　　http://www.cnblogs.com/ka200812/archive/2011/09/02/2164404.html

扩展欧几里德的递归代码：

//扩展欧几里德的递归代码：
int exgcd(int a,int b,int &x,int &y)
{
    if(b==0)
    {
    //b==0，也就是说gcd(a,0)==a。
    //原式变为ax+by==a --> x==1,y==0。
        x=1;
        y=0;
        return a;
    }
    int d=exgcd(b,a%b,x,y);
    int t=x;
    x=y;   //更新当前的x
    y=t-a/b*y;   //更新当前的y
    return d; //返回最大公约数
}

　　扩展欧几里德的应用：

　　（一）：求解不定方程：

　　利用扩展欧几里德，可以求出

　　ax+by=gcd(a,b)的一组解，x1、y1

　　则通解为：

　　x = x1+b/gcd*i

　　y = y1-a/gcd*i

　　那么对于一般的不定式ax+by=c，

　　若gcd(a,b)整除c，则方程有解

　　一组解为：x0=x1*(c/gcd(a,b)),y0=y1*(c/gcd(a,b))。

　　通解为：x=x0+b/gcd*i,y=y0-a/gcd*i。

　　否则无解。

　　若ax+by=c有解，且gcd(a,b)=d，x的特解为x*(c/d)。

　　令ans=x*(c/d),s=b/d;

　　则方程的x最小整数解为：(ans%s+s)%s;

　　（二）：求解同余方程

　　ax+by=c，实际上它等价于同余方程，ax =c mod b,这样同余方程就可以用扩展欧几里德算法求解。

　　（三）：求解模的逆元

　　另外如果a b互质那么ax+by=gcd(a,b) = 1,转换为同余方程ax = 1 mod b这样x就变成了a的关于mod b的逆元。

　　求解扩展欧几里德过后，只需x=(x%b+b)%b，求出最小正整数解x，即为a的逆元。

　　也就是说求逆元也不过是扩展欧几里德算法的一个特殊情况。

　　两个简单的求解不定方程的例子：

　　POJ 2115  C Looooops

#include <iostream>
#include <stdio.h>
using namespace std;
long long A,B,C;
int k;
long long exgcd(long long a,long long b,long long &x,long long &y) {
    if(b==0) {
        x=1;
        y=0;
        return a;
    }
    long long d=exgcd(b,a%b,x,y);
    long long tmp=x;
    x=y;
    y=tmp-(a/b)*y;
    return d;
}
int main() {
    while(scanf("%I64d%I64d%I64d%d",&A,&B,&C,&k)!=EOF) {
        if(A==0 && B==0 && C==0 && k==0)
            break;
        long long a,b,c,x,y,d,s;
        a=C;
        b=(long long)1<<k;
        c=B-A;
        d=exgcd(a,b,x,y);
        if(c%d==0) {
            x=x*c/d;
            s=b/d;
            x=(x%s+s)%s;
            printf("%I64d
",x);
        } else
            printf("FOREVER
");

    }
    return 0;
}

　　POJ 2142  The Balance

http://www.cnblogs.com/chenxiwenruo/p/3557315.html

二：解方程x^2+y^2=z^2

　　设不定方程:x^2+y^2=z^2
　　若正整数三元组(x,y,z)满足上述方程，则称为毕达哥拉斯三元组。
　　若gcd(x,y,z)=1，则称为本原的毕达哥拉斯三元组。

　　定理：
　　正整数x,y,z构成一个本原的毕达哥拉斯三元组且y为偶数,当且仅当存在互素的正整数m,n(m>n),其中m,n的奇偶性不同,
　　并且满足
　　　　x=m^2-n^2,y=2*m*n, z=m^2+n^2

　　例题：

　　POJ 1305 Fermat vs. Pythagoras

　　http://www.cnblogs.com/chenxiwenruo/p/3558297.html

三：佩尔方程x^2-dy^2=1

　　形如x^2-dy^2=1 (d>1且d不为完全平方数)的不定方程称为佩尔方程。

 　　若d是完全平方数，那么x^2-dy^2=1是无解的。

　　

 

　　求解佩尔方程，可以用暴力求解法。

　　还有一种更高效的解法，连分数法，有兴趣的可以自己看看。

　　暴力解法代码：

void solve(int d){
    y=1;
    while(1){
        x=(long long)sqrt(double(d*y*y+1));
        if(x*x-d*y*y==1){
            break;
        }
        y++;
    }

}

　　给出两个简单例题：

　　POJ    1320  Street Numbers

　　求解两个不相等的正整数n和m (n<M)，使得1+2+3+…+n=n+(n+1)+…+m。

　　思路：实际上求解(2m+1)^2-8n^2=1。令x=2m+1，y=n。可以利用迭代公式求解出前十组的(n,m)。

 

　　HDU 3292  No more tricks, Mr Nanguo

　　求解X^2-NY^2=1的第K大的满足式子的X，如果无解则输出相应的句子。利用矩阵迭代+快速幂即可求出。

#include <iostream>
#include <cstdio>
#include <string.h>
#include <algorithm>
#include <math.h>
using namespace std;
const int mod=8191;
int n,k;
int vis[30];
long long x,y;
void solve(int d){
    y=1;
    while(1){
        x=(long long)sqrt(double(d*y*y+1));
        if(x*x-d*y*y==1){
            break;
        }
        y++;
    }

}
struct Matrix{
    long long m[2][2];
    void init(long long x,long long y){
        m[0][0]=x;m[0][1]=n*y%mod;
        m[1][0]=y;m[1][1]=x;
    }
    void initE(){
        memset(m,0,sizeof(m));
        m[0][0]=m[1][1]=1;
    }
}A;
Matrix operator*(Matrix a,Matrix b){
    Matrix c;
    for(int i=0;i<2;i++){
        for(int j=0;j<2;j++){
            c.m[i][j]=0;
            for(int k=0;k<2;k++){
                c.m[i][j]=(c.m[i][j]+a.m[i][k]*b.m[k][j]%mod)%mod;
            }
        }
    }
    return c;
}
Matrix quickPow(Matrix a,int b){
    Matrix ret;
    ret.initE();
    while(b){
        if(b&1)
            ret=ret*a;
        a=a*a;
        b=b>>1;
    }
    return ret;
}

int main()
{
    memset(vis,0,sizeof(vis));
    for(int i=1;i*i<30;i++){
        vis[i*i]=1;  //若是完全平方数，则无解，标记为1
    }
    while(scanf("%d%d",&n,&k)!=EOF){
        if(!vis[n]){
            solve(n);
            x=x%mod;y=y%mod;
            A.init(x,y);
            Matrix B;
            B=quickPow(A,k-1);
            long long xk,yk;
            xk=(B.m[0][0]*x%mod+B.m[0][1]*y%mod)%mod;
            printf("%I64d
",xk);
        }
        else{
            printf("No answers can meet such conditions
");
        }

    }
    return 0;
}

四：解高次同余方程A^x=B(mod C)

　　Baby Step Giant Step算法：复杂度O( sqrt(C) )

　　请戳：

　　http://www.cnblogs.com/chenxiwenruo/p/3554885.html