CASE (1sim 3)
(n)组测试数据,每次输入一个数(x),求(19^x)。
测试点(1),(x=0,1,dots n-1),可以直接递推。
测试点(2)要开long long
并用快速幂。
测试点(3)中(x)超出了long long
范围。根据欧拉定理,当(a,p)互质时,(a^bequiv a^{bmodvarphi(p)}pmod p)。因为模数(p=998244353)是一个质数,所以(varphi(p)=p-1=998244352)。将数字一位一位读进来然后对(998244352)取模即可。
代码片段:
typedef long long ll;
ll pow_mod(...){
//快速幂略
}
ll read_mod(){
ll x=0;
char ch=getchar();
while(!isdigit(ch))ch=getchar();
while(isdigit(ch))x=(x*10+ch-'0')%Math::PHI_MOD,ch=getchar();
return x;
}
void solve1(ll x=998244353){
Init(x);//MOD=x
int n=read();while(n--){
ll x=read_mod();
printf("%lld
",pow_mod(19,x));
}
}
CASE (4sim 5)
测试点(4),功能编号是1?
,容易猜到意思是说模数不确定。于是考虑枚举模数。先找出输出文件中最大的数(mx),则模数一定大于(mx)。我们从(mx+1)开始枚举模数。容易找到模数为(1145141)。
找模数代码:
注意这里我们只判断了模数是质数的情况,因为我们相信出题人足够良心。
ll read_mod(string s){
ll x=0;
for(int i=0;i<(int)s.length();++i)x=(x*10+s[i]-'0')%Math::PHI_MOD;
return x;
}
bool test_mod(ll m,vector<pair<string,ll> >&vec){
if(!is_prime(m))return false;
Init(m);
for(int i=0;i<(int)vec.size();++i){
ll x=read_mod(vec[i].fst);
if(pow_mod(19,x)!=vec[i].scd)return false;
}
return true;
}
void get_mod(){
freopen("software4.in","r",stdin);
int n;
string str;cin>>str>>n;
for(int i=1;i<=3;++i)cin>>str;
vector<pair<string,ll> >vec;
for(int i=0;i<5;++i){
cin>>str;
vec.pb(mk(str,233));
}
freopen("software4.ans","r",stdin);
for(int i=1;i<=3;++i)cin>>str;
ll mx=0;
for(int i=0;i<n;++i){
ll x;cin>>x;
if(i<(int)vec.size())vec[i].scd=x;
mx=max(mx,x);
}
// for(int i=0;i<(int)vec.size();++i){
// cerr<<vec[i].fst<<" "<<vec[i].scd<<endl;
// }
cerr<<"max = "<<mx<<endl;
for(ll i=mx+1;i<=LLONG_MAX;++i){
if(test_mod(i,vec)){
cerr<<i<<endl;
return;
}
}
}
找到模数后:
void solve2(){
solve1(1145141);
}
测试点(5),发现答案都很大,说明模数也很大,用刚刚的方法大概是别想在考试结束前找到模数了。
考虑更巧妙的数论技巧。我们在输入文件中找两个接近的数(x,y),满足(x<y)且(ans_x>ans_y)。假设模数为(p),则可以知道(ans_xcdot19^{y-x}equiv ans_ypmod p)。
所以模数(p)就一定是(ans_xcdot 19^{y-x}-ans_y)的约数。注意这个数本身很大,要用__int128
存。暴力枚举其所有约数复杂度是(O(sqrt{n}))的,依然无法接受。于是做如下剪枝:
-
从输出文件中找出最大的数(mx),模数一定大于(mx)。
-
因为我们相信出题人足够良心,所以模数应该是质数。
-
因为我们相信出题人足够良心,所以模数应该在
long long
范围内。
判断大质数的方法会在CASE 8~10部分详细介绍。
用如下的程序可以很快地求出,模数是(5211600617818708273)。
pair<int,ll>data[10004];
void print128(__int128 x){
const ll base=1e18;
if(x<base)cout<<(ull)x<<endl;
else{
cout<<(ll)(x/base);
ll t=x%base;int len=0;
while(t)len++,t/=10;
for(int i=1;i<=18-len;++i)cout<<"0";
if(x%base)cout<<(ll)(x%base)<<endl;
}
}
void get_big_mod(){
freopen("software5.in","r",stdin);
string s;cin>>s;
int n;cin>>n;
for(int i=1;i<=n;++i)cin>>data[i].fst;
freopen("software5.ans","r",stdin);
ll mx=0;
for(int i=1;i<=n;++i){cin>>data[i].scd;mx=max(mx,data[i].scd);}
sort(data+1,data+n+1);
int x=0,y=0;
for(int i=1;i+1<=n;++i){
for(int j=i+1;j<=n;++j){
if(data[j].fst==data[i].fst)continue;
if(data[j].scd>data[i].scd)continue;
// i.f<j.f && i.s>j.s
if(!x||data[j].fst-data[i].fst<data[y].fst-data[x].fst)x=i,y=j;
}
}
cerr<<data[x].fst<<" "<<data[x].scd<<endl;
cerr<<data[y].fst<<" "<<data[y].scd<<endl;
cerr<<"--------------------"<<endl;
__int128 X=(__int128)data[x].scd*(__int128)19*(__int128)19-(__int128)data[y].scd;
print128(X);
cerr<<"--------------------"<<endl;
vector<pair<string,ll> >vec;
for(int i=n;i>=n-5+1;--i){
string str="";
ull x=data[i].fst;
while(x)str=(char)('0'+x%10)+str,x/=10;
vec.pb(mk(str,data[i].scd));
}
for(__int128 i=2;i*i<=X;++i){
if(X%i!=0)continue;
if(i>mx&&i<LLONG_MAX){
if(test_mod(i,vec)){
cerr<<(ll)i<<endl;
}
}
if(X/i>mx&&X/i<LLONG_MAX){
if(test_mod(X/i,vec)){
cerr<<(ll)(X/i)<<endl;
}
}
}
}
找到模数后:
void solve3(){
solve1(5211600617818708273LL);
}
CASE (6sim 7)
功能编号是wa
,注意到前几个数和第一个测试点相同,后面开始出现负数。结合题目下方关于自然溢出的提示,不难猜到这个两个测试点应该是(19^x)发生了自然溢出。
注意,自然溢出的过程非常奇特,不能用快速幂,否则就溢不出我们想要的效果。最保险的方法是按照题目的提示做加法。
对于测试点(7),输入的数值都非常大,直接做加法递推显然是不合适的。打表可以发现,(19^x)的自然溢出有循环节。这个循环节开始于(55245),长度是(45699)。
找循环节的代码:
//#include <windows.h>
map<int,int>mp;
void get_repetend(){
//找自然溢出循环节
Init(998244353);
int x=1;
for(int i=0;i<10000000;++i){
if(mp.count(x)==0)mp[x]=i;
else{
//cerr<<mp[x]<<" "<<i<<endl;//Sleep(100);
cerr<<"start: "<<mp[x]<<endl<<"term: "<<i-mp[x]<<endl;
return;
}
int new_x=x;
for(int j=1;j<=18;++j)new_x=((int)((uint)new_x+(uint)x));
x=new_x%Math::MOD;
}
}
找出循环节之后的实现:
const int ST=55245,LEN=45699;
int ans[ST+LEN];
void solve4(){
Init(998244353);//MOD=998244353
int x=1;
for(int i=0;i<ST+LEN;++i){
ans[i]=x;
int new_x=x;
for(int j=1;j<=18;++j)new_x=((int)((uint)new_x+(uint)x));
x=new_x%Math::MOD;
}
int n=read();
for(int i=1;i<=n;++i){
ll x=readll();
if(x<ST+LEN)printf("%d
",ans[x]);
else printf("%d
",ans[(x-ST)%LEN+ST]);
}
}
CASE (8sim 10)
首先观察到,输入的每一行,对应地输出了一个字符串。设输入的数为(x,y),则输出字符串的长度为(y-x+1)。不难猜到输出的是([x,y])这一区间内每个数的信息。由功能编号p
可以联想到prime number,质数。由此发现,输出字符串的每一位,表示([x,y])区间里的每个数是否是质数,若为质数则该位为p
,否则为.
。
测试点(8)可以使用线性筛法。
对于测试点(9,10),这里介绍一种快速的素数判定算法:米勒-罗宾法。
前文我们已经介绍了欧拉定理,对于它的一种特殊情况:模数为质数时,我们得到费马小定理:
若(p)为质数,(forall a<p),有:
[a^{p-1}equiv 1pmod p ]
我们选取一些(a),然后把(p)带入,进行判定。若对于其中某个(a_i),(p)不满足费马小定理,则可以确定(p)不是质数。
但这个条件是必要不充分条件。也就是说,存在一类伪素数,它们虽然不是素数,但满足费马小定理。
为避免这种情况,我们的算法要进行二次探测。它基于如下定理:
对一个素数(p),若有(x^2equiv 1pmod p),则(xequiv 1pmod p)或(xequiv p-1pmod p)。
证明:移项得(x^2-1equiv 0pmod p)。用平方差公式展开得((x+1)(x-1)equiv 0pmod p)。所以(xequiv 1pmod p)或(xequiv p-1pmod p)。
很遗憾,这个定理同样是必要不充分条件。但我们可以利用它来加强我们的条件,提高素数判定的准确性。
请注意,接下来我们讨论的是一个递归的过程,在递归过程中,指数(p)是在不断变化的,而模数(P)是始终不变的。初始时,(P=p),都是我们要判定的那个数。
首先,经过费马小定理的验证,我们已经确定(a^{p-1}equiv 1pmod P)。如果(p-1)是偶数,则有(left(a^{frac{p-1}{2}} ight)^2equiv 1pmod P)。那么可以将(a^{frac{p-1}{2}})视作“二次探测定理”里的(x),考虑(a^{frac{p-1}{2}})在(mod P)意义下的值:
-
如果既不是(1)也不是(P-1),说明(P)不是质数。直接反回
false
。 -
如果是(P-1),那么不能利用二次探测继续递归,说明目前无法验证(P)为合数,返回
true
。 -
如果是(1),则我们可以递归下去,令(p:=frac{p-1}{2}+1)。继续判断(a^{frac{p-1}{2}})。
如果每一层递归都做一次快速幂,复杂度会多个(log)。事实上,我们可以先把(p-1)每次除以(2)直到不能除为止,这相当于递归的最底层。然后在最底层做一次快速幂。再从底向上逐层推出当前层要判定的数的值。每一层的值都是在它下一层的基础上平方一下。求出每一层的数后,把它们存在一个vector
里,从后往前扫一遍,就能模拟出递归的过程了。
复杂度(O(log p))。
参考代码:
inline ll mul_mod(ll a,ll b,ll m=MOD) {
a%=m;b%=m;
ll c=(long double)a*b/m;
ll ans=a*b-c*m;
if(ans<0)ans+=m;
if(ans>=m)ans-=m;
return ans;
}
inline ll pow_mod(ll x,ll i,ll m=MOD){
ll y=1;
while(i){
if(i&1LL)y=mul_mod(y,x,m);
x=mul_mod(x,x,m);
i>>=1;
}
return y;
}
const int B[5]={2,3,5,7,11};
bool check(ll x,int b){
if(!(x&1))return false;
ll k=x-1;
while(!(k&1LL))k>>=1;
ll t=pow_mod(b,k,x);
vector<ll>vec;
while(k!=x-1)k<<=1,t=mul_mod(t,t,x),vec.pb(t);
if(t!=1)return false;
for(int i=vec.size()-2;i>=0;--i){
if(vec[i]!=1&&vec[i]!=x-1)return false;
if(vec[i]==x-1)return true;
}
return true;
}
bool is_prime(ll x,int RC=2){
/*
RC 容错
是一个1~5的整数,默认为2
数字越大,准确度越高,效率越低
*/
if(x==1)return false;
for(int i=0;i<RC;++i)if(x==B[i])return true;
for(int i=0;i<RC;++i)if(!check(x,B[i]))return false;
return true;
}
解决了素数判定,CASE 8~10的代码就相当简单了:
void solve5(){
int n=read();while(n--){
ll l=readll(),r=readll();
for(ll i=l;i<=r;++i){
if(is_prime(i,4))putchar('p');
else putchar('.');
}
putchar('
');
}
}
CASE (11sim 13)
套路和CASE 8~10差不多,只不过之前是对区间每个数做素数判定,现在是求区间每个数的莫比乌斯函数。
设(x=prod_{i=1}^{m}p_i^{c_i}),则莫比乌斯函数(mu(x))定义为:
与测试点(8)类似,测试点(11)我们可以用线性筛法轻松解决。
对于测试点(12,13),(l,rleq 10^{18}),需要一点数论技巧。
我们先求出([1,10^6])的所有质数,用它们筛掉([l,r])的数的所有(leq 10^6)的质因子,筛的过程中顺便维护一下每个数的(mu)值:每新增一个质因子,(mu)值要乘以(-1),若一个质因子出现了大于(1)次,则(mu)值直接等于(0)。
每个数被筛去这些小质因子后,留下的部分的质因子数最多不超过2。具体来讲,分如下几种情况:
-
没有其他质因子,即留下的数(=1)。
-
留下的数是一个质数。此时(mu)值要乘以(-1)。
-
留下的数是一个完全平方数,那它肯定是某个质数的平方。此时(mu)的值直接(=0)。
-
留下的数是两个不同的质数的积,(mu)值不变。
枚举每个(leq 10^6)的质因子,再枚举其倍数,这样做的复杂度是调和级数,即(O(nlog n))级别。输出答案的部分,对每个数可能要使用大素数判定,复杂度(O(nlog p)),(p)是值的大小,在(10^{18})级别。常数很大,我朴素地实现了一下,在本地要跑6秒。做了如下优化后勉强跑到2秒左右:
-
筛质数时筛到某个(mu)值已经为(0)的数,直接跳过。
-
最后的输出答案部分中,如果(mu)值已经为(0),则不要再进行素数判定。
参考代码:
using Math::is_prime;
using Math::sieve_mu;
ll val[1000006];
int mu[1000006];
void solve6(){
#define C(x) ((x==0)?('0'):((x==1)?'+':'-'))
sieve_mu();
int n=read();
//cout<<Math::cnt<<" "<<Math::p[78499]<<endl;
for(int i=1;i<=n;++i){
ll l=readll(),r=readll();
if(l<=1000000){
for(ll j=l;j<=min(r,1000000LL);++j)putchar(C(Math::mu[j]));
}
if(r<=1000000){
putchar('
');
continue;
}
l=max(l,1000001LL);
for(ll j=l;j<=r;++j)val[j-l]=j,mu[j-l]=1;
for(int j=1;j<=Math::cnt;++j){
ll x=l/Math::p[j],y=r/Math::p[j];
for(ll k=x;k<=y;++k){
if(k*Math::p[j]<l)continue;
//assert(k*Math::p[j]>=l);
//assert(k*Math::p[j]<=r);
if(!mu[k*Math::p[j]-l])continue;
int cnt=0;
while(val[k*Math::p[j]-l]%Math::p[j]==0){
val[k*Math::p[j]-l]/=(ll)Math::p[j];
cnt++;
if(cnt>1){mu[k*Math::p[j]-l]=0;break;}
}
if(cnt==1)mu[k*Math::p[j]-l]*=-1;
}
}
for(ll j=l;j<=r;++j){
//cout<<val[j-l]<<" "<<mu[j-l]<<endl;
if(val[j-l]==1||mu[j-l]==0)putchar(C(mu[j-l]));
else if(is_prime(val[j-l]))putchar(C(-mu[j-l]));
else{
ll s=sqrt(val[j-l]);
if(s*s==val[j-l])putchar(C(0));
else putchar(C(mu[j-l]));
}
}
putchar('
');
}
#undef C
}
CASE (14sim 16)
根据前面的经验,结合功能编号g
,不难想到这部分数据是要让大家判断区间([l,r])里每个数是不是(m)的原根。
简单讲一下原根的定义和性质。
根据欧拉定理我们知道,(a,m)互质时,(a^{varphi(m)}equiv 1pmod m)。所以显然存在一类(x)使得(a^xequiv 1pmod m)(如(x=varphi(m)))。但(varphi(m))不一定是最小的(x)。我们希望找到满足条件的最小的正整数(x),即(min{x|x>0,a^xequiv 1pmod m}),并将这样的(x)命名为(a)对(m)的阶。
对于一个(a),若(a)对(m)的阶为(varphi(m)),则称(a)是(m)的原根。
结论1:设(a)对(m)的阶为(y),那么({x|x>0,a^xequiv 1pmod m})这个集合一定等于({ky|kinmathbb{N_+}}),也即,(a^xequiv 1pmod mLeftrightarrow y|x)。
证明:
首先所有(y)的倍数(x)一定满足(a^xequiv 1pmod m)。因为根据定义,在(mod m)意义下,每乘以(y)个(a)就会形成一个周期。
考虑如果存在一个满足(a^xequiv 1pmod m)的(x)不被(y)整除,设(x=qy+r(0<r<y))。因为(a^{qy}equiv 1pmod m),所以(a^requiv 1pmod m),这与(y)是(a)对(m)的阶矛盾。故原命题不成立。即:所有满足(a^xequiv 1pmod m)的(x)都能被(y)整除。
结论2:一个数(m)有原根,当且仅当(m)可以表示为如下四种形式的任意一种:(2,4,p^k,2p^k),其中(p)是一个质数。
证明略。
结论3:如果(m)有一个原根(g),那么(g,g^2,dots,g^{varphi(m)})在(mod m)意义下的值一定恰好通过(m)的简化剩余系。
证明:
首先我们要证明,对任意(1leq j<kleqvarphi(m)),(g^j eq g^kpmod m)。反证法。假设(g^kequiv g^jpmod m),则(frac{g^k}{g^j}=g^{k-j}equiv1 pmod m)。即:存在一个(1<k-j<varphi(m))使得(g^{k-j}equiv 1pmod m),这与(g)是(m)的原根矛盾。
于是我们证明了我们取的数在(mod m)意义下两两不同。又因为我们取了(varphi(m))个数,所以这些数恰好通过(m)的简化剩余系。
回到本题,我们要判断一个数(a)是否是(m)的原根。朴素的算法我们要枚举(a)的(1simvarphi(m))次幂,复杂度是(O(varphi(m)))的,无法承受。
但是根据结论1,任何一个满足(a^xequiv 1pmod m)的(x)都是(a)对(m)的阶的倍数。(varphi(m))显然是其中一个满足条件的(x)。要找(a)对(m)的阶,我们只需要枚举(varphi(m))的约数即可。
更进一步,如果一个约数(x)满足(a^xequiv 1pmod m),那么(x)的倍数也满足。于是我们只需要枚举(varphi(m))的每个质因数(p_i),看(frac{varphi(m)}{p_i})是否满足条件。如果存在(frac{varphi(m)}{p_i})满足条件,说明(a)对(m)的阶一定小于(varphi(m)),(a)必不是(m)的原根。否则,可以说明(a)必是(m)的原根。
对于测试点(14),可以预处理出(varphi(m)=998244352)的所有质因数,判断每个(ain[l,r])是否是原根时都扫一遍所有质因数。复杂度(O(sqrt{m}+ntlog m)),其中(n=r-l+1),(t)是(varphi(m))的质因数个数,(log m)是判断时做快速幂的复杂度。
测试点(15)的区间长度很大,但模数(m)较小,考虑换一种思路。先求出(m)的一个原根(g)。因为(m)是质数,那么根据结论3,([1,m-1])的每个数都能被表示为(g^i(1leq i<m))。设当前要判断的数(a=g^x)。若(gcd(x,m-1) eq 1),则(a)不是原根。因为如果(gcd(x,m-1)=d),则(a^{frac{m-1}{d}}=g^{frac{x(m-1)}{d}}equiv 1pmod m)。否则不存在(m-1)的约数(d),使得(a^{frac{m-1}{d}}equiv 1pmod m),根据前面基于结论1的讨论,此时(a)一定是原根。
如果每次暴力求(gcd)显然还是太慢了,我们可以用(m-1)的所有质因数依次筛去它们的所有倍数。具体实现见代码。
测试点(16)里的?
,我们从(10^9)开始暴力求。注意,根据结论2以及基于对出题人的信任,此处我们还是只枚举质数就好。用下面的代码大约10分钟左右可以求出来。
void get_pr(){
int L=1e9,R=2e9,l=233333333,r=234133333;
freopen("software16.ans","r",stdin);
string str;cin>>str;cin>>str;cin>>str;
for(int P=L,c=0;P<=R;++P){
if(!is_prime(P))continue;
++c;
if(c%10000==0)cerr<<P<<endl;
div(P-1);
bool fail=0;
for(int i=l;i<=r;++i){
char ch=(is_pr(i,P)?'g':'.');
if(ch!=str[i-l]){
fail=1;break;
}
}
if(fail)continue;
cerr<<"! "<<P<<endl;
return;
}
}
CASE 14~16的代码:
int p[100005],cnt;
void div(int x){
cnt=0;
for(int i=2;i*i<=x;++i)if(x%i==0){
while(x%i==0)x/=i;
p[++cnt]=i;
}
if(x!=1)p[++cnt]=x;
}
inline int is_pr(int x,int P){
for(int i=1;i<=cnt;++i){
if(pow_mod(x,(P-1)/p[i],P)==1)return false;
}
return true;
}
bool vis[20000007];
int ind[20000007];
void solve7(){
div(998244353-1);
int n=read();while(n--){
int l=read(),r=read(),P;
if(l==233333333)P=1515343657,div(P-1);
else P=read();
if(P==998244353||P==1515343657){
for(int i=l;i<=r;++i){
if(is_pr(i,P))putchar('g');
else putchar('.');
}
}
else{
div(P-1);
for(int i=1;i<=cnt;++i){
for(int j=1;j*p[i]<=P-1;++j){
vis[j*p[i]]=1;
}
}
int g=0;
for(int i=1;i<P;++i)if(is_pr(i,P)){g=i;break;}
for(int i=1,t=g;i<=P-1;++i,t=(ll)t*g%P)ind[t]=i;//对数
for(int i=l;i<=r;++i){
if(vis[ind[i]])putchar('.');
else putchar('g');
}
}
putchar('
');
}
}