(root)/
gmp-6.3.0/
mpn/
generic/
sec_pi1_div.c
       1  /* mpn_sec_pi1_div_qr, mpn_sec_pi1_div_r -- Compute Q = floor(U / V), U = U
       2     mod V.  Side-channel silent under the assumption that the used instructions
       3     are side-channel silent.
       4  
       5     Contributed to the GNU project by Torbjörn Granlund.
       6  
       7     THE FUNCTIONS IN THIS FILE ARE INTERNAL WITH MUTABLE INTERFACES.  IT IS ONLY
       8     SAFE TO REACH THEM THROUGH DOCUMENTED INTERFACES.  IN FACT, IT IS ALMOST
       9     GUARANTEED THAT THEY WILL CHANGE OR DISAPPEAR IN A FUTURE GNU MP RELEASE.
      10  
      11  Copyright 2011-2013 Free Software Foundation, Inc.
      12  
      13  This file is part of the GNU MP Library.
      14  
      15  The GNU MP Library is free software; you can redistribute it and/or modify
      16  it under the terms of either:
      17  
      18    * the GNU Lesser General Public License as published by the Free
      19      Software Foundation; either version 3 of the License, or (at your
      20      option) any later version.
      21  
      22  or
      23  
      24    * the GNU General Public License as published by the Free Software
      25      Foundation; either version 2 of the License, or (at your option) any
      26      later version.
      27  
      28  or both in parallel, as here.
      29  
      30  The GNU MP Library is distributed in the hope that it will be useful, but
      31  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
      32  or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
      33  for more details.
      34  
      35  You should have received copies of the GNU General Public License and the
      36  GNU Lesser General Public License along with the GNU MP Library.  If not,
      37  see https://www.gnu.org/licenses/.  */
      38  
      39  #include "gmp-impl.h"
      40  #include "longlong.h"
      41  
      42  /* This side-channel silent division algorithm reduces the partial remainder by
      43     GMP_NUMB_BITS/2 bits at a time, compared to GMP_NUMB_BITS for the main
      44     division algorithm.  We actually do not insist on reducing by exactly
      45     GMP_NUMB_BITS/2, but may leave a partial remainder that is D*B^i to 3D*B^i
      46     too large (B is the limb base, D is the divisor, and i is the induction
      47     variable); the subsequent step will handle the extra partial remainder bits.
      48  
      49     With that partial remainder reduction, each step generates a quotient "half
      50     limb".  The outer loop generates two quotient half limbs, an upper (q1h) and
      51     a lower (q0h) which are stored sparsely in separate limb arrays.  These
      52     arrays are added at the end; using separate arrays avoids data-dependent
      53     carry propagation which could else pose a side-channel leakage problem.
      54  
      55     The quotient half limbs may be between -3 to 0 from the accurate value
      56     ("accurate" being the one which corresponds to a reduction to a principal
      57     partial remainder).  Too small quotient half limbs correspond to too large
      58     remainders, which we reduce later, as described above.
      59  
      60     In order to keep quotients from getting too big, corresponding to a negative
      61     partial remainder, we use an inverse which is slightly smaller than usually.
      62  */
      63  
      64  #if OPERATION_sec_pi1_div_qr
      65  /* Needs (dn + 1) + (nn - dn) + (nn - dn) = 2nn - dn + 1 limbs at tp. */
      66  #define FNAME mpn_sec_pi1_div_qr
      67  #define Q(q) q,
      68  #define RETTYPE mp_limb_t
      69  #endif
      70  #if OPERATION_sec_pi1_div_r
      71  /* Needs (dn + 1) limbs at tp.  */
      72  #define FNAME mpn_sec_pi1_div_r
      73  #define Q(q)
      74  #define RETTYPE void
      75  #endif
      76  
      77  RETTYPE
      78  FNAME (Q(mp_ptr qp)
      79         mp_ptr np, mp_size_t nn,
      80         mp_srcptr dp, mp_size_t dn,
      81         mp_limb_t dinv,
      82         mp_ptr tp)
      83  {
      84    mp_limb_t nh, cy, q1h, q0h, dummy, cnd;
      85    mp_size_t i;
      86    mp_ptr hp;
      87  #if OPERATION_sec_pi1_div_qr
      88    mp_limb_t qh;
      89    mp_ptr qlp, qhp;
      90  #endif
      91  
      92    ASSERT (dn >= 1);
      93    ASSERT (nn >= dn);
      94    ASSERT ((dp[dn - 1] & GMP_NUMB_HIGHBIT) != 0);
      95  
      96    if (nn == dn)
      97      {
      98        cy = mpn_sub_n (np, np, dp, dn);
      99        mpn_cnd_add_n (cy, np, np, dp, dn);
     100  #if OPERATION_sec_pi1_div_qr
     101        return 1 - cy;
     102  #else
     103        return;
     104  #endif
     105      }
     106  
     107    /* Create a divisor copy shifted half a limb.  */
     108    hp = tp;					/* (dn + 1) limbs */
     109    hp[dn] = mpn_lshift (hp, dp, dn, GMP_NUMB_BITS / 2);
     110  
     111  #if OPERATION_sec_pi1_div_qr
     112    qlp = tp + (dn + 1);				/* (nn - dn) limbs */
     113    qhp = tp + (nn + 1);				/* (nn - dn) limbs */
     114  #endif
     115  
     116    np += nn - dn;
     117    nh = 0;
     118  
     119    for (i = nn - dn - 1; i >= 0; i--)
     120      {
     121        np--;
     122  
     123        nh = (nh << GMP_NUMB_BITS/2) + (np[dn] >> GMP_NUMB_BITS/2);
     124        umul_ppmm (q1h, dummy, nh, dinv);
     125        q1h += nh;
     126  #if OPERATION_sec_pi1_div_qr
     127        qhp[i] = q1h;
     128  #endif
     129        mpn_submul_1 (np, hp, dn + 1, q1h);
     130  
     131        nh = np[dn];
     132        umul_ppmm (q0h, dummy, nh, dinv);
     133        q0h += nh;
     134  #if OPERATION_sec_pi1_div_qr
     135        qlp[i] = q0h;
     136  #endif
     137        nh -= mpn_submul_1 (np, dp, dn, q0h);
     138      }
     139  
     140    /* 1st adjustment depends on extra high remainder limb.  */
     141    cnd = nh != 0;				/* FIXME: cmp-to-int */
     142  #if OPERATION_sec_pi1_div_qr
     143    qlp[0] += cnd;
     144  #endif
     145    nh -= mpn_cnd_sub_n (cnd, np, np, dp, dn);
     146  
     147    /* 2nd adjustment depends on remainder/divisor comparison as well as whether
     148       extra remainder limb was nullified by previous subtract.  */
     149    cy = mpn_sub_n (np, np, dp, dn);
     150    cy = cy - nh;
     151  #if OPERATION_sec_pi1_div_qr
     152    qlp[0] += 1 - cy;
     153  #endif
     154    mpn_cnd_add_n (cy, np, np, dp, dn);
     155  
     156    /* 3rd adjustment depends on remainder/divisor comparison.  */
     157    cy = mpn_sub_n (np, np, dp, dn);
     158  #if OPERATION_sec_pi1_div_qr
     159    qlp[0] += 1 - cy;
     160  #endif
     161    mpn_cnd_add_n (cy, np, np, dp, dn);
     162  
     163  #if OPERATION_sec_pi1_div_qr
     164    /* Combine quotient halves into final quotient.  */
     165    qh = mpn_lshift (qhp, qhp, nn - dn, GMP_NUMB_BITS/2);
     166    qh += mpn_add_n (qp, qhp, qlp, nn - dn);
     167  
     168    return qh;
     169  #else
     170    return;
     171  #endif
     172  }