Line data    Source code

       1             : #include "tommath_private.h"
       2             : #ifdef BN_MP_KRONECKER_C
       3             : 
       4             : /* LibTomMath, multiple-precision integer library -- Tom St Denis */
       5             : /* SPDX-License-Identifier: Unlicense */
       6             : 
       7             : /*
       8             :    Kronecker symbol (a|p)
       9             :    Straightforward implementation of algorithm 1.4.10 in
      10             :    Henri Cohen: "A Course in Computational Algebraic Number Theory"
      11             : 
      12             :    @book{cohen2013course,
      13             :      title={A course in computational algebraic number theory},
      14             :      author={Cohen, Henri},
      15             :      volume={138},
      16             :      year={2013},
      17             :      publisher={Springer Science \& Business Media}
      18             :     }
      19             :  */
      20           0 : mp_err mp_kronecker(const mp_int *a, const mp_int *p, int *c)
      21             : {
      22           0 :    mp_int a1, p1, r;
      23           0 :    mp_err err;
      24           0 :    int v, k;
      25             : 
      26           0 :    static const int table[8] = {0, 1, 0, -1, 0, -1, 0, 1};
      27             : 
      28           0 :    if (MP_IS_ZERO(p)) {
      29           0 :       if ((a->used == 1) && (a->dp[0] == 1u)) {
      30           0 :          *c = 1;
      31             :       } else {
      32           0 :          *c = 0;
      33             :       }
      34           0 :       return MP_OKAY;
      35             :    }
      36             : 
      37           0 :    if (MP_IS_EVEN(a) && MP_IS_EVEN(p)) {
      38           0 :       *c = 0;
      39           0 :       return MP_OKAY;
      40             :    }
      41             : 
      42           0 :    if ((err = mp_init_copy(&a1, a)) != MP_OKAY) {
      43           0 :       return err;
      44             :    }
      45           0 :    if ((err = mp_init_copy(&p1, p)) != MP_OKAY) {
      46           0 :       goto LBL_KRON_0;
      47             :    }
      48             : 
      49           0 :    v = mp_cnt_lsb(&p1);
      50           0 :    if ((err = mp_div_2d(&p1, v, &p1, NULL)) != MP_OKAY) {
      51           0 :       goto LBL_KRON_1;
      52             :    }
      53             : 
      54           0 :    if ((v & 1) == 0) {
      55           0 :       k = 1;
      56             :    } else {
      57           0 :       k = table[a->dp[0] & 7u];
      58             :    }
      59             : 
      60           0 :    if (p1.sign == MP_NEG) {
      61           0 :       p1.sign = MP_ZPOS;
      62           0 :       if (a1.sign == MP_NEG) {
      63           0 :          k = -k;
      64             :       }
      65             :    }
      66             : 
      67           0 :    if ((err = mp_init(&r)) != MP_OKAY) {
      68           0 :       goto LBL_KRON_1;
      69             :    }
      70             : 
      71           0 :    for (;;) {
      72           0 :       if (MP_IS_ZERO(&a1)) {
      73           0 :          if (mp_cmp_d(&p1, 1uL) == MP_EQ) {
      74           0 :             *c = k;
      75           0 :             goto LBL_KRON;
      76             :          } else {
      77           0 :             *c = 0;
      78           0 :             goto LBL_KRON;
      79             :          }
      80             :       }
      81             : 
      82           0 :       v = mp_cnt_lsb(&a1);
      83           0 :       if ((err = mp_div_2d(&a1, v, &a1, NULL)) != MP_OKAY) {
      84           0 :          goto LBL_KRON;
      85             :       }
      86             : 
      87           0 :       if ((v & 1) == 1) {
      88           0 :          k = k * table[p1.dp[0] & 7u];
      89             :       }
      90             : 
      91           0 :       if (a1.sign == MP_NEG) {
      92             :          /*
      93             :           * Compute k = (-1)^((a1)*(p1-1)/4) * k
      94             :           * a1.dp[0] + 1 cannot overflow because the MSB
      95             :           * of the type mp_digit is not set by definition
      96             :           */
      97           0 :          if (((a1.dp[0] + 1u) & p1.dp[0] & 2u) != 0u) {
      98           0 :             k = -k;
      99             :          }
     100             :       } else {
     101             :          /* compute k = (-1)^((a1-1)*(p1-1)/4) * k */
     102           0 :          if ((a1.dp[0] & p1.dp[0] & 2u) != 0u) {
     103           0 :             k = -k;
     104             :          }
     105             :       }
     106             : 
     107           0 :       if ((err = mp_copy(&a1, &r)) != MP_OKAY) {
     108           0 :          goto LBL_KRON;
     109             :       }
     110           0 :       r.sign = MP_ZPOS;
     111           0 :       if ((err = mp_mod(&p1, &r, &a1)) != MP_OKAY) {
     112           0 :          goto LBL_KRON;
     113             :       }
     114           0 :       if ((err = mp_copy(&r, &p1)) != MP_OKAY) {
     115           0 :          goto LBL_KRON;
     116             :       }
     117             :    }
     118             : 
     119           0 : LBL_KRON:
     120           0 :    mp_clear(&r);
     121           0 : LBL_KRON_1:
     122           0 :    mp_clear(&p1);
     123           0 : LBL_KRON_0:
     124           0 :    mp_clear(&a1);
     125             : 
     126           0 :    return err;
     127             : }
     128             : 
     129             : #endif