[HDU 4953 Periodic Binary String] Euclid的应用

来源:互联网 发布:php销售管理系统 编辑:程序博客网 时间:2024/06/05 06:23

题目

http://acm.hdu.edu.cn/showproblem.php?pid=4953

分析

首先,由于循环,这题的l,r  其实可等价于1,(r-l)+1。设n=(r-l+1) 然后设循环外多出来的一部分值为x, 另一部分为y。

原问题等价为   ax+by=c(mod p)  y=a‘y+b’(mod p)    0<=x<2^l    0<=y<2^r  

[ay+b%p<k]=[(ay+b)/p]-[(ay+b-k)/p]                 则答案为  sigma  [(ay+b)/p]-[(ay+b-k)/p] 

运用Euclid  求[(ay+b)/p]      sigma(0-n) [(ay+b)/p]   =   sigma (1-[(an+b)/p]) n- [(py-b-1)/a]

代码

#include<stdio.h>#include<stdlib.h>#include<string.h>#include<iostream>#include<algorithm>#define INF 2147483647#define S(a) scanf("%d",&a)#define cl(a,b) memset(a,b,sizeof(a))#define foru(i,a,b) for(int i=a;i<=b;++i)#define ford(i,a,b) for(int i=a;i>=b;--i)#define fore(i,a) for(int i=g[a];i;i=next[i])using namespace std;typedef long long int64;const int64 MOD=1000000007;int64 p,x,l,r,k,n,a,b,c,sa,sb,k1,k2,ans,ny;int test;int64 qadd(int64 a,int64 b,int64 MOD){    int64 res=0;    a=(a%MOD+MOD)%MOD; b=(b%MOD+MOD)%MOD;    for(;b;b>>=1,a=(a+a)%MOD)      if(b&1) res=(res+a)%MOD;    return (res+MOD)%MOD;}int64 qpow(int64 a,int64 b,int64 MOD){    int64 res=1;    for(;b;b>>=1,a=qadd(a,a,MOD))      if(b&1) res=qadd(res,a,MOD);    return (res+MOD)%MOD;}int64 div(int64 &a,int64 b,int64 c,int64 d,int64 MOD){    int64 res=(d/c)%MOD,rst=d%c,p=(a/c)%MOD,pr=a%c;    for(;b;b>>=1,p=(p*2+pr*2/c)%MOD,pr=pr*2%c)      if(b&1) res=(res+p+(rst+pr)/c)%MOD,rst=(rst+pr)%c;    a=rst;    return (res+MOD)%MOD;}int64 qdiv(int64 a,int64 b,int64 c,int64 MOD)//???{    int64 res=0,rst=1,p=a/c,pr=a%c;    for(;b;b>>=1,p=(qadd(qadd(p,p,MOD),c,MOD)+qadd(p,pr,MOD)*2)%MOD,p=(p+div(pr,pr,c,0,MOD))%MOD)      if(b&1) res=(qadd(qadd(res,p,MOD),c,MOD)+qadd(res,pr,MOD)+qadd(p,rst,MOD))%MOD,res=(res+div(rst,pr,c,0,MOD))%MOD;    return (res+MOD)%MOD;}int64 work2(int64 a,int64 b,int64 n,int64 MOD){    int64 res;    if(n&1) res=qadd(n,(n+1)/2,MOD);    else res=qadd(n/2,n+1,MOD);    res=(qadd(a,res,MOD)+qadd(b,n+1,MOD))%MOD;    return (res+MOD)%MOD;}int64 work(int64 a,int64 b,int64 c,int64 n,int64 MOD){    int64 aa=a/c,bb=(b%c+c)%c;    int64 res=work2(aa,(b-bb)/c,n,MOD);    a%=c; if(!a) return res;    int64 t=a,m=div(t,n,c,bb,p);    res=(res+qadd(n,m,MOD)-((bb+1)/a+((bb+1)%a>0))-work(c,-bb-1,a,m,MOD)+2*MOD)%MOD;    return res;}int main(){    freopen("input.txt","r",stdin);freopen("output.txt","w",stdout);    while(scanf("%I64d%I64d%I64d%I64d%I64d",&p,&x,&l,&r,&k)!=EOF)    {        if(p+x+l+r+k==0) break;        n=(r-l+1); l=n%k; r=k-l; n=n/k; ans=0;        ny=(qpow(2,k,p)+p-1)%p; c=x%p;        if(ny)        {            a=(qpow(2,qadd(k,n+1,p-1),p)+p-1)%p; b=(qpow(2,qadd(k,n,p-1),p)+p-1)%p; b=qadd(b,qpow(2,l,p),p);            ny=qpow(ny,p-2,p); b=qadd(ny,b,p); a=qadd(ny,a,p);        }        else        {            a=(n+1)%p; b=qadd(n,qpow(2,l,p),p);        }        sa=qdiv(2,l,p,MOD); sb=qdiv(2,r,p,MOD);        k1=qpow(2,l,p); k2=qpow(2,r,p);        if(!a)        {            if(!b)            {                if(!c) printf("Case #%d: %I64d\n",++test,qpow(2,l+r,MOD));                else printf("Case #%d: 0\n",++test);                continue;            }            c=qadd(c,qpow(b,p-2,p),p); ans=qadd(qpow(2,l,MOD),sb+(c<k2),MOD); //????            printf("Case #%d: %I64d\n",++test,ans);            continue;        }        ny=qpow(a,p-2,p); b=(p-qadd(b,ny,p))%p; c=qadd(x,ny,p);        if(!b)        {            ans=qpow(2,r,MOD); ans=qadd(ans,sa+(k1>c),MOD);            printf("Case #%d: %I64d\n",++test,ans);            continue;        }        ans=qadd(qadd(sa,sb,MOD),p,MOD);        ans=(ans+qadd(k1,sb,MOD)+qadd(k2,sa,MOD))%MOD;        ans=(ans+work(b,c,p,k2-1,MOD)-work(b,c-k1,p,k2-1,MOD)+MOD)%MOD;        printf("Case #%d: %I64d\n",++test,ans);    }    return 0;}


0 0
原创粉丝点击