第一章 数论(基础知识 1.9)

1.9Baby-Step-Giant-Step及扩展算法

BSGS算法用于解决形如A^xB(mod C)的问题

以下分析完全复制此博客：http://www.cnblogs.com/hchlqlz-oj-mrj/p/4769180.html

                                                https://blog.csdn.net/shengtao96/article/details/51362592

学这个知识需要具备以下知识：快速幂取模，扩展欧几里得，同余知识，哈希表

一，普通的Baby-Step-Giant-Step

对于普通的BSGS，使用时有个限制条件就是：C只能是一个素数。；当C不是素数的时候，便要用到扩展BSGS，先来看普通的：

做法：

step1:令   x=a*m+b;

                          m=ceil(sqrt(C));(对上取整)             （0<=a<m,0<=b<m;）？？？？？？？？？？？？？？？？？？？？？？？

           所求等式等价于：  ( ( A ^ a ) ^ m ) * A ^ b ≡ B (mod C )

step2:枚举a,从0到m,对每个A算出A^a,放到Hash表中，如h=A^a,则mp[h]=a;

step3:用D来表示((A^a)^m),用Y来表示A^b，原式便可化为D*YB(mod C);

step4:再次枚举a从0到m，可以逐个算出D的值，在有了D的值的基础上，运用扩展欧几里得算法可以算出Y，现在查找原本记录的              表，查看这个Y值是否存在，如果存在可返回a^m+b;如果找完0到m，依然没有发现可行的解，表示无解，退出。

 

先把代码放在这吧，好难~   :(

#include<stdio.h>
#include<math.h>
#include<map>
using namespace std;
/*
  Baby-Step-Giant-Step:   用于解决形如: A^x = B (mod C) 的问题
  <><> 普通的解法只能解决C是素数的情况
  首先，令 x= a*m + b ; -- m= Ceil(sqrt(c));  --  Ceil 表示上取整
  then 原式等价于 ( (A^a )^m )*(A^b)=B (mod c);
  then 那么我们可以枚举 a: 0-m ,记录到 (a,A^a) ，用H 来表示: (A^a)^m ;
  那么 原式等价于 H * A^b = B (mod c) ;
  then 求出 A^m ; 枚举 a: 0-m ，对于每个值判断有没有 A^b 存在，如果有便退出...
*/
// 这个代码不敢保证正确无误，因为没有找到题来A，不过思路应该是没错
//了解了普通的BSGS，可以去看下面那份代码，注释得比较详细...
#define EPS 0.0000001
typedef long long LL;
LL pow(LL a,LL b,LL c)    //快速幂取模
{
  LL ans=1;
  while(b>0)
  {
    if(b&1)
    {
      ans=ans*a%c;
    }
    a=a*a%c;
    b>>=1;
  }
  return ans;
}
void exGcd(LL a,LL b,LL &d,LL &x,LL &y)
{
  if(b==0)
  {
    d=a;
    x=1;
    y=0;
    return ;
  }
  exGcd(b,a%b,d,x,y);
  LL t=x;
  x=y;
  y=t-a/b*y;
}
LL BSGS(LL A,LL B,LL C)
{
  map<LL,int> mp;
  // LL m=(LL)ceil(sqrt(C));  据英明神武的某朱说Gcc不支持这个ceil
  LL m=(LL)(sqrt(C)+1-EPS);
  for(int i=m-1;i>=0;i--)
  {
    LL h=pow(A,i,C);
    mp[h]=i;
  }
  LL h=pow(A,m,C);
  LL d,x,y;
  for(int i=0;i<m;i++)
  {
    LL k=pow(h,i,C);
    exGcd(k,C,d,x,y);
    if(B%d!=0) continue;
    x=((x*B/d%C)+C)%C;
    if(mp.find(x)!=mp.end()) return mp[x]+i*m;
  }
  return -1;
}
int main()
{
  LL A,B,C;
  while(scanf("%I64d%I64d%I64d",&A,&B,&C)!=EOF)
  {
    printf("%I64d\n",BSGS(A,B,C));
  }
  return 0;
}

二，扩展Baby-Step-Giant-Step

对于C不是素数的情况，不能照搬上面的做法，但也大同小异。

给出3个数，A,C,B表示A^xB(mod C),要求求最小的x：

#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#include<conio.h>
#define Mod 40007       //取模的大小，哈希表的大小...
#define Max 100000     //存放的总数
#define EPS 0.000000001//精度
typedef long long LL;
class Hash             //手写哈希
{
  public:
    int hs[Mod];       //哈希值  设定的哈希函数为 原值 % Mod ，所以哈希值有可能是 0 ~ Mod-1
    int next[Max];     //链表    存放哈希值相等的一条链，他的大小取决于所有原值的数量
    LL S[Max];         //存放原值
    int H[Max];        //存放所有哈希值
    int sn;            //不同原值的数量
    int hn;            //不同哈希值的数量
    Hash()             //构造函数: 定义Hash类变量时初始化
    {
      sn=0;
      hn=0;
      for(int i=0;i<Mod;i++)
        hs[i]=0;
    }
    void clear()       //清空函数
    {
      sn=0;
      for(int i=0;i<hn;i++)
        hs[H[i]]=0;
      hn=0;
    }
    void add(LL s)           //加入
    {
      int ha=abs(s)%Mod;     //计算哈希值
      if(hs[ha]==0)          //如果该哈希值还未出现过
      {
        H[hn++]=ha;          //将该哈希值记录起来，同时哈希值数量加 1
      }
      sn++;                  //0 表示结尾，所以从1 开始存，原值数量加 1，特别针对 hs数组
      S[sn]=s;               //将原值记录起来
      next[sn]=hs[ha];       //原本原值记录位置
      hs[ha]=sn;             //最新原值记录位置，如果从0 开始存，就无法判断此时是空还是1个值
      //比如：5 和 10 有一样的哈希值 ，并且 5 和 10 先后加入 那么有：
      //5 加入: next[1] = 0; hs[5] = 1; hs[5] 是哈希值为5 的头，表示第一个原值在1的位置
      //10加入: next[2] = 1; hs[5] = 2; 表示第一个哈希值为5的在2，第二个在1，第三个不存在
    }
    int find(LL s)           //查找
    {
      int ha=abs(s)%Mod;     //计算哈希值
      int k=hs[ha];          //头
      while(k!=0)
      {
        if(S[k]==s) return k;//找到
        k=next[k];           //下一个节点
      }
      return 0;              //表示没找到
    }
};
int pow(int a,int b,int c)   //快速幂取模
{
  LL ans=1;
  a%=c;
  while(b>0)
  {
    if(b&1)
    {
      ans=(LL)ans*a%c;
    }
    b>>=1;
    a=(LL)a*a%c;
  }
  return ans;
}
void exGcd(int a,int b,int &d,int &x,int &y)     //扩展欧几里得
{
  if(b==0)
  {
    d=a;
    x=1;
    y=0;
    return ;
  }
  exGcd(b,a%b,d,x,y);
  LL t=x;
  x=y;
  y=t-a/b*y;
}
int Gcd(int a,int b)            //最大公约数
{
  return b?Gcd(b,a%b):a;
}
/*
  对于 A^x = B (mod C) 这个式子，如果无法保证C是素数，就得使用扩展Baby-Step-Giant-Step
  首先: 有 A^x = C * k + B ; 然后，我们一步步地消去 A 和 C 的公因子
  令 D 从 1 开始， 原式就可以化为 D*A^x = C*k+B;
                              --> D*A/(A,C)*A^(x-1) = C/(A,C)*k +B/(A,C);
                              --> ...
                              --> D*(A/(A,C))^(co))*A^(x-co) = C/((A,C)^co)*k + B/((A,C)^co)
  上面的那个循环知道 A、C互质时结束，没进行一次，D便等于 D*A/(A,C);
  最后，原式化为 D * A^(x-co) = B (mod C) 如果上面循环中，出现B无法整除的情况就表示无解退出
  这样，我们便保证了 D 和 A 都和 C 互质，这时候，便可以利用普通的Baby-Step-Giant-Step解决
  同时，要小心的一点就是，因为这样解出来的最小解是肯定大于co的，但如果真实解小于co，就会出
  错，所以我们可以事先枚举0 - p 的i，进行特判，这里的p=log(C)，因为每次消因子最少也要消去2，
  所以最多消log(C)次，co最大也就是这个了。
*/
int BSGS(int A,int B,int C)     //扩展 Baby-Step-Giant-Step
{
  Hash hp;
  for(int i=0;i<50;i++)         //枚举1 - log(C) ,可以稍大一点...
  {
    if(pow(A,i,C)==B) return i;
  }
  int tmp,co=0;
  LL D=1,k=1;
  while((tmp=Gcd(A,C))!=1)      //消因子循环
  {
    co++;
    if(B%tmp!=0) return -1;
    C/=tmp;
    B/=tmp;
    D=D*A/tmp%C;
  }
  hp.clear();
  int m=(int)(sqrt((double)C)-EPS);   //m=ceil(sqrt(C)),这里的ceil表示上取整
  for(int i=0;i<m;i++,k=k*A%C)
    if(hp.find(k)==0) hp.add(k);
  int d,x,y;
  for(int i=0;i<m;i++)
  {
    exGcd(D,C,d,x,y);
    x=((LL)x*B%C+C)%C;                //不用进行B是否能整除d的判断，因为D和C永远互质
    int w=hp.find(x);
    if(w!=0) return w-1+co+i*m;
    D=D*k%C;                          //D最开始的时候和C互质，同时k也和C互质
  }
  return -1;                          //找不到
}
int main()
{
  int A,B,C;
  while(scanf("%d%d%d",&A,&C,&B)!=EOF)
  {
    if(B>=C)
    {
      printf("Orz,I can’t find D!\n");
      continue;
    }
    int t=BSGS(A,B,C);
    if(t==-1)
      printf("Orz,I can’t find D!\n");
    else printf("%d\n",t);
  }
  return 0;
}

恐惧的一批~~~