• 二次同余方程的解


    算法

    问题是解方程(x^2 equiv n (mod p)),其中(p)是奇质数。

    引理(n^{frac{p-1}2}equiv pm 1 (mod p))

    证明:由费马小定理,(n^{p-1}-1equiv (n^frac{p-1}2-1)(n^frac{p-1}2+1)equiv 0 (mod p) Rightarrow n^frac{p-1}2 equiv pm 1 (mod p))

    引理:方程有解当且仅当(n^frac{p-1}2 equiv 1)

    定理:设(a)满足(omega=a^2-n)不是模(p)的二次剩余,那么(x=(a+sqrt{omega})^frac{p+1}2)是二次剩余方程(x^2 equiv n (mod p))的解。

    证明:由((a+sqrt{omega})^p=a^p+omega^frac{p-1}2sqrt{w}=a-sqrt{omega}),前面的等号用二项式定理和(inom pi=0 (mod p),i eq 0,p)得到,后面的等号用费马小定理和(omega)是模(p)的二次非剩余。然后

    [(a+sqrt{omega})^{p+1}\ equiv(a+sqrtomega)^p(a+sqrtomega)\ equiv(a-sqrtomega)(a+sqrtomega)\ equiv a^2-omegaequiv n (mod p) ]

    在算法实现的时候,对(a​)的选择可以随机,因为大约有一半数是模(p​)的二次非剩余。然后快速幂即可。

    Timus Online Judge 例题

    1132. Square Root

    Time limit: 1.0 second
    Memory limit: 64 MB
    The number x is called a square root of a modulo n (root(a,n)) if x*x = a (mod n). Write the program to find the square root of number a by given modulo n.

    Input

    One number K in the first line is an amount of tests (K ≤ 100000). Each next line represents separate test, which contains integers a and n (1 ≤ a, n ≤ 32767, n is prime, a and n are relatively prime).

    Output

    For each input test the program must evaluate all possible values root(a,n) in the range from 1 to n − 1 and output them in increasing order in one separate line using spaces. If there is no square root for current test, the program must print in separate line: ‘No root’.

    Sample

    inputoutput
    5
    4 17
    3 7
    2 7
    14 31
    10007 20011
    
    2 15
    No root
    3 4
    13 18
    5382 14629
    
    Problem Author: Michael Medvedev

    p=2时特判输出1,不然会WA。龟速乘会TLE。

    #include<bits/stdc++.h>
    #define rg register
    #define il inline
    #define co const
    template<class T>il T read(){
        rg T data=0,w=1;
        rg char ch=getchar();
        while(!isdigit(ch)){
            if(ch=='-') w=-1;
            ch=getchar();
        }
        while(isdigit(ch))
            data=data*10+ch-'0',ch=getchar();
        return data*w;
    }
    template<class T>il T read(rg T&x){
        return x=read<T>();
    }
    typedef long long ll;
    
    ll n,mod;
    ll add(ll x,ll y){
    	return (x+y)%mod;
    }
    //ll mul(ll x,ll y){
    //	x%=mod,y%=mod;
    //	ll re=0;
    //	for(;y;y>>=1,x=add(x,x))
    //		if(y&1) re=add(re,x);
    //	return re;
    //}
    ll mul(ll x,ll y){
    	return x*y%mod;
    }
    ll pow(ll x,ll k){
    	x%=mod,k%=(mod-1);
    	ll re=1;
    	for(;k;k>>=1,x=mul(x,x))
    		if(k&1) re=mul(re,x);
    	return re;
    }
    ll omega;
    struct complex{ll a,b;};
    complex operator*(co complex&x,co complex&y){
    	return (complex){add(mul(x.a,y.a),mul(x.b,mul(y.b,omega))),add(mul(x.a,y.b),mul(x.b,y.a))};
    }
    complex operator^(complex x,ll k){
    	complex re=(complex){1,0};
    	for(;k;k>>=1,x=x*x)
    		if(k&1) re=re*x;
    	return re;
    }
    ll sqrt(ll n){
    	if(mod==2) return 1;
    	n%=mod;
    	if(!n) return 0;
    	ll a=rand()%mod;
    	while(pow(add(mul(a,a),mod-n),(mod-1)/2)!=mod-1) a=rand()%mod;
    	omega=add(mul(a,a),mod-n);
    	return ((complex){a,1}^(mod+1)/2).a;
    }
    int main(){
    //	freopen(".in","r",stdin),freopen(".out","w",stdout);
    	int kase=read<int>();
    	while(kase--){
    		read(n),read(mod);
    		if(pow(n,(mod-1)/2)!=1) {puts("No root");continue;}
    		ll ans1=sqrt(n),ans2=mod-ans1;
    		if(ans1>ans2) std::swap(ans1,ans2);
    		if(ans2!=ans1) printf("%lld %lld
    ",ans1,ans2);
    		else printf("%lld
    ",ans1);
    	}
    	return 0;
    }
    
  • 相关阅读:
    1085 PAT单位排行
    安装MongoDB并且添加用户
    同源政策,发送请求时携带cookie信息
    博客园文章编辑时实现语法高亮
    template中的时间格式如何修改
    Node模块下载路径的更改设置
    JavaScript--遍历
    JavaScript--作用域
    JavaScript--arguments
    JavaScript--apply&call
  • 原文地址:https://www.cnblogs.com/autoint/p/quadratic_residue.html
Copyright © 2020-2023  润新知