這個是求擴展離散對數問題。XY mod Z = K,給出X,Z,K,求Y。
當Z是素數的時候直接用baby-step算法即可了。但是,模數不是素數的情況怎麼辦了。
方程a^X = b % c,可以進行一系列的轉化。假設d = gcd(a,c),由a^(x-1)*a = b % c,知道a^(x-1)要存在必須滿足
gcd(a,c)|b,如果滿足這個條件,那麼我們可以在方程2邊同時除以d,方程是不變的。因為a^x = b + k * c,再除以公約數
d,得到方程a^(x-1)*a/d = b / d + k * c / d。根據以上推論,我們可以不斷的除以d,直到gcd(a,c)=1。
假設我們除了k次,那麼方程轉化為a^(x-k) * a^k/d^k = b / d^k + k * c / d^k。令d = a^k/d^k,b' = b / d^k,
c' = c / d^k,x' = x - k,方程轉化為a^x' * d = b' % c',得到a^x' = b' * d^-1 % c'。
現在直接用baby-step解方程a^x' = b' * (d^-1) % c'即可。注意到x=x'+k,如果存在x小於k的解,那麼x'小於0,但是
由baby-step是不會求負的次數的,所以需要先枚舉一下是否存在小於k的解,由於輸入的數據不會超過10^9的,假設k不超過50
進行枚舉即可了。
代碼如下:
#include <stdio.h>
#include <math.h>
#include <algorithm>
using namespace std;
typedef long long INT;
#define MAX (1000000)
INT nData[MAX];
INT nKey[MAX];
INT HashPos(INT key)
{
return ((unsigned)(key ^ 0xA5A5A5A5)) % MAX;
}
void HashAdd(INT key, INT data)
{
INT nPos = HashPos(key);
while (nData[nPos] != -1)
{
nPos = (nPos + 1) % MAX;
}
nData[nPos] = data;
nKey[nPos] = key;
}
INT HashQuery(INT key)
{
INT nPos = HashPos(key);
while (nData[nPos] != -1)
{
if (nKey[nPos] == key)
{
return nData[nPos];
}
nPos = (nPos + 1) % MAX;
}
return -1;
}
INT MultMod(INT nA, INT nB, INT nC)
{
INT nAns = 0;
while (nB)
{
if (nB & 1)
{
nAns = (nAns + nA) % nC;
}
nA = (2 * nA) % nC;
nB >>= 1;
}
return nAns;
}
INT PowerMod(INT nA, INT nX, INT nC)
{
INT nAns = 1;
nA %= nC;
while (nX)
{
if (nX & 1)
{
nAns = MultMod(nAns, nA, nC);
}
nA = MultMod(nA, nA, nC);
nX >>= 1;
}
return nAns;
}
INT gcd(INT nA, INT nB)
{
if (nA < nB)swap(nA, nB);
while (nB)
{
INT nT = nA;
nA = nB;
nB = nT % nB;
}
return nA;
}
//d = nA * nX + nB * nY(nA > nB, nA是模數)
INT egcd(INT nA, INT nB, INT& nX, INT& nY)
{
if (nA < nB)swap(nA, nB);
if (nB == 0)
{
nX = 1;
nY = 0;
return nA;
}
INT nRet = egcd(nB, nA % nB, nX, nY);
INT nT = nX;
nX = nY;
nY = nT - (nA / nB) * nY;
return nRet;
}
INT GetAns(INT nA, INT nB, INT nC)
{
if (nC == 0)return -1;
//先枚舉0-50,擴展baby-step的過程可能會漏掉這些解
INT nTemp = 1;
nB %= nC;
for (INT i = 0; i <= 50; ++i)
{
if (nTemp == nB)
{
return i;
}
nTemp = MultMod(nTemp, nA, nC);
}
//如果nC不是素數,那麼方程nA^x = nB + k*nC
//可以不到除以gcd(nC,nA)
//如果gcd(nC,nA)|nB不成立,方程無解,
//這個由a*x=b%c有解必須滿足gcd(a,c)|b一樣
INT d;
INT nD = 1;//nD最後是A^k次,k是nC中因子d的次數
INT k = 0;
while ((d = gcd(nC, nA)) != 1)
{
k++;
nC /= d;
if (nB % d)return -1;
nB /= d;
nD = MultMod(nD, nA / d, nC);
}
//現在方程轉化為nA^(x-k) * nA^k/d^k = nB/d^k % nC/d^k
//其實就是方程2側除以d^k次而已,這樣的做法與原方程是等價的
//令nD = nA^k/d^k,則nA^x'*nD = nB' % nC',
//解該方程,那麼x=x'+k
//注意,如果x<k,那麼x'為負數,baby-step無法求出,故在函數開頭進行枚舉
memset(nKey, -1, sizeof(nKey));
memset(nData, -1, sizeof(nData));
INT nM = ceil(sqrt(1.0 * nC));
nTemp = 1;
for (INT j = 0; j <= nM; ++j)
{
HashAdd(nTemp, j);
nTemp = MultMod(nTemp, nA, nC);
}
INT nK = PowerMod(nA, nM, nC);
for (int i = 0; i <= nM; ++i)
{
INT x, y;
egcd(nC, nD, x, y);//y = nD^-1,nD = nD*(nA^m)^i
y = (y + nC) % nC;//這句話是必須的,y很可能就是負數
INT nR = MultMod(y, nB, nC);//nR=nB*nD^-1
int j = HashQuery(nR);
if (j != -1)
{
return nM * i + j + k;
}
nD = MultMod(nD, nK, nC);
}
return -1;
}
int main()
{
INT nA, nB, nC;
while (scanf("%I64d%I64d%I64d", &nA, &nC, &nB), nA + nB + nC)
{
INT nAns = GetAns(nA, nB, nC);
if (nAns == -1)
{
printf("No Solution\n");
}
else
{
printf("%I64d\n", nAns);
}
}
return 0;
}