2013年11月18日 星期一

STEP5 0035 蘿莉不同色問題

[數論]
先備:1. 指數原根 2. Baby step-Giant step

這題真不是普通的好玩XDD

總之,照題意,要解x:
先抓p的原根g
然後對兩邊取index:
然後解Linear方程。


/*      LANG:C++
**  http://web2.ck.tp.edu.tw/~step5/probdisp.php?pid=0035
**      NUMBER THEORY
**      2013-11-8
*/
#include<cstring>
#include<stdio.h>
#include<stdlib.h>
#include<algorithm>
#include<math.h>
#include<assert.h>
#include<vector>
#define PMAX 100000
#define Mid(a,b) ((a+b)>>1)
#define LL long long int
#define PS printf("PS[%d]\n",step++);
using namespace std;
int step=0;
///Set up Prime Table
vector<LL> PRIME;
bool seive[PMAX+1];

void BuildPrimeTable()
{
    memset(seive,true,sizeof(seive));
    seive[1]=false;
    PRIME.clear();
    for(LL lx=2;lx<=PMAX;lx++)
    {
        if(seive[lx])
        {
            PRIME.push_back(lx);
            for(LL ly=2;ly*lx<=PMAX;ly++)
                seive[lx*ly]=false;
        }
    }
}

///Get the primitive root of P
LL P;
vector<LL> Drc;
LL _Pow(LL a,LL n,LL _PP)
{
    if(n==0) return 1;
    if(n==1) return a%_PP;
    LL k=_Pow(a%_PP,n/2,_PP);
    if(n%2)
        return (((k*k)%_PP)*a)%_PP;
    else
        return (k*k)%_PP;
}
LL _Pow(LL a,LL n)
{
    return _Pow(a,n,P);
}
//536877277
void DrcPrime()
{

    LL prc=P-1;
    Drc.clear();
    //printf("a1");
    for(int lx=0;(PRIME.size()>lx)&&(PRIME[lx]<prc);lx++)
    {
        if(prc%PRIME[lx]==0)
        {
            Drc.push_back(PRIME[lx]);
            while(prc%PRIME[lx]==0)
                prc/=PRIME[lx];
            //printf("%I64d\n",PRIME[lx]);
        }

    }

    if(prc>1)
        Drc.push_back(prc);
    //printf("a1");
}
LL GetPrimitiveRoot()
{
    for(LL lx=2;lx<=P;lx++)
    {
        bool isRoot=true;
        for(LL ly=0;ly<Drc.size();ly++)
            if(_Pow(lx,(P-1)/Drc[ly])==1){isRoot=false;break;}
        if(isRoot)
            return lx;
    }
}
///Babystep  prt^a=r
class Data
{
public:
    LL index;
    LL cont;
    Data(){}
    void set(const LL _i,const LL _c){index=_i,cont=_c;}
};
bool operator<(Data _aa,Data _bb){return (_aa.cont<_bb.cont);}
Data* brysch(Data *I,LL c,LL m)
{
    //prLLf("SCH:%d\n",c);
    LL start=0,end=m-1;
    while((end-start)>1)
    {
        //prLLf("END=%d START=%d\n",end,start);
        if(I[start].cont==c)
            return &I[start];
        else if(I[end].cont==c)
            return &I[end];
        else if(I[Mid(start,end)].cont==c)
            return &I[Mid(start,end)];
        if(c>I[Mid(start,end)].cont)
            start=Mid(start,end);
        else
            end=Mid(start,end);
    }
    if(I[start].cont==c)
        return &I[start];
    else if(I[end].cont==c)
        return &I[end];
    else
        return NULL;
}
Data sqrtb[PMAX];
LL GetInd(LL prt,LL r)
{

    r=(r%P+P)%P;//prLLf("GetInd(%d,%d)\n",prt,r);
    LL m=ceil(sqrt(P));

    for(LL lx=0,prc=1;lx<m;lx++,prc=(prc*prt)%P)
        sqrtb[lx].set(lx,prc);
    sort(sqrtb,sqrtb+m);
    /*for(LL lx=0;lx<m;lx++)
        prLLf("%d:%d\n",sqrtb[lx].index,sqrtb[lx].cont);*/
    LL base=_Pow(prt,m*(P-2));
    //printf("base=%d\n",base);
    for(LL lx=0,prc=1;lx<m;lx++,prc=(prc*base)%P)
    {
        //printf("prc=%d\n",prc);
        Data *RES=brysch(sqrtb,(r*prc)%P,m);
        if(RES==NULL) continue;
        return (*RES).index+m*lx;
    }
    return -1;
}

///ax=b mod P-1 ,return x
LL phi(LL m)
{
    LL prc=P;
    LL ret=1;
    for(LL lx=0,prc=m;(lx<PRIME.size())&&(PRIME[lx]<=prc);lx++)
    {
        if(prc%PRIME[lx]==0)
        {
            LL nn=1;prc/=PRIME[lx];
            while(prc%PRIME[lx]==0)
                prc/=PRIME[lx],nn*=PRIME[lx];
            ret*=nn*(PRIME[lx]-1);
        }
    }
    return ret;
}
LL Linear(LL a,LL b)
{
    LL MM=P-1;
    LL gg=__gcd(a,MM);
    //printf("gg=%I64d\n",gg);
    //assert(gg);
    MM/=gg;
    a/=gg;
    b/=gg;
    //printf("{%I64d}}\n",phi(MM));
    return (b*_Pow(a,phi(MM)-1,MM))%MM;
}
int main()
{
    BuildPrimeTable();
    LL X0,X1,K,N;
    LL proot;
    long long int Kqq;
    //scanf("%d",&P);
    scanf("%I64d %I64d %I64d %I64d %I64d",&X0,&X1,&Kqq,&N,&P);
    //PS
    if(Kqq>=0)
        K=Kqq;
    else
        K=P+Kqq;
    DrcPrime();
    proot=GetPrimitiveRoot();
    LL aa,bb;//ax=b;
    LL qq=0;
    if(X1>=X0)
        qq=(X1-X0);
    else
        qq=P+X1-X0;
    aa=N-1;bb=GetInd(proot,(K*_Pow(qq,P-2))%P);
    //printf("a=%d b=%d\n",aa,bb);
    //printf("proot=%I64d\n",proot);
    //PS
    LL a,b;a=_Pow(proot,Linear(aa,bb));
    if(X1>=X0*a)
        qq=X1-X0*a;
    else
        qq=P+X1-X0*a;
    b=(qq%P+P)%P;
    //PS
    printf("%I64d %I64d\n",a,b);
    for(LL lx=0,prc=X0;lx<=N;lx++,prc=(prc*a+b)%P)
        printf("%I64d ",prc);
    printf("\n");
    return 0;
}

沒有留言:

張貼留言