Open In App

Find Square Root under Modulo p | Set 2 (Shanks Tonelli algorithm)

Last Updated : 20 Dec, 2022
Comments
Improve
Suggest changes
Like Article
Like
Report

Given a number ‘n’ and a prime ‘p’, find square root of n under modulo p if it exists. 
Examples: 
 

Input: n = 2, p = 113
Output: 62
62^2 = 3844  and  3844 % 113 = 2

Input:  n = 2, p = 7
Output: 3 or 4
3 and 4 both are square roots of 2 under modulo
7 because (3*3) % 7 = 2 and (4*4) % 7 = 2

Input:  n = 2, p = 5
Output: Square root doesn't exist


We have discussed Euler’s criterion to check if square root exists or not. We have also discussed a solution that works only when p is in form of 4*i + 3 
In this post, Shank Tonelli’s algorithm is discussed that works for all types of inputs.
Algorithm steps to find modular square root using shank Tonelli’s algorithm :
1) Calculate n ^ ((p - 1) / 2) (mod p), it must be 1 or p-1, if it is p-1, then modular square root is not possible.
2) Then after write p-1 as (s * 2^e) for some integer s and e, where s must be an odd number and both s and e should be positive.
3) Then find a number q such that q ^ ((p - 1) / 2) (mod p) = -1 
4) Initialize variable x, b, g and r by following values 

   x = n ^ ((s + 1) / 2 (first guess of square root)
   b = n ^ s                
   g = q ^ s
   r = e   (exponent e will decrease after each updation) 


5) Now loop until m > 0 and update value of x, which will be our final answer. 

   Find least integer m such that b^(2^m) = 1(mod p)  and  0 <= m <= r – 1 
   If m = 0, then we found correct answer and return x as result
   Else update x, b, g, r as below
       x = x * g ^ (2 ^ (r – m - 1))
       b = b * g ^(2 ^ (r - m))
       g = g ^ (2 ^ (r - m))
       r = m 


so if m becomes 0 or b becomes 1, we terminate and print the result. This loop guarantees to terminate because value of m is decreased each time after updation.
Following is the implementation of above algorithm. 
 

C++
// C++ program to implement Shanks Tonelli algorithm for
// finding Modular  Square Roots
#include <bits/stdc++.h>
using namespace std;

//  utility function to find pow(base, exponent) % modulus
int pow(int base, int exponent, int modulus)
{
    int result = 1;
    base = base % modulus;
    while (exponent > 0)
    {
        if (exponent % 2 == 1)
           result = (result * base)% modulus;
        exponent = exponent >> 1;
        base = (base * base) % modulus;
    }
    return result;
}

//  utility function to find gcd
int gcd(int a, int b)
{
    if (b == 0)
        return a;
    else
        return gcd(b, a % b);
}

//  Returns k such that b^k = 1 (mod p)
int order(int p, int b)
{
    if (gcd(p, b) != 1)
    {
        printf("p and b are not co-prime.\n");
        return -1;
    }

    //  Initializing k with first odd prime number
    int k = 3;
    while (1)
    {
        if (pow(b, k, p) == 1)
            return k;
        k++;
    }
}

//  function return  p - 1 (= x argument) as  x * 2^e,
//  where x will be odd  sending e as reference because
//  updation is needed in actual e
int convertx2e(int x, int& e)
{
    e = 0;
    while (x % 2 == 0)
    {
        x /= 2;
        e++;
    }
    return x;
}

//  Main function for finding the modular square root
int STonelli(int n, int p)
{
    //  a and p should be coprime for finding the modular
    // square root
    if (gcd(n, p) != 1)
    {
        printf("a and p are not coprime\n");
        return -1;
    }

    //  If below expression return (p - 1)  then modular
    // square root is not possible
    if (pow(n, (p - 1) / 2, p) == (p - 1))
    {
        printf("no sqrt possible\n");
        return -1;
    }

    //  expressing p - 1, in terms of s * 2^e,  where s
    // is odd number
    int s, e;
    s = convertx2e(p - 1, e);

    //  finding smallest q such that q ^ ((p - 1) / 2)
    //  (mod p) = p - 1
    int q;
    for (q = 2; ; q++)
    {
        // q - 1 is in place of  (-1 % p)
        if (pow(q, (p - 1) / 2, p) == (p - 1))
            break;
    }

    //  Initializing variable x, b and g
    int x = pow(n, (s + 1) / 2, p);
    int b = pow(n, s, p);
    int g = pow(q, s, p);

    int r = e;

    // keep looping until b become 1 or m becomes 0
    while (1)
    {
        int m;
        for (m = 0; m < r; m++)
        {
            if (order(p, b) == -1)
                return -1;

            //  finding m such that b^ (2^m) = 1
            if (order(p, b) == pow(2, m))
                break;
        }
        if (m == 0)
            return x;

        // updating value of x, g and b according to
        // algorithm
        x = (x * pow(g, pow(2, r - m - 1), p)) % p;
        g = pow(g, pow(2, r - m), p);
        b = (b * g) % p;

        if (b == 1)
            return x;
        r = m;
    }
}

//  driver program to test above function
int main()
{
    int n = 2;

    // p should be prime
    int p = 113;

    int x = STonelli(n, p);

    if (x == -1)
        printf("Modular square root is not exist\n");
    else
        printf("Modular square root of %d and %d is %d\n",
                n, p, x);
}
Java
// Java program to implement Shanks
// Tonelli algorithm for finding 
// Modular Square Roots 
import java.util.*;
class GFG
{
    static int z = 0;
    
// utility function to find
// pow(base, exponent) % modulus 
static int pow1(int base1, 
    int exponent, int modulus) 
{ 
    int result = 1; 
    base1 = base1 % modulus; 
    while (exponent > 0) 
    { 
        if (exponent % 2 == 1) 
            result = (result * base1) % modulus; 
        exponent = exponent >> 1; 
        base1 = (base1 * base1) % modulus; 
    } 
    return result; 
} 

// utility function to find gcd 
static int gcd(int a, int b) 
{ 
    if (b == 0) 
        return a; 
    else
        return gcd(b, a % b); 
} 

// Returns k such that b^k = 1 (mod p) 
static int order(int p, int b) 
{ 
    if (gcd(p, b) != 1) 
    { 
        System.out.println("p and b are" + 
                            "not co-prime."); 
        return -1; 
    } 

    // Initializing k with first
    // odd prime number 
    int k = 3; 
    while (true) 
    { 
        if (pow1(b, k, p) == 1) 
            return k; 
        k++; 
    } 
} 

// function return p - 1 (= x argument)
// as x * 2^e, where x will be odd 
// sending e as reference because 
// updation is needed in actual e 
static int convertx2e(int x) 
{ 
    z = 0;
    while (x % 2 == 0) 
    { 
        x /= 2; 
        z++; 
    } 
    return x; 
} 

// Main function for finding 
// the modular square root 
static int STonelli(int n, int p) 
{ 
    // a and p should be coprime for  
    // finding the modular square root 
    if (gcd(n, p) != 1) 
    { 
        System.out.println("a and p are not coprime"); 
        return -1; 
    } 

    // If below expression return (p - 1) then modular 
    // square root is not possible 
    if (pow1(n, (p - 1) / 2, p) == (p - 1)) 
    { 
        System.out.println("no sqrt possible"); 
        return -1; 
    } 

    // expressing p - 1, in terms of  
    // s * 2^e, where s is odd number 
    int s, e; 
    s = convertx2e(p - 1);
    e = z;

    // finding smallest q such that q ^ ((p - 1) / 2) 
    // (mod p) = p - 1 
    int q; 
    for (q = 2; ; q++) 
    { 
        // q - 1 is in place of (-1 % p) 
        if (pow1(q, (p - 1) / 2, p) == (p - 1)) 
            break; 
    } 

    // Initializing variable x, b and g 
    int x = pow1(n, (s + 1) / 2, p); 
    int b = pow1(n, s, p); 
    int g = pow1(q, s, p); 

    int r = e; 

    // keep looping until b 
    // become 1 or m becomes 0 
    while (true) 
    { 
        int m; 
        for (m = 0; m < r; m++) 
        { 
            if (order(p, b) == -1) 
                return -1; 

            // finding m such that b^ (2^m) = 1 
            if (order(p, b) == Math.pow(2, m)) 
                break; 
        } 
        if (m == 0) 
            return x; 

        // updating value of x, g and b
        // according to algorithm 
        x = (x * pow1(g, (int)Math.pow(2, 
                            r - m - 1), p)) % p; 
        g = pow1(g, (int)Math.pow(2, r - m), p); 
        b = (b * g) % p; 

        if (b == 1) 
            return x; 
        r = m; 
    } 
} 

// Driver code 
public static void main (String[] args) 
{

    int n = 2; 

    // p should be prime 
    int p = 113; 

    int x = STonelli(n, p); 

    if (x == -1)
        System.out.println("Modular square" + 
                        "root is not exist\n"); 
    else
        System.out.println("Modular square root of " +
                            n + " and " + p + " is " +
                            x + "\n"); 
} 
}

// This code is contributed by mits
Python3
# Python3 program to implement Shanks Tonelli
# algorithm for finding Modular Square Roots 

# utility function to find pow(base, 
# exponent) % modulus 
def pow1(base, exponent, modulus): 

    result = 1; 
    base = base % modulus; 
    while (exponent > 0): 
        if (exponent % 2 == 1):
            result = (result * base) % modulus; 
        exponent = int(exponent) >> 1; 
        base = (base * base) % modulus; 

    return result; 

# utility function to find gcd 
def gcd(a, b): 
    if (b == 0): 
        return a; 
    else:
        return gcd(b, a % b); 

# Returns k such that b^k = 1 (mod p) 
def order(p, b): 

    if (gcd(p, b) != 1):
        print("p and b are not co-prime.\n"); 
        return -1; 

    # Initializing k with first 
    # odd prime number 
    k = 3; 
    while (True): 
        if (pow1(b, k, p) == 1): 
            return k; 
        k += 1; 

# function return p - 1 (= x argument) as 
# x * 2^e, where x will be odd sending e 
# as reference because updation is needed
# in actual e 
def convertx2e(x):
    z = 0; 
    while (x % 2 == 0):
        x = x / 2; 
        z += 1; 
        
    return [x, z]; 

# Main function for finding the
# modular square root 
def STonelli(n, p): 

    # a and p should be coprime for
    # finding the modular square root 
    if (gcd(n, p) != 1):
        print("a and p are not coprime\n"); 
        return -1; 

    # If below expression return (p - 1) then
    # modular square root is not possible 
    if (pow1(n, (p - 1) / 2, p) == (p - 1)):
        print("no sqrt possible\n"); 
        return -1; 

    # expressing p - 1, in terms of s * 2^e, 
    # where s is odd number 
    ar = convertx2e(p - 1);
    s = ar[0];
    e = ar[1];

    # finding smallest q such that 
    # q ^ ((p - 1) / 2) (mod p) = p - 1 
    q = 2; 
    while (True):
        
        # q - 1 is in place of (-1 % p) 
        if (pow1(q, (p - 1) / 2, p) == (p - 1)): 
            break;
        q += 1;

    # Initializing variable x, b and g 
    x = pow1(n, (s + 1) / 2, p); 
    b = pow1(n, s, p); 
    g = pow1(q, s, p); 

    r = e; 

    # keep looping until b become 
    # 1 or m becomes 0 
    while (True):
        m = 0; 
        while (m < r):
            if (order(p, b) == -1): 
                return -1; 

            # finding m such that b^ (2^m) = 1 
            if (order(p, b) == pow(2, m)): 
                break;
            m += 1;

        if (m == 0): 
            return x; 

        # updating value of x, g and b 
        # according to algorithm 
        x = (x * pow1(g, pow(2, r - m - 1), p)) % p; 
        g = pow1(g, pow(2, r - m), p); 
        b = (b * g) % p; 

        if (b == 1): 
            return x; 
        r = m; 

# Driver Code
n = 2; 

# p should be prime 
p = 113; 

x = STonelli(n, p); 

if (x == -1): 
    print("Modular square root is not exist\n"); 
else:
    print("Modular square root of", n, 
          "and", p, "is", x); 
        
# This code is contributed by mits
C#
// C# program to implement Shanks
// Tonelli algorithm for finding 
// Modular Square Roots 
using System;

class GFG
{
    
static int z=0;
    
// utility function to find 
// pow(base, exponent) % modulus 
static int pow1(int base1, 
        int exponent, int modulus) 
{ 
    int result = 1; 
    base1 = base1 % modulus; 
    while (exponent > 0) 
    { 
        if (exponent % 2 == 1) 
            result = (result * base1) % modulus; 
        exponent = exponent >> 1; 
        base1 = (base1 * base1) % modulus; 
    } 
    return result; 
} 

// utility function to find gcd 
static int gcd(int a, int b) 
{ 
    if (b == 0) 
        return a; 
    else
        return gcd(b, a % b); 
} 

// Returns k such that b^k = 1 (mod p) 
static int order(int p, int b) 
{ 
    if (gcd(p, b) != 1) 
    { 
        Console.WriteLine("p and b are" +
                            "not co-prime."); 
        return -1; 
    } 

    // Initializing k with 
    // first odd prime number 
    int k = 3; 
    while (true) 
    { 
        if (pow1(b, k, p) == 1) 
            return k; 
        k++; 
    } 
} 

// function return p - 1 (= x argument) 
// as x * 2^e, where x will be odd sending
// e as reference because updation is 
// needed in actual e 
static int convertx2e(int x) 
{ 
    z = 0;
    while (x % 2 == 0) 
    { 
        x /= 2; 
        z++; 
    } 
    return x; 
} 

// Main function for finding
// the modular square root 
static int STonelli(int n, int p) 
{ 
    // a and p should be coprime for 
    // finding the modular square root 
    if (gcd(n, p) != 1) 
    { 
        Console.WriteLine("a and p are not coprime"); 
        return -1; 
    } 

    // If below expression return (p - 1) then 
    // modular square root is not possible 
    if (pow1(n, (p - 1) / 2, p) == (p - 1)) 
    { 
        Console.WriteLine("no sqrt possible"); 
        return -1; 
    } 

    // expressing p - 1, in terms of s * 2^e,
    //  where s is odd number 
    int s, e; 
    s = convertx2e(p - 1);
    e=z;

    // finding smallest q such that q ^ ((p - 1) / 2) 
    // (mod p) = p - 1 
    int q; 
    for (q = 2; ; q++) 
    { 
        // q - 1 is in place of (-1 % p) 
        if (pow1(q, (p - 1) / 2, p) == (p - 1)) 
            break; 
    } 

    // Initializing variable x, b and g 
    int x = pow1(n, (s + 1) / 2, p); 
    int b = pow1(n, s, p); 
    int g = pow1(q, s, p); 

    int r = e; 

    // keep looping until b become 
    // 1 or m becomes 0 
    while (true) 
    { 
        int m; 
        for (m = 0; m < r; m++) 
        { 
            if (order(p, b) == -1) 
                return -1; 

            // finding m such that b^ (2^m) = 1 
            if (order(p, b) == Math.Pow(2, m)) 
                break; 
        } 
        if (m == 0) 
            return x; 

        // updating value of x, g and b  
        // according to algorithm 
        x = (x * pow1(g, (int)Math.Pow(2, r - m - 1), p)) % p; 
        g = pow1(g, (int)Math.Pow(2, r - m), p); 
        b = (b * g) % p; 

        if (b == 1) 
            return x; 
        r = m; 
    } 
} 

// Driver code 
static void Main() 
{ 
    int n = 2; 

    // p should be prime 
    int p = 113; 

    int x = STonelli(n, p); 

    if (x == -1)
        Console.WriteLine("Modular square root" +
                            "is not exist\n"); 
    else
        Console.WriteLine("Modular square root of" + 
                        "{0} and {1} is {2}\n", n, p, x); 
} 
}

// This code is contributed by mits
PHP
<?php
// PHP program to implement Shanks Tonelli 
// algorithm for finding Modular Square Roots 

// utility function to find pow(base,
// exponent) % modulus 
function pow1($base, $exponent, $modulus) 
{ 
    $result = 1; 
    $base = $base % $modulus; 
    while ($exponent > 0) 
    { 
        if ($exponent % 2 == 1) 
        $result = ($result * $base) % $modulus; 
        $exponent = $exponent >> 1; 
        $base = ($base * $base) % $modulus; 
    } 
    return $result; 
} 

// utility function to find gcd 
function gcd($a, $b) 
{ 
    if ($b == 0) 
        return $a; 
    else
        return gcd($b, $a % $b); 
} 

// Returns k such that b^k = 1 (mod p) 
function order($p, $b) 
{ 
    if (gcd($p, $b) != 1) 
    { 
        print("p and b are not co-prime.\n"); 
        return -1; 
    } 

    // Initializing k with first
    // odd prime number 
    $k = 3; 
    while (1) 
    { 
        if (pow1($b, $k, $p) == 1) 
            return $k; 
        $k++; 
    } 
} 

// function return p - 1 (= x argument) as
// x * 2^e, where x will be odd sending e 
// as reference because updation is needed
// in actual e 
function convertx2e($x, &$e) 
{ 
    $e = 0; 
    while ($x % 2 == 0) 
    { 
        $x = (int)($x / 2); 
        $e++; 
    } 
    return $x; 
} 

// Main function for finding the 
// modular square root 
function STonelli($n, $p) 
{ 
    // a and p should be coprime for 
    // finding the modular square root 
    if (gcd($n, $p) != 1) 
    { 
        print("a and p are not coprime\n"); 
        return -1; 
    } 

    // If below expression return (p - 1) then 
    // modular square root is not possible 
    if (pow1($n, ($p - 1) / 2, $p) == ($p - 1)) 
    { 
        printf("no sqrt possible\n"); 
        return -1; 
    } 

    // expressing p - 1, in terms of s * 2^e, 
    // where s is odd number 
    $e = 0; 
    $s = convertx2e($p - 1, $e); 

    // finding smallest q such that 
    // q ^ ((p - 1) / 2) (mod p) = p - 1 
    $q = 2; 
    for (; ; $q++) 
    { 
        // q - 1 is in place of (-1 % p) 
        if (pow1($q, ($p - 1) / 2, $p) == ($p - 1)) 
            break; 
    } 

    // Initializing variable x, b and g 
    $x = pow1($n, ($s + 1) / 2, $p); 
    $b = pow1($n, $s, $p); 
    $g = pow1($q, $s, $p); 

    $r = $e; 

    // keep looping until b become
    // 1 or m becomes 0 
    while (1) 
    { 
        $m = 0; 
        for (; $m < $r; $m++) 
        { 
            if (order($p, $b) == -1) 
                return -1; 

            // finding m such that b^ (2^m) = 1 
            if (order($p, $b) == pow(2, $m)) 
                break; 
        } 
        if ($m == 0) 
            return $x; 

        // updating value of x, g and b 
        // according to algorithm 
        $x = ($x * pow1($g, pow(2, $r - $m - 1), $p)) % $p; 
        $g = pow1($g, pow(2, $r - $m), $p); 
        $b = ($b * $g) % $p; 

        if ($b == 1) 
            return $x; 
        $r = $m; 
    } 
} 

// Driver Code
$n = 2; 

// p should be prime 
$p = 113; 

$x = STonelli($n, $p); 

if ($x == -1) 
    print("Modular square root is not exist\n"); 
else
    print("Modular square root of " . 
          "$n and $p is $x\n"); 
    
// This code is contributed by mits
?>
JavaScript
<script>

// JavaScript program to implement Shanks
// Tonelli algorithm for finding 
// Modular Square Roots 

    let z = 0;
      
// utility function to find
// pow(base, exponent) % modulus 
function pow1(base1, 
    exponent, modulus) 
{ 
    let result = 1; 
    base1 = base1 % modulus; 
    while (exponent > 0) 
    { 
        if (exponent % 2 == 1) 
            result = (result * base1) % modulus; 
        exponent = exponent >> 1; 
        base1 = (base1 * base1) % modulus; 
    } 
    return result; 
} 
  
// utility function to find gcd 
function gcd(a, b) 
{ 
    if (b == 0) 
        return a; 
    else
        return gcd(b, a % b); 
} 
  
// Returns k such that b^k = 1 (mod p) 
function order(p, b) 
{ 
    if (gcd(p, b) != 1) 
    { 
        document.write("p and b are" + 
                            "not co-prime." + "<br/>"); 
        return -1; 
    } 
  
    // Initializing k with first
    // odd prime number 
    let k = 3; 
    while (true) 
    { 
        if (pow1(b, k, p) == 1) 
            return k; 
        k++; 
    } 
} 
  
// function return p - 1 (= x argument)
// as x * 2^e, where x will be odd 
// sending e as reference because 
// updation is needed in actual e 
function convertx2e(x) 
{ 
    z = 0;
    while (x % 2 == 0) 
    { 
        x /= 2; 
        z++; 
    } 
    return x; 
} 
  
// Main function for finding 
// the modular square root 
function STonelli(n, p) 
{ 
    // a and p should be coprime for  
    // finding the modular square root 
    if (gcd(n, p) != 1) 
    { 
        System.out.prletln("a and p are not coprime"); 
        return -1; 
    } 
  
    // If below expression return (p - 1) then modular 
    // square root is not possible 
    if (pow1(n, (p - 1) / 2, p) == (p - 1)) 
    { 
        document.write("no sqrt possible" + "<br/>"); 
        return -1; 
    } 
  
    // expressing p - 1, in terms of  
    // s * 2^e, where s is odd number 
    let s, e; 
    s = convertx2e(p - 1);
    e = z;
  
    // finding smallest q such that q ^ ((p - 1) / 2) 
    // (mod p) = p - 1 
    let q; 
    for (q = 2; ; q++) 
    { 
        // q - 1 is in place of (-1 % p) 
        if (pow1(q, (p - 1) / 2, p) == (p - 1)) 
            break; 
    } 
  
    // Initializing variable x, b and g 
    let x = pow1(n, (s + 1) / 2, p); 
    let b = pow1(n, s, p); 
    let g = pow1(q, s, p); 
  
    let r = e; 
  
    // keep looping until b 
    // become 1 or m becomes 0 
    while (true) 
    { 
        let m; 
        for (m = 0; m < r; m++) 
        { 
            if (order(p, b) == -1) 
                return -1; 
  
            // finding m such that b^ (2^m) = 1 
            if (order(p, b) == Math.pow(2, m)) 
                break; 
        } 
        if (m == 0) 
            return x; 
  
        // updating value of x, g and b
        // according to algorithm 
        x = (x * pow1(g, Math.pow(2, 
                            r - m - 1), p)) % p; 
        g = pow1(g, Math.pow(2, r - m), p); 
        b = (b * g) % p; 
  
        if (b == 1) 
            return x; 
        r = m; 
    } 
}

// Driver Code

     let n = 2; 
  
    // p should be prime 
    let p = 113; 
  
    let x = STonelli(n, p); 
  
    if (x == -1)
        document.write("Modular square" + 
                        "root is not exist\n"); 
    else
        document.write("Modular square root of " +
                            n + " and " + p + " is " +
                            x + "\n"); 

</script>

Output
Modular square root of 2 and 113 is 62

For more detail about above algorithm please visit : 
http://cs.indstate.edu/~sgali1/Shanks_Tonelli.pdf
For detail of example (2, 113) see : 
http://www.math.vt.edu/people/brown/class_homepages/shanks_tonelli.pdf
 


Next Article

Similar Reads