Tonelli-Shanks algorithm Quadratic Residue System Solution (Ural 1132. square root)

Explanation and proof on Wikipedia: http://en.wikipedia.org/wiki/Tonelli%E2%80%93Shanks_algorithm

TheTonelli-shanksAlgorithm (referred to by shanks as the ressol algorithm) is used within modular arithmetic to solve a congruence of the form

WhereNIs a quadratic residue (modP), AndPIs an odd prime.

Tonelli-Shanks cannot be used for composite moduli; finding square roots modulo composite numbers is a computational problem equivalent to integer factorization.^[1]

An equivalent, but slightly more redundant version of this algorithm was developed by Albert to Tonelli in 1891. The version discussed here was developed independently by Daniel shanksin 1973, who explained:

"My tardiness in learning of these historical references was because I had lent Volume 1 of Dickson's history to a friend and it was never returned ."^[2]

(Note: all are taken to mean, unless indicated otherwise). [edit] The Algorithm

Inputs:P, An odd prime.N, An integer which is a quadratic residue (modP), Meaning that the legendh symbol.

Outputs:R, An integer satisfying.

1. Factor out powers of 2 fromP−1, definingQAndSAs:QOdd. Note that if,I. e., Then solutions are given directly.
2. SelectZSuch that the legendh symbol (that is,ZShocould be a quadratic non-residue moduloP), And set.
3. Let
4. Loop:
   1. If, returnR.
   2. Otherwise, find the lowestI, Such that;E.g.Via repeated squaring.
   3. Let, and set and.

Once you have solved the congruenceRThe second solution isP−R.

Example

Solving the congruence. It is clear that is odd, and since, 10 is a quadratic residue (by Euler's criterion ).

• Step 1: Observe so ,.

• Step 2: Take as the quadratic nonresidue (2 is a quadratic nonresidue since (again, Euler's criterion). Set

• Step 3:

• Step 4: Now we start the loop: So;I. e. 
  • Let, so.
  • Set. Set, and
  • We restart the loop, and since we are done, returning

Indeed, observe that and naturally also. So the algorithm yields two solutions to our congruence.

Proof

First write. Now write and, observing that. This latter congruence will be true after every iteration of the algorithm's main loop. If at any point, then and the algorithm terminates.

If, then consider, a quadratic non-residue of. Let. Then and, which shows that the order of is.

Similarly we have, so the order of divides. Suppose the order of is. Since is a square modulo, is also a square, and hence.

Now we set and with this, And. As before, holds; however with this construction both and have order. This implies that has order.

If then, and the algorithm stops, returning. Else, we restart the loop with analogous definitions of, and until we arrive at an that equals 0. Since the sequenceSIs strictly decreasing the algorithm terminates.

The general process of the algorithm has been clearly stated, and then the process can be simulated. However, for this question on Ural, we need to determine N = 2. For details, see the code:

//#pragma comment(linker,"/STACK:327680000,327680000")#include <iostream>#include <cstdio>#include <cmath>#include <vector>#include <cstring>#include <algorithm>#include <string>#include <set>#include <functional>#include <numeric>#include <sstream>#include <stack>#include <map>#include <queue>#define CL(arr, val)    memset(arr, val, sizeof(arr))#define REP(i, n)       for((i) = 0; (i) < (n); ++(i))#define FOR(i, l, h)    for((i) = (l); (i) <= (h); ++(i))#define FORD(i, h, l)   for((i) = (h); (i) >= (l); --(i))#define L(x)    (x) << 1#define R(x)    (x) << 1 | 1#define MID(l, r)   (l + r) >> 1#define Min(x, y)   (x) < (y) ? (x) : (y)#define Max(x, y)   (x) < (y) ? (y) : (x)#define E(x)        (1 << (x))#define iabs(x)     (x) < 0 ? -(x) : (x)#define OUT(x)  printf("%I64d\n", x)#define Read()  freopen("data.in", "r", stdin)#define Write() freopen("data.out", "w", stdout);typedef long long LL;const double eps = 1e-8;const double pi = acos(-1.0);const double inf = ~0u>>2;using namespace std;LL MOD;LL mod_exp(LL a, LL b) {    LL res = 1;    while(b > 0) {        if(b&1)    res = (res*a)%MOD;        a = (a*a)%MOD;        b >>= 1;    }    return res;}LL solve(LL n, LL p) {    int Q = p - 1, S = 0;    while(Q%2 == 0) {        Q >>= 1;        S++;    }    if(S == 1) {return mod_exp(n, (p + 1)/4);}    int z;    while(1) {        z = 1 + rand()%(p - 1);        if(mod_exp(z, (p - 1)/2) != 1)   break;    }    LL c = mod_exp(z, Q);    LL R = mod_exp(n, (Q + 1)/2);    LL t = mod_exp(n, Q);    LL M = S, b, i;    while(1) {        if(t%p == 1)  break;        for(i = 1; i < M; ++i) {            if(mod_exp(t, 1<<i) == 1)    break;        }        b = mod_exp(c, 1<<(M-i-1));        R = (R*b)%p;        t = (t*b*b)%p;        c = (b*b)%p;        M = i;    }    return (R%p + p)%p;}int main() {    //Read();    int T, a, n, i;    scanf("%d", &T);    while(T--) {        scanf("%d%d", &a, &n);        if(n == 2) {    //***            if(a%n == 1)    printf("1\n");            else    puts("No root");            continue;        }        MOD = n;        if(mod_exp(a, (n-1)/2) != 1)    {puts("No root"); continue; }        i = solve(a, n);        if(i == n - i)  printf("%d\n", i);        else    printf("%d %d\n", min(i, n - i), max(i, n - i));    }    return 0;}

 

Poj 1808 legendh symbol determines view code

// # Pragma comment (linker, "/Stack: 327680000,327680000 ") # include <iostream> # include <cstdio> # include <cmath> # include <vector> # include <cstring> # include <algorithm> # include <string> # include <set> # include <functional> # include <numeric> # include <sstream> # include <stack> # include <map> # include <queue> # define Cl (ARR, val) memset (ARR, Val, sizeof (ARR) # define rep (I, n) for (I) = 0; (I) <(N ); ++ (I) # define For (I, L, H) for (I) = (l); (I) <= (h); ++ (I) # define Ford (I, h, L) for (I) = (h); (I)> = (l); -- (I) # define L (x) (X) <1 # define R (x) <1 | 1 # define mid (L, R) (L + r)> 1 # define min (X, y) (x) <(y )? (X): (y) # define max (x, y) (x) <(y )? (Y): (x) # define e (x) (1 <(x) # define iabs (x) <0? -(X): (x) # define out (x) printf ("% i64d \ n", x) # define read () freopen ("data. in "," r ", stdin) # define write () freopen (" data. out "," W ", stdout); typedef long ll; const double EPS = 1e-8; const double Pi = ACOs (-1.0); const double INF = ~ 0u> 2; using namespace STD; ll MOD; ll mod_exp (ll a, LL B) {ll res = 1; while (B> 0) {If (B & 1) res = (RES * A) % MOD; A = (A * A) % MOD; B >>= 1;} return res ;} int main () {// read (); ll a, P; int T, CAS = 0; scanf ("% d", & T); While (t --) {scanf ("% LLD", & A, & P); mod = P; printf ("Scenario # % d: \ n", ++ CAS ); if (mod_exp (a % P + p) % P, (p-1)/2) = 1) puts ("1 "); // note that a may be a negative else puts ("-1"); cout <Endl ;}}