码迷，mamicode.com
首页 > 其他好文 > 详细

数论快速入门（同余、扩展欧几里德、中国剩余定理、大素数测定和整数分解、素数三种筛法、欧拉函数以及各种模板）

时间：2016-07-09 20:56:07      阅读：516      评论：0      收藏：0      [点我收藏+]

标签：

数学渣渣愉快的玩了一把数论，来总结一下几种常用的算法入门，不过鶸也是刚刚入门，

 所以也只是粗略的记录下原理，贴下模板，以及入门题目（感受下模板怎么用的）

（PS:文中蓝色字体都可以点进去查看百度原文）

附赠数论入门训练专题：点我打开专题（题目顺序基本正常，用以配套数论入门）

一、同余定理

简单粗暴的说就是：若 a-b == m 那么 a%m == b%m

这个模运算性质一眼看出。。。直接上入门水题：

Reduced ID Numbers

附AC代码（这个也没啥模板。。。。知道就好）

#include<iostream>
#include<cstdio>
#include<cstring>
#include<stdlib.h>
#define mem(a,x) memset(a,x,sizeof(a))
using namespace std;
typedef long long ll;
/*
同余：
if
    a-b == m
then
    a%m == b%m

*/
const int N = 1000000;
bool vis[N+5];
int a[N+5];
bool ol[N+5];
int main()
{
    int T;cin>>T;
    while (T--)
    {
        int n;
        scanf("%d",&n);
        mem(vis,0);
        for (int i = 1;i <= n;++i) scanf("%d",a+i);
        for (int i = 1;i <= n;++i)
        {
            for (int j = i+1;j <= n;++j)
            {
                int d = abs(a[i]-a[j]);
                vis[d] = 1;
            }
        }
        bool fd = 0;
        for (int m = 1;!fd;++m)
        {
            if (!vis[m])
            {
                mem(ol,0);bool ok = 1;
                for (int i = 1;ok&&i <= n;++i)
                {
                    if (ol[a[i]%m]) ok = 0;
                    else ol[a[i]%m] = 1;
                }
                if (ok)
                {
                    printf("%d\n",m);
                    fd = 1;
                }
            }
        }
    }
    return 0;
}

二、扩展欧几里德算法

这个扩展是从原欧几里德算法扩展而来，这个算法真心非常有用！非常有用！非常有用！（中国剩余定理也要用到它）

首先说欧几里德算法(其实就是我们小时候数学课上就学过的辗转相除法求gcd)

欧几里德说：gcd(a,b) = gcd(b,a%b)

于是得到欧几里德算法：

int gcd(int a,int b)
{
    return b==0?a:gcd(b,a%b);
}
该算法可以在几乎是log的时间算a,b的最大公因数gcd

PS:__gcd(a,b)库函数可以直接调用，但是有些OJ上提交会CE

现在，我们令a,b的最大公因数为gcd,那么我们一定可以找到一组x,y满足：

a*x + b*y = gcd;

此时，我们设x0,y0是上述不定方程的特解，那么就可以用x0,y0表示出整个不定方程的通解 

x = x0 + (b/gcd)*t;
y = y0 - (a/gcd)*t;

然后是怎么求解通解x,y?

倒过去看看欧几里德算法，显然，它的结束条件是 b = 0的时候返回a

这就意味着，终止状态是 ：a = gcd ,b = 0;

将这组a,b代会不定方程ax+by=gcd,可以得到一组特解：x = 1,y = 0

找到终止状态，我们再看看递归的过程：

gcd = b*x1 + (a%b)*y1 
    = b*x1 + (a-(a/b)*b)*y1 // a%b = a-(a/b)*b
    = b*x1 + a*y1 - (a/b)*b*y1
    = a*y1 + b*(x1-a/b*y1)

最后得到gcd = a*y1 + b*(x1-a/b*y1)和原来的式子gcd = a*x + b*y一对比就得到: 

x = y1,y = x1 - a/b*y1;

上面两个等式很重要！！（因为模板里面会直接用到，如果是直接背模板就得把它背熟）

下面就是扩展欧几里德算法：

int e_gcd(int a,int b,int &x,int &y)
{
    if (b == 0)
    {
        x = 1,y = 0;
        return a;
    }
    int ans = e_gcd(b,a%b,x,y);
    int tmp = x;
    x = y;
    y = tmp - a/b*y;
}

这个算法本质上还是求a,b的最大公因数，只是在计算的过程中顺带计算x,y(同时体验了一把引用的神奇)

扩展欧几里德算法有什么用？反正用法很多，它可以求解形如 ax+by=c的通解

但是实际应用中通常指要求在通解中选一些特殊的解，

比如一个数对于另一个数的乘法逆元。那么什么是乘法逆元？

ax ≡ 1(mod b) // ax % b = 1
这里，我们称x是a关于b的乘法逆元

ax%b = 1,可以化成 ax = by + 1

显然，它等价于这样的表达式：ax + by = 1             

这个式子很眼熟有木有！！！如果等式右边是gcd(a,b)就好了！！！

然后这里copy欧几里德的三个定理：

定理一:如果d = gcd(a,b),则必能找到正的或负的整数x和y使d = ax + by 
定理二:若gcd(a,b) = 1,则方程ax≡c(mod b)在[0,b-1]上有唯一解
定理三:若gcd(a,b) = d,则方程ax≡c(mod b)在[0,b|d01]上有唯一解

于是，对上述方程 ax + by = 1, 当gcd(a,b) != 1的时候无解

故方程ax + by = c 有解的条件是 ： c%gcd(a,b) = 0;

按乘法逆元讲，一般，我们能找到无数组解满足条件，但一般题目是求解最小的那组解

假设我们求出了特解x0,那么 ，只需要用x0 % b就是最小解了

为什么？（这个我没管，反正知道是这个。。。。后面贴的参考资料链接里面有。。。）

另外，有的时候求得的特解可能是个负数，或者说b是负数，怎么办？

如果b是负数，取b的绝对值

如果x0是负数，让x0对|b|取模后再加上|b|

然后，直接上模板代码：

int Cal(int a,int b)//求最小的x使ax+ by = 1
{
    int x,y;
    int gcd = e_gcd(a,b,x,y);
    if (1%gcd) return -1;//无解
    x*=1/gcd;
    b = abs(b);
    int ans = x%b;
    if (ans <= 0) ans += b;
    return ans;
}


下面给出求最小逆元的代码： 

typedef long long ll;
ll e_gcd (ll a, ll b, ll& x, ll& y)
{
    if (b == 0)
    {
        x = 1, y = 0;
        return a;
    }
    ll ans = e_gcd (b, a % b, y, x);
    y -= a / b * x; //这个和前面用的方法不一样，不过是对的，写起来更快、
    return ans;
}
ll Cal(ll a,ll b,ll c)//求最小的x使ax+ by = c
{
    ll x,y;
    ll gcd = e_gcd(a,b,x,y);
    if (c%gcd) return -1;//无解
    x*=c/gcd;
    b /= gcd;
    if (b < 0) b = -b;
    ll ans = x%b;
    if (ans <= 0) ans += b;
    return ans;
}

来一道入门级的裸的扩展欧几里德的题感受下模板怎么用：

PS:方程还是要自己列，然后用扩展欧几里德求解

青蛙的约会

AC代码供参考：

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#define mem(a,x) memset(a,x,sizeof(a))
using namespace std;
typedef long long ll;
ll e_gcd(ll a,ll b,ll &x,ll &y)
{
    if (b == 0)
    {
        x = 1,y = 0;
        return a;
    }
    ll ans = e_gcd(b,a%b,x,y);
    ll tmp = x;
    x = y;
    y = tmp - a/b*y;
    return ans;
}
ll cal(ll a,ll b,ll c)//求最小的x使ax+by=c
{
    ll x,y;
    ll gcd = e_gcd(a,b,x,y);
    if (c%gcd != 0) return -1;
    x *= c/gcd;
    b/=gcd;
    if (b < 0) b = -b;
    ll ans = x%b;
    if (ans <= 0) ans += b;
    return ans;
}
int main()
{
    ll xa,xb,va,vb,L;
    while (~scanf("%lld %lld %lld %lld %lld",&xa,&xb,&va,&vb,&L))
    {
        ll ans = cal(vb-va,L,xa-xb);
        if (ans == -1) puts("Impossible");
        else printf("%lld\n",ans);
    }
    return 0;
}

三、中国剩余定理

基本看到一个算法都是用外国人的名字命名，而中国剩余定理又称孙子定理，没错，就是中国的孙子！（这话怎么怪怪的(⊙o⊙)......）
《孙子算经》里面好像有个求同余方程组的问题，不过古文就不看了。
用现代数学的语言来说，中国剩余定理给出了以下一元线性同余方程组：

    x ≡ a1 (mod m1)    x%m1 = a1     
    x ≡ a2 (mod m2)    x%m2 = a2
(S):x ≡ a3 (mod m3) -> x%m3 = a3
    ...                ...
    x ≡ an (mod mn)    x%mn = an

有解的判定条件，并用构造法给出了在有解的情况下解的具体形式
中国剩余定理说明：假设整数m1,m2, ... ,mn两两互质，则对任意的整数：a1,a2, ... ,an，方程组（S）有解，并且通解可以用如下方式构造得到：
设
 M = m1*m2*m3......*mn = ∏mi(i从1-n)
是整数m1,m2, ... ,mn的乘积，并设
 Mi = M/mi(i从1到n) 
是除了mi以外的n- 1个整数的乘积。
设
 ti = Mi^-1
为
 Mi
模
 mi
的数论倒数
 ti*Mi ≡1（mod mi）(i从1到n) 
：
方程组
 （S）
的通解形式为
 x = a1t1M1 + a2t2M2 + ... + antnMn + kM = kM + ∑aitiMi,(i = 1,2,3,...,n,k ∈ Z) 
：在模
 M 
的意义下，方程组
 (S) 
只有一个解：
 x = ∑aitiMi,(i = 1,2,3,...,n)
人生苦短，暂不证明（其实鶸不会证明。。。）
直接上模板代码：
typedef long long ll;
ll e_gcd (ll a, ll b, ll& x, ll& y)
{
    if (b == 0)
    {
        x = 1, y = 0;
        return a;
    }
    ll ans = e_gcd (b, a % b, y, x);
    y -= a / b * x; //这个和前面用的方法不一样，不过是对的，写起来更快、
    return ans;
}
ll CR(int a[],int m[],int n)
{
    ll M = 1;
    for (int i = 1;i <= n;++i) M*=m[i];
    ll ans = 0;
    for (int i = 1;i <= n;++i)
    {
        ll Mi = M/m[i];ll x,y;
        ll t = e_gcd(m[i],Mi,x,y);
        ans = (ans + y*Mi*a[i])%M;
    }
    return (M+ans%M)%M;
}

=-=这个是我最初的模板，当然他非常渣，在经历了WA无数之后，模板得到了很好的改进，
就不直接贴了，直接上入门题目和AC代码（AC代码用的改进版模板=-=）
题目：Hello Kiki 
代码：
#include<iostream>
#include<algorithm>
#include<cstdio>
#include<cstring>
#include<cmath>
#define mem(a,x) memset(a,x,sizeof(a))
using namespace std;
typedef long long ll;
ll e_gcd (ll a, ll b, ll& x, ll& y)
{
    if (b == 0)
    {
        x = 1, y = 0;
        return a;
    }
    ll ans = e_gcd (b, a % b, y, x);
    y -= a / b * x;
    return ans;
}
ll CR (int a[], int m[], int n)
{
    ll Mi = m[1], ans = a[1];
    for (int i = 2; i <= n; ++i)
    {
        ll mi = m[i], ai = a[i];
        ll x, y;
        ll gcd = e_gcd (Mi, mi, x, y);
        ll c = ai - ans;
        if (c % gcd != 0) return -1;
        ll M = mi / gcd;
        ans += Mi * ( ( (c /gcd*x) % M + M) % M);
        Mi *= M;
    }
    if (ans == 0) //当余数都为0
    {
        ans = 1;
        for (int i = 1; i <= n; ++i)
        {
            ans = ans*m[i]/__gcd(ans,(ll)m[i]);
        }
    }
    return ans;
}
int main()
{
    int T;cin>>T;int kas = 0;
    while (T--)
    {
        int n,a[100],m[100];
        scanf("%d",&n);
        for (int i = 1;i <= n;++i) scanf("%d",m+i);
        for (int i = 1;i <= n;++i) scanf("%d",a+i);
        printf("Case %d: %lld\n",++kas,CR(a,m,n));
    }
    return 0;
}

不爽再来一题？Biorhythms
这题要想清楚同余方程组怎么用的，本鶸开始没想清楚直接套模板使劲WA。。。
AC代码：
#include<iostream>
#include<cstdio>
#include<cstring>
#include<stdlib.h>
#include<cmath>
#include<algorithm>
#define mem(a,x) memset(a,x,sizeof(a))
using namespace std;
typedef long long ll;
ll e_gcd(ll a,ll b,ll &x,ll &y)
{
    if (b == 0)
    {
        x = 1,y = 0;
        return a;
    }
    ll ans = e_gcd(b,a%b,y,x);
    y -= a/b*x;
    return ans;
}
ll CR(int a[],int m[],int n)
{
    ll M = 1;
    for (int i = 1;i <= n;++i) M*=m[i];
    ll ans = 0;
    for (int i = 1;i <= n;++i)
    {
        ll Mi = M/m[i]; ll x,y;
        ll t = e_gcd(m[i],Mi,x,y);
        ans = (ans+y*Mi*a[i])%M;
    }
    ans = (M+ans%M)%M;
    if (ans == 0) //当余数都为0
    {
        ans = 1;
        for (int i = 1; i <= n; ++i)
        {
            ans = ans*m[i]/__gcd(ans,(ll)m[i]);
        }
    }
    return ans;
}
int main()
{
    int a[4];
    int m[4] = {0,23,28,33};int d;int kas = 0;
    while (~scanf("%d %d %d %d",a+1,a+2,a+3,&d))
    {
        if (a[1]==-1&&a[2]==-1&&a[3]==-1&&d==-1) break;
        for (int i = 1;i <= 3;++i)
        {
            a[i] = a[i] - d;
            while (a[i]<0) a[i]+=m[i];
        }
        ll ans = CR(a,m,3);
        printf("Case %d: the next triple peak occurs in %lld days.\n",++kas,ans);
    }
    return 0;
}

四、素数筛法
很多数学题里面都有用到，主要就是卡超时的问题。
不多说，直接上三种模板：
素数筛法一：（最简单的）
#include<iostream>
#include<cstring>
#include<cstdio>
#define mem(a,x) memset(a,x,sizeof(a))
#define inf (1<<29)
using namespace std;
typedef long long ll;
const int N = 100;
bool p[N+5];
void init()
{
    mem(p,1);
    for (int i = 2;i*i <= N;++i)
    {
        if (p[i])
        {
            for (int j = i*i;j <= N;j+=i)
            {
                p[j] = 0;
            }
        }
    }
}
int main()
{
    init();
    for (int i = 2;i <= N;++i)
    {
        if (p[i]) cout<<i<<endl;
    }
    return 0;
}
素数筛法二：
#include<iostream>
#include<cmath>
#include<cstdio>
#define mem(a,x) memset(a,x,sizeof(a))
#define inf (1<<29)
using namespace std;
typedef long long ll;
const int N = 100;
int p[N+5];
void init()
{
    int sq = (int)sqrt(N*2+1);
    for (int i = 3;i <= sq;i+=2)
    {
        if (p[i>>1]) continue;
        for (int j = i*i;j <= N<<1;j+=i<<1)
        {
            p[j>>1] = 1;
        }
    }
}
int main()
{
    init();
    puts("2");
    for (int i = 1;i < N;++i)
    {
        if (p[i]) continue;
        printf("%d\n",(i<<1)+1);
    }
    return 0;
}

素数筛法三：
#include<iostream>
#include<cstring>
#include<cstdio>
#define mem(a,x) memset(a,x,sizeof(a))
#define inf (1<<29)
using namespace std;
typedef long long ll;
const int N = 100000;
int a[N+5];
int p[N+5];
void init()
{
    int t = 0;
    for (int i = 2;i <= N;++i)
    {
        if (a[i] == 0)
        {
            p[++t] = i;
        }
        for (int j = 1,k;(j<=t)&&(k=i*p[j])<=N;++j)
        {
            a[k] = 1;
            if (i%p[j]==0) break;
        }
    }
}
int main()
{
    init();
    for (int i = 1;p[i]>1;++i)
    {
        printf("%d\n",p[i]);
    }
    return 0;
}

本来素数筛完应该直接上欧拉函数，不过觉得大素数测定蛮好玩的，于是插一句大素数测定与整数分解
五、大素数测定与整数分解：
上面的素数筛法很常用，但是如果问你一个很大的数n（n < 2^63）是不是素数，怎么办？
于是引入了Miller Rabin算法
该算法是一种基于概率的素数测试算法，特点是速度快，能判断<2^63的数是不是素数
虽然是基于概率，但是其实该算法还是蛮可靠的，首先是可以通过增加测试次数提高测试的正确率
另外，我自己玩的时候，将测试次数调为1，一直测，没有出现判断错误的情况（可能我人品太好(⊙o⊙)？）
原理大致看了下，表示不是很懂，直接上代码
#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<algorithm>
#include<queue>
#include<cmath>
#include<stdlib.h>
#include<cctype>
#include<time.h>
#define mem(a,x) memset(a,x,sizeof(a))
using namespace std;
typedef long long ll;
const int S = 8;//测试次数
ll mult_mod (ll a,ll b, ll c)
{
    a%=c,b%=c;
    ll ret = 0;
    ll tmp = a;
    while (b)
    {
        if (b&1)
        {
            ret += tmp;
            if (ret > c) ret -= c;
        }
        tmp<<=1;
        if (tmp>c) tmp-=c;
        b>>=1;
    }
    return ret;
}
ll pow_mod(ll a,ll n,ll mod)
{
    ll ret = 1;
    ll temp = a%mod;
    while (n)
    {
        if (n&1) ret = mult_mod(ret,temp,mod);
        temp = mult_mod(temp,temp,mod);
        n>>=1;
    }
    return ret;
}
bool check(ll a,ll n,ll x,ll t)
{
    ll ret = pow_mod(a,x,n);
    ll last = ret;
    for (int i = 1;i <= t;++i)
    {
        ret = mult_mod(ret,ret,n);
        if (ret == 1&&last!=1&&last!=n-1) return 1;
        last = ret;
    }
    if (ret != 1) return 1;
    else return 0;
}
bool MR(ll n)
{
    if(n < 2) return 0;
    if (n == 2) return 1;
    if ((n&1)==0) return 0;
    ll x = n - 1;
    ll t = 0;
    while ((x&1)==0)
    {
        x>>=1;++t;
    }
    srand(time(NULL));
    for (int i = 0;i < S;++i)//做S次测试
    {
        ll a = rand()%(n-1) + 1;
        if (check(a,n,x,t)) return 0;//只要其中有一次判定是合数就可以确定一定是合数
    }
    return 1;
}
int main()
{
    ll n;
    while (cin>>n)
    {
        if (MR(n)) puts("YES");
        else puts("NO");
    }
    return 0;
}
水题：Prime Test
这个直接裸模板，就不贴AC代码了
题外话：表示上面算法蛮好玩，于是某次比赛用该算法打素数表，然后就。。。。。呵呵。。。。。所以觉得这个算法主要就是针对很大的数的素数判断的，如果没有这个条件，请贪玩的童鞋不要在比赛随意使用。。。。。。。。
然后整数分解用的是Pollard rho算法,具体原理我也没懂，反正就是可以求2^63以内的数分解成素因子
模板是直接copy的kuangbin的模板（蛮长的。。。。）
#include<stdio.h>
#include<string.h>
#include<stdlib.h>
#include<time.h>
#include<iostream>
#include<algorithm>
using namespace std;


//****************************************************************
// Miller_Rabin 算法进行素数测试
//速度快，而且可以判断 <2^63的数
//****************************************************************
const int S=20;//随机算法判定次数，S越大，判错概率越小
//计算 (a*b)%c.   a,b都是long long的数，直接相乘可能溢出的
//  a,b,c <2^63
long long mult_mod(long long a,long long b,long long c)
{
    a%=c;
    b%=c;
    long long ret=0;
    while(b)
    {
        if(b&1){ret+=a;ret%=c;}
        a<<=1;
        if(a>=c)a%=c;
        b>>=1;
    }
    return ret;
}
//计算  x^n %c
long long pow_mod(long long x,long long n,long long mod)//x^n%c
{
    if(n==1)return x%mod;
    x%=mod;
    long long tmp=x;
    long long ret=1;
    while(n)
    {
        if(n&1) ret=mult_mod(ret,tmp,mod);
        tmp=mult_mod(tmp,tmp,mod);
        n>>=1;
    }
    return ret;
}

//以a为基,n-1=x*2^t      a^(n-1)=1(mod n)  验证n是不是合数
//一定是合数返回true,不一定返回false
bool check(long long a,long long n,long long x,long long t)
{
    long long ret=pow_mod(a,x,n);
    long long last=ret;
    for(int i=1;i<=t;i++)
    {
        ret=mult_mod(ret,ret,n);
        if(ret==1&&last!=1&&last!=n-1) return true;//合数
        last=ret;
    }
    if(ret!=1) return true;
    return false;
}

// Miller_Rabin()算法素数判定
//是素数返回true.(可能是伪素数，但概率极小)
//合数返回false;

bool Miller_Rabin(long long n)
{
    if(n<2)return false;
    if(n==2)return true;
    if((n&1)==0) return false;//偶数
    long long x=n-1;
    long long t=0;
    while((x&1)==0){x>>=1;t++;}
    for(int i=0;i<S;i++)
    {
        long long a=rand()%(n-1)+1;//rand()需要stdlib.h头文件
        if(check(a,n,x,t))
            return false;//合数
    }
    return true;
}


//************************************************
//pollard_rho 算法进行质因数分解
//************************************************
long long factor[100];//质因数分解结果（刚返回时是无序的）
int tol;//质因数的个数。数组小标从0开始

long long gcd(long long a,long long b)
{
    if(a==0)return 1;//???????
    if(a<0) return gcd(-a,b);
    while(b)
    {
        long long t=a%b;
        a=b;
        b=t;
    }
    return a;
}

long long Pollard_rho(long long x,long long c)
{
    long long i=1,k=2;
    long long x0=rand()%x;
    long long y=x0;
    while(1)
    {
        i++;
        x0=(mult_mod(x0,x0,x)+c)%x;
        long long d=gcd(y-x0,x);
        if(d!=1&&d!=x) return d;
        if(y==x0) return x;
        if(i==k){y=x0;k+=k;}
    }
}
//对n进行素因子分解
void findfac(long long n)
{
    if(Miller_Rabin(n))//素数
    {
        factor[tol++]=n;
        return;
    }
    long long p=n;
    while(p>=n)p=Pollard_rho(p,rand()%(n-1)+1);
    findfac(p);
    findfac(n/p);
}

int main()
{
    //srand(time(NULL));//需要time.h头文件//POJ上G++不能加这句话
    long long n;
    while(scanf("%I64d",&n)!=EOF)
    {
        tol=0;
        findfac(n);
        for(int i=0;i<tol;i++)printf("%I64d ",factor[i]);
        printf("\n");
        if(Miller_Rabin(n))printf("Yes\n");
        else printf("No\n");
    }
    return 0;
}

整数分解也有个题：GCD & LCM Inverse（用dfs去找答案）
AC代码：
#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<algorithm>
#include<queue>
#include<cmath>
#include<stdlib.h>
#include<cctype>
#define mem(a,x) memset(a,x,sizeof(a))
using namespace std;
typedef unsigned long long ll;
const int S = 20;
ll mult_mod (ll a, ll b, ll c)
{
    a %= c, b %= c;
    ll ret = 0;
    ll tmp = a;
    while (b)
    {
        if (b & 1)
        {
            ret += tmp;
            if (ret > c) ret -= c;
        }
        tmp <<= 1;
        if (tmp > c) tmp -= c;
        b >>= 1;
    }
    return ret;
}
ll pow_mod (ll a, ll n, ll mod)
{
    ll ret = 1;
    ll temp = a % mod;
    while (n)
    {
        if (n & 1) ret = mult_mod (ret, temp, mod);
        temp = mult_mod (temp, temp, mod);
        n >>= 1;
    }
    return ret;
}
bool check (ll a, ll n, ll x, ll t)
{
    ll ret = pow_mod (a, x, n);
    ll last = ret;
    for (int i = 1; i <= t; ++i)
    {
        ret = mult_mod (ret, ret, n);
        if (ret == 1 && last != 1 && last != n - 1) return 1;
        last = ret;
    }
    if (ret != 1) return 1;
    else return 0;
}
bool MR (ll n)
{
    if (n < 2) return 0;
    if (n == 2) return 1;
    if ( (n & 1) == 0) return 0;
    ll x = n - 1;
    ll t = 0;
    while ( (x & 1) == 0)
    {
        x >>= 1;
        ++t;
    }
//    srand(time(NULL));
    for (int i = 0; i < S; ++i) //做S次测试
    {
        ll a = rand() % (n - 1) + 1;
        if (check (a, n, x, t) ) return 0; //只要其中有一次判定是合数就可以确定一定是合数
    }
    return 1;
}
ll fac[11111];//质因数分解结果(刚返回时时无序的)
int tot;//质因数的个数，数组下标从0开始
ll gcd (ll a, ll b)
{
    if (a == 0) return 1;
    if (a < 0) return gcd (-a, b);
    while (b)
    {
        ll t = a % b;
        a = b, b = t;
    }
    return a;
}
ll PR (ll x, ll c)
{
    ll i = 1, k = 2;
    ll x0 = rand() % x;
    ll y = x0;
    while (1) //美丽的循环，不会死
    {
        ++i;
        x0 = (mult_mod (x0, x0, x) + c) % x;
        ll d = gcd (y - x0, x);
        if (d != 1 && d != x) return d;
        if (y == x0) return x;
        if (i == k)
        {
            y = x0;
            k += k;
        }
    }
}
void findfac (ll n) //对n进行素因子分解
{
    if (MR (n) ) //如果n是素数
    {
        fac[tot++] = n;
        return;
    }
    ll p = n;
    while (p >= n) p = PR (p, rand() % (n - 1) + 1);
    findfac (p);
    findfac (n / p);
}
ll x[111];
ll k;
ll a, b, mn;
ll ans;
ll g, lcm;
const ll inf = 1LL << 62LL;
bool fd;
void dfs (ll cur, ll p)
{
    ll q = ans / p;
    if (gcd (p, q) == 1)
    {
        fd = 1;
        ll tmp = p * g + q * g;
        if (mn > tmp)
        {
            mn = tmp;
            a = p*g;
            b = q*g;
        }
    }
    if (cur > k) return;
    dfs(cur+1,p);p*=x[cur];
    if (p > mn ) return ;
    dfs(cur+1,p);
}
int main()
{
    while (~scanf ("%llu %llu", &g, &lcm) )
    {
        if (g == lcm)
        {
            printf("%llu %llu\n",g,lcm);
            continue;
        }
        ans = lcm / g;
        tot = 0;
//        cout << ans << endl;
        findfac (ans);
        sort (fac, fac + tot); //fac保存的是ans的素因子，比如3 60对应的fac数组是2 2 5
//        for (int i = 0; i < tot; ++i) cout << fac[i] << " ";
//        cout << "---------------------------------------------------" << endl;
        k = 0;
        x[0] = fac[0];
        for (int i = 1; i < tot; ++i)
        {
            if (fac[i] == fac[i - 1]) x[k] *= fac[i];
            else x[++k] = fac[i];
        }
        sort (x, x + k + 1); //x保存的是所有不相同的素因子，比如3 60 对应的x数组是4 5
//        for (int i = 0; i <= k; ++i) cout << x[i] << " ";
//        cout << endl;
        //用dfs将数组x分成2部分p*q,满足p,q互质,找到所有p,q中使a+b最小的情况，其中a = p*lcm,b = q*lcm
        //比如3 60 ans = 20 4 5 ,a = 4*3 b = 5*3
        mn = inf;
        fd = 0;
        dfs (0, 1);
        if (a > b) swap (a, b);
//        if (!fd) puts ("???");
        printf ("%llu %llu\n", a, b);
    }
    return 0;
}
/*
7 635040

*/


六、欧拉函数
定义：
在数论，对正整数n，欧拉函数是小于等于n的数中与n互质的数的数目。
此函数以其首名研究者欧拉命名(Euler‘so totient function)，它又称为Euler‘s totient function、φ函数、欧拉商数等。
 例如φ(8)=4，因为1,3,5,7均和8互质。 
通式：
 φ(x) = x(1 - 1/p1)(1 - 1/p2)(1 - 1/p3) ...(1 - 1/pn)
,其中p1, p2……pn为x的所有质因数，x是不为0的整数。φ(1)=1（唯一和1互质的数(小于等于1)就是1本身）。 (注意：每种质因数只一个。比如12=2*2*3那么φ（12）=12*（1-1/2）*(1-1/3)=4
若n是质数p的k次幂，
 φ(n) = p^k - p^(k-1) = (p-1)p^(k-1)
，因为除了p的倍数外，其他数都跟n互质。
设n为正整数，以 φ(n)表示不超过n且与n互
素的正整数的个数，称为n的欧拉函数值，这里函数
φ：N→N，n→φ(n)称为欧拉函数。
欧拉函数是积性函数——若m,n互质，
 φ(mn) = φ(m)φ(n);
特殊性质：当n为奇数时，
 φ(2n) = φ(n); 
, 证明与上述类似。
若n为质数则
 φ(n) = n - 1;
数学渣渣，再次直接上模板：
#include<iostream>
#include<cmath>
#include<algorithm>
using namespace std;

const int MAXN = 1e5;

int isprime[MAXN];
int prime[MAXN];
int cnt;
void getP()
{
    cnt = 0;
    for(int i = 1; i < MAXN; i++)
        isprime[i] = 1;
    for(int i = 2; i < MAXN; i++)
    {
        if(!isprime[i])continue;
        prime[cnt++] = i;
        for(int j = 2 * i; j < MAXN; j += i)
        {
            isprime[j] = 0;
        }
    }
}

int euler( int n )  //求小于n且与n互质的数的个数
{
    int ans = n;

    for(int i = 0; prime[i] * prime[i] <= n; i++)
    {
        if(n%prime[i] == 0)
        {
            ans = ans - ans / prime[i]; //ans=ans*(1-1/pi)
            while(n%prime[i] == 0)
            {
                n /= prime[i];
            }
        }
    }
    if(n > 1)  ans = ans - ans / n;
    return ans;
}
int main()
{
    getP();
    int n;
    while(cin >> n&&n)
    {
        cout << euler( n ) << endl;
    }
    return 0;
}

七、其他
入门就先到这里，虽然基本都是套模板，但是还是要把模板运用熟练，
接下来就是多做些题，多体会各种模板在各种题目中的灵活运用，
可能很多现在不是很懂的地方等做的题多了，自然而然懂了，加油↖(^ω^)↗~~~
错误的地方望指出~~~
参考资料：
1.各种百度百科，阅读过程点击蓝色字跳转

数论快速入门（同余、扩展欧几里德、中国剩余定理、大素数测定和整数分解、素数三种筛法、欧拉函数以及各种模板）

标签：

原文地址：http://blog.csdn.net/tomorrowtodie/article/details/51865496
踩
(0)
赞
(0)
   
举报
评论 一句话评论（0）
登录后才能评论！
分享档案
更多>
2021年07月29日 (22)
2021年07月28日 (40)
2021年07月27日 (32)
2021年07月26日 (79)
2021年07月23日 (29)
2021年07月22日 (30)
2021年07月21日 (42)
2021年07月20日 (16)
2021年07月19日 (90)
2021年07月16日 (35)
周排行
mamicode.com排行更多图片
更多
友情链接
兰亭集智   国之画   百度统计   站长统计   阿里云   chrome插件   新版天听网
关于我们 - 联系我们 - 留言反馈
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码！