FZU code carnival-GCDLCM (最大公约数与最小公倍数）（Miller Rabin质数判定）

题目描述

In FZU ACM team, BroterJ and Silchen are good friends, and they often play some interesting games.

在福建大学ACM队伍里，J大哥和SilchenS是好朋友，他们经常一块玩有趣的游戏。
One day they play a game about GCD and LCM. firstly BrotherJ writes an integer A and Silchen writes an integer B on the paper. Then BrotherJ gives Silchen an integer X. Silchen will win if he can find two integers Y1 and Y2 that satisfy the following conditions:

有一天他们玩一个关于最大公约数与最小公倍数的游戏。首先J大哥写了一个整数A，Silchen写了一个整数B在纸上。接着J大哥给了Silchen一个整数X。Silchen 如果可以找到两个满足下列条件的整数Y1，Y2的话就能获胜。
• GCD(X, Y1) = A

X与Y1的最大公约数为A
• LCM(X, Y2) = B

X与Y2的最小公倍数为B
• Fuction GCD(X, Y ) means greatest common divisor between X and Y .

。。。不作赘述
• Fuction LCM(X, Y ) means lowest common multiple between X and Y .
BrotherJ loves Silchen so much that he wants Silchen to win the game. Now he wants to calculate how many number of X he can give to Silchen.

J大哥爱Silchen爱的不行以至于他想让Silchen获胜。现在他想计算一下他有多少个X可以给Silchen。

输入

Input is given from Standard Input in the following format:

Emergency Caution：多组输入，否则直接wrong answer到死（血的教训）。
A B
Constraints
1 ≤ A, B ≤ 1018
Both A and B are integers.

输出

Print one integer denotes the number of X.

input

3 12

output

根据题给条件，我们可以得出以下两条结论

已知Y1，Y2是任意取的，那也就是说X是不可能小于A（不能小于最大公约数），而且X也不可能大于B（不能大于最小公倍数），否则A，B都会改变。
如果X的最小公倍数是B，那么B一定是X的整数倍，即X必定为B的因子，又因为X不可能与A互质，那么X必为A的整数倍。

既然X同时为A的倍数，B的因子，那么就可以得到B为A的倍数，A为B的因子。

那么这个问题就变成了：A有多少个X可以使之乘到B，即B除了A之外还有什么因子——>找出(B/A)的因子个数。

找因子个数我们自然知道算数基本定理（唯一分解定理）这一个好用的工具，但是如果在这里直接使用很可能会超时，不过它仍然扮演很重要的角色。我们用一个特殊的质数判定算法——Miller Rabin算法来优化程序。

这个算法基于费马小定理的逆定理：

$\large a^{p-1}\equiv 1(mod p)$

若a为质数，且满足上式关系，则P大概率为质数。为什么说大概率呢？因为有一些反例，所以这个定理并不是一般关系，只是大部分的“特殊关系”...为了尽可能提高准确性，我们采用二次检测来进一步验证。

二次检测基于以下的一个式子：

$\large x^2\equiv 1(mod p)$

若p为质数，可以解得 $\large x=-1 || x=1$ ，在mod p的情况下就是 $\large x=p-1||x=1$

在验证过程中我们令 $\large a^{p-1}$ 为 $\large x^2$ 进行二次检测，首先可以确定是p-1一定是个偶数或者2，不然的话p+1一定就是个偶数，大于2的偶数是不可能为质数的。那么可以把p-1变成 $\large k*2^n$ 的形式（n和k为整数）。

最终二次检测式变为：

$\large a^{k*2^{n-1}}*a^{k*2^{n-1}}\equiv 1(mod p)$ ，通过同底数幂的乘法我们很容易就知道这是费马小定理的原式变形得到的。

在代码中我们怎么操作这样的检测呢？

*~*第一大块*~*
对于一个待测数p，首先把p-1里边所有的2提取出来，并记录提取了多少个（牢记 $\large k*2^n$ 这个形式），最后留存的就是k。
*~*第二大块*~*
通过一个循环，每一次更换一个底数a（可以从数组里预存几个数字，也可以直接rand（）几个），并对它进行k次的幂运算（mod p）,这个过程相当于同余方程模拟。由此得到同余方程右侧数字res
接着就是把k还原为 $\large k*2^n$ 的过程了。顺便，如果res=1的话直接符合条件，跳过这一步。由于还原的过程相当于两边同乘2，那么为了等价，对res进行快速乘（这是为了防止爆long long）的操作，期间如果得到res=1或者res=p-1的话跳出循环。
进行如下判断：如果经过step(5)处理后仍然得不到res=p-1，那么这个数的质数可能性就可以直接否定了。（由于k在处理前一定是一个奇数,而且既然要处理，那么res肯定不是1，所以这就可以直接下结论）
如果5次检测完全符合，那么断定他是质数。

值得一说的是，这里的质数判定指的是判定算数基本定理筛剩下的那个数的性质。因为在一开始的时候我们预处理了1e7的质数，所以余下的那个数必然很大，这不是朴素质数判定法则就能解决的。

以这个思路实现所有功能，我们需要5个模块：

质数筛
快速乘（基于乘法分配律）
快速幂
Miller rabin质数判定
算数基本定理

详情见代码。

#include <iostream>
#include <cstdio>
#include <bits/stdc++.h>
#include <algorithm>
#include <iomanip>
#include <cstring>
#include <cmath>
#define DETERMINATION main
#define lldin(a) scanf("%lld",&a)
#define din(a) scanf("%d",&a)
#define printlnlld(a) printf("%lld\n",a)
#define printlnd(a) printf("%d\n",a)
#define printlld(a) printf("%lld",a)
#define printd(a) printf("%d",a)
#define reset(a,b) memset(a,b,sizeof(a))
const long long int INF=0x3f3f3f3f;
using namespace std;
const double PI=acos(-1);
typedef long long ll;
typedef unsigned long long ull;
typedef long double ld;
const int mod=1000000007;
const int tool_const=1999112620000907;
const int tool_const2=33;
inline ll Unnamed_Scanner()
{
    ll tmp=0,si=1;
    char c=getchar();
    while(c>'9'||c<'0')
    {
        if(c=='-')
            si=-1;
        c=getchar();
    }
    while(c>='0'&&c<='9')
    {
        tmp=tmp*10+c-'0';
        c=getchar();
    }
    return si*tmp;
}
///Schlacht von Stalingrad
/**Although there will be many obstructs ahead,
the desire for victory still fills you with determination..**/
const int limit=1e7+200;
bool vis[limit];
ll TmpPrimes[9]= {2,3,5,211,317};//预处理的几个质数
ll primes[limit];
ll cnt=0;
void filter()//质数筛
{
    for(int i=2; i<=limit; i++)
    {
        if(vis[i]==0)
            primes[cnt++]=i;
        for(int j=0; j<cnt&&i*primes[j]<limit; j++)
        {
            vis[primes[j]*i]=true;
            if(i%primes[j]==0)
                break;
        }
    }
}
ll quick_multiply(ll a,ll b,ll mod)//快速乘
{
    ll ans=0;
    while(b)
    {
        if(b&1)
            ans=(a+ans)%mod;
        b>>=1;
        a=(a+a)%mod;
    }
    return ans;
}
ll quick_power(ll a,ll b,ll mod)//快速幂
{
    ll ans=1;
    while(b)
    {
        if(b&1)
            ans=quick_multiply(a,ans,mod);
        b>>=1;
        a=quick_multiply(a,a,mod);
    }
    return ans%mod;
}
bool miller_rabin(ll num)//质数判定
{
    if(num<2)
        return false;
    else if(num==2)
        return true;
    else if(num%2==0)
        return false;
    ll tmp2=num-1;
    while(!(tmp2&1))//拆除2的幂次
    {
        tmp2>>=1;
    }
    for(int i=0;i<5;i++)
    {
        if(num==TmpPrimes[i])
            return true;
        ll tmp=tmp2;
        ll mesa=quick_power(TmpPrimes[i],tmp2,num);
        while(tmp!=num-1&&mesa!=num-1&&mesa!=1)//恢复2的幂次
        {
            mesa=quick_multiply(mesa,mesa,num);
            tmp<<=1;
        }
        if(mesa!=num-1&&!(tmp&1))//如果就算这样都没找着num-1这样的结果
            return false;
    }
    return true;
}
ll judge(ll uxx)//算数基本定理
{
    ll totality=0;
    ll sum[5000]= {0};
    //cout<<miller_rabin(uxx)<<endl;
    for(int i=0; i<cnt&&primes[i]*primes[i]<=uxx; i++)//筛质因子
    {
        if(uxx%primes[i]==0)
        {
            totality++;
            while(uxx%primes[i]==0)
            {
                uxx/=primes[i];
                sum[totality]++;
            }
      
        }
        if(uxx==1)
            break;
    }
    //cout<<uxx<<endl;
    ll tmpres=1;
    //cout<<sum[1]<<endl;
    for(int i=1; i<=totality; i++)//计算因子数
        tmpres*=(sum[i]+1);
    //cout<<uxx<<endl;
    if(uxx>1)
    {
        if(miller_rabin(uxx)==false)
        {
            ll tmp=sqrt(uxx);
            if(tmp*tmp==uxx)//如果这是个完全平方数，那么肯定由一个数的平方构成。
                tmpres*=3;
            else
                tmpres*=4;//如果这是个普通合数，那么必定存在两个不同的质因子
        }
        else
            tmpres*=2;//质数本身自己就是一个质因子
    }
    return tmpres;
}
int DETERMINATION()
{
    filter();
    ll a,b;
    while(cin>>a>>b)
    {
        if(b%a!=0)
    {
        cout<<0<<endl;
    }
    else
    {
        ll tmp=b/a;
        cout<<judge(tmp)<<endl;
    }
    }
    return 0;
}

FZU code carnival-GCDLCM (最大公约数与最小公倍数）（Miller Rabin质数判定）

猜你喜欢