「余姚中学 2019 联测 Day 6」解码

「余姚中学 2019 联测 Day 6」解码

先不考虑求(p,q)

根据人人都知道的欧拉定理(x^cequiv x^{cmod varphi(n)} (mod n))

那么(varphi(n)=(p-1)(q-1))，而((c,varphi(n))=1)

所以求出(frac{1}{c} pmod {varphi(n)})

带入原来的式子((x^c)^{frac{1}{c}pmod {varphi(n)}}equiv x^1pmod n)

即(x=m^{frac{1}{c}pmod {varphi(n)}}pmod n)

求逆最好用扩展欧几里得算法，复杂度为(O(log n))

那么直接快速幂即可，但是如果快速幂套快速乘复杂度为(O(log ^2))，实际常数极大，很有可能超时(如果用long double O(1)快速乘另谈。。。)

由于知道(p,q)可以分别求出(m^{frac{1}{c}pmod {varphi(n)}}pmod p)，(m^{frac{1}{c}pmod {varphi(n)}}pmod q)

然后中国剩余定理合并，即(O(log n))

---

现在问题是求(p,q)

对于前3档分，由于素数密度是(O(log n))的，所以(sqrt n -p)期望只有(log n)

而对于最后一档分，考虑更好表示，枚举(p+q)解出答案，发现

(4nleq (p+q)^2= (q-p)^2+4pqleq lambda^2+4n)

即(2sqrt{n}leq (p+q)le sqrt{lambda^2+4n})

由于(nge p^2ge 10^{18},lambda leq 3cdot 10^5)，这个范围实际非常非常非常非常小，大概只有(22)

tips:实际上(4n)可能爆long long

枚举(p+q)之后，(O(1))解出(p,q)即可

竟然有人问怎么解(p,q)，我震惊了

(q-p=sqrt{(p+q)^2-4n})，然后判一下是不是整数就好了

写的时候害怕sqrt炸精度，很多奇怪的语句请忽略

#include<bits/stdc++.h>
using namespace std;

#define reg register
typedef long long ll;
typedef unsigned long long ull;
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)

#define pb push_back
#define Mod1(x) ((x>=P)&&(x-=P))
#define Mod2(x) ((x<0)&&(x+=P))

template <class T> inline void cmin(T &a,T b){ if(a>b) a=b; }
template <class T> inline void cmax(T &a,T b){ if(a<b) a=b; }

char IO;
template<class T=int> T rd(){
	T s=0; int f=0;
	while(!isdigit(IO=getchar())) if(IO=='-') f=1;
	do s=(s<<1)+(s<<3)+(IO^'0');
	while(isdigit(IO=getchar()));
	return f?-s:s;
}

const int N=3e5+10;

ll n,m,c;

ll qmul(ll x,ll k,ll P){
	k=(k%P+P)%P;
	ll res=0;
	for(;k;k>>=1,x=(x+x)%P) if(k&1) res=(res+x)%P;
	return res;
}

ll qpow(ll x,ll k,ll P) {
	ll res=1;
	for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
	return res;
}

void Exgcd(ll a,ll b,ll &x,ll &y){
	if(b==0) x=1,y=0;
	else Exgcd(b,a%b,y,x),y-=a/b*x;
}

ll Inv(ll a,ll P){
	ll x,y;
	Exgcd(a,P,x,y);
	return (x%P+P)%P;
}

int main(){
	freopen("rsa.in","r",stdin),freopen("rsa.out","w",stdout);
	rep(kase,1,rd()) {
		n=rd<ll>(),m=rd<ll>(),c=rd<ll>();
		ll T=sqrt(n),p=-1,q;
		for(ll i=T;i>=T-100;--i) if(n%i==0) {
			p=i,q=n/i;
			break;
		}
		if(p==-1) {
			// 2*sqrt(n) <= p+q <= sqrt(4*n+lambda * lambda)
			//ll R=ceil(sqrt((long double)4*n+9e10)+1);
			ll L=ceil(sqrt(n+0.5));
			if((L-1)*(L-1)>=n) L--;

			for(ll x=2*L;;++x) {
				ull y=x*x-4*(ull)n; 
				// x=p+q
				// t=q-p
				ull t=sqrt(y);
				if((t+1)*(t+1)==y) t++;
				if((t-1)*(t-1)==y) t--;
				if(t*t==y) {
					p=(x-t)/2;
					q=(x+t)/2;
					break;
				}
			}
		}

		ll t=Inv(c,(p-1)*(q-1));
		// t= 1/c (mod phi(n))
		 
		ll k1=p,b1=qpow(m%p,t,p);
		ll k2=q,b2=qpow(m%q,t,q);
		// k1 x + b1 = k2 y + b2
		// k1 x = b2-b1 (mod k2)
		// x= (b2-b1)/k1;
		// x' = (b2-b1)/k1 (mod k2) * k1 + b1
		ll inv=qpow(k1,k2-2,k2);
		b1=(b2-b1)%k2*inv%k2  * k1+b1;
		k1*=k2;
		b1=(b1%k1+k1)%k1;
		printf("%lld
",b1);
	}
}