github.com/aergoio/aergo@v1.3.1/libtool/src/gmp-6.1.2/tests/mpn/t-bdiv.c (about)

     1  /* Copyright 2006, 2007, 2009, 2010 Free Software Foundation, Inc.
     2  
     3  This file is part of the GNU MP Library test suite.
     4  
     5  The GNU MP Library test suite is free software; you can redistribute it
     6  and/or modify it under the terms of the GNU General Public License as
     7  published by the Free Software Foundation; either version 3 of the License,
     8  or (at your option) any later version.
     9  
    10  The GNU MP Library test suite is distributed in the hope that it will be
    11  useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
    12  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General
    13  Public License for more details.
    14  
    15  You should have received a copy of the GNU General Public License along with
    16  the GNU MP Library test suite.  If not, see https://www.gnu.org/licenses/.  */
    17  
    18  
    19  #include <stdlib.h>		/* for strtol */
    20  #include <stdio.h>		/* for printf */
    21  
    22  #include "gmp.h"
    23  #include "gmp-impl.h"
    24  #include "longlong.h"
    25  #include "tests/tests.h"
    26  
    27  
    28  static void
    29  dumpy (mp_srcptr p, mp_size_t n)
    30  {
    31    mp_size_t i;
    32    if (n > 20)
    33      {
    34        for (i = n - 1; i >= n - 4; i--)
    35  	{
    36  	  printf ("%0*lx", (int) (2 * sizeof (mp_limb_t)), p[i]);
    37  	  printf (" ");
    38  	}
    39        printf ("... ");
    40        for (i = 3; i >= 0; i--)
    41  	{
    42  	  printf ("%0*lx", (int) (2 * sizeof (mp_limb_t)), p[i]);
    43  	  printf (i == 0 ? "" : " ");
    44  	}
    45      }
    46    else
    47      {
    48        for (i = n - 1; i >= 0; i--)
    49  	{
    50  	  printf ("%0*lx", (int) (2 * sizeof (mp_limb_t)), p[i]);
    51  	  printf (i == 0 ? "" : " ");
    52  	}
    53      }
    54    puts ("");
    55  }
    56  
    57  static unsigned long test;
    58  
    59  void
    60  check_one (mp_ptr qp, mp_srcptr rp, mp_limb_t rh,
    61  	   mp_srcptr np, mp_size_t nn, mp_srcptr dp, mp_size_t dn, const char *fname)
    62  {
    63    mp_size_t qn;
    64    int cmp;
    65    mp_ptr tp;
    66    mp_limb_t cy = 4711;		/* silence warnings */
    67    TMP_DECL;
    68  
    69    qn = nn - dn;
    70  
    71    if (qn == 0)
    72      return;
    73  
    74    TMP_MARK;
    75  
    76    tp = TMP_ALLOC_LIMBS (nn + 1);
    77  
    78    if (dn >= qn)
    79      mpn_mul (tp, dp, dn, qp, qn);
    80    else
    81      mpn_mul (tp, qp, qn, dp, dn);
    82  
    83    if (rp != NULL)
    84      {
    85        cy = mpn_add_n (tp + qn, tp + qn, rp, dn);
    86        cmp = cy != rh || mpn_cmp (tp, np, nn) != 0;
    87      }
    88    else
    89      cmp = mpn_cmp (tp, np, nn - dn) != 0;
    90  
    91    if (cmp != 0)
    92      {
    93        printf ("\r*******************************************************************************\n");
    94        printf ("%s inconsistent in test %lu\n", fname, test);
    95        printf ("N=   "); dumpy (np, nn);
    96        printf ("D=   "); dumpy (dp, dn);
    97        printf ("Q=   "); dumpy (qp, qn);
    98        if (rp != NULL)
    99  	{
   100  	  printf ("R=   "); dumpy (rp, dn);
   101  	  printf ("Rb=  %d, Cy=%d\n", (int) cy, (int) rh);
   102  	}
   103        printf ("T=   "); dumpy (tp, nn);
   104        printf ("nn = %ld, dn = %ld, qn = %ld", nn, dn, qn);
   105        printf ("\n*******************************************************************************\n");
   106        abort ();
   107      }
   108  
   109    TMP_FREE;
   110  }
   111  
   112  
   113  /* These are *bit* sizes. */
   114  #define SIZE_LOG 16
   115  #define MAX_DN (1L << SIZE_LOG)
   116  #define MAX_NN (1L << (SIZE_LOG + 1))
   117  
   118  #define COUNT 500
   119  
   120  mp_limb_t
   121  random_word (gmp_randstate_ptr rs)
   122  {
   123    mpz_t x;
   124    mp_limb_t r;
   125    TMP_DECL;
   126    TMP_MARK;
   127  
   128    MPZ_TMP_INIT (x, 2);
   129    mpz_urandomb (x, rs, 32);
   130    r = mpz_get_ui (x);
   131    TMP_FREE;
   132    return r;
   133  }
   134  
   135  int
   136  main (int argc, char **argv)
   137  {
   138    gmp_randstate_ptr rands;
   139    unsigned long maxnbits, maxdbits, nbits, dbits;
   140    mpz_t n, d, tz;
   141    mp_size_t maxnn, maxdn, nn, dn, clearn, i;
   142    mp_ptr np, dp, qp, rp;
   143    mp_limb_t rh;
   144    mp_limb_t t;
   145    mp_limb_t dinv;
   146    int count = COUNT;
   147    mp_ptr scratch;
   148    mp_limb_t ran;
   149    mp_size_t alloc, itch;
   150    mp_limb_t rran0, rran1, qran0, qran1;
   151    TMP_DECL;
   152  
   153    if (argc > 1)
   154      {
   155        char *end;
   156        count = strtol (argv[1], &end, 0);
   157        if (*end || count <= 0)
   158  	{
   159  	  fprintf (stderr, "Invalid test count: %s.\n", argv[1]);
   160  	  return 1;
   161  	}
   162      }
   163  
   164  
   165    maxdbits = MAX_DN;
   166    maxnbits = MAX_NN;
   167  
   168    tests_start ();
   169    rands = RANDS;
   170  
   171    mpz_init (n);
   172    mpz_init (d);
   173    mpz_init (tz);
   174  
   175    maxnn = maxnbits / GMP_NUMB_BITS + 1;
   176    maxdn = maxdbits / GMP_NUMB_BITS + 1;
   177  
   178    TMP_MARK;
   179  
   180    qp = TMP_ALLOC_LIMBS (maxnn + 2) + 1;
   181    rp = TMP_ALLOC_LIMBS (maxnn + 2) + 1;
   182  
   183    alloc = 1;
   184    scratch = __GMP_ALLOCATE_FUNC_LIMBS (alloc);
   185  
   186    for (test = 0; test < count;)
   187      {
   188        nbits = random_word (rands) % (maxnbits - GMP_NUMB_BITS) + 2 * GMP_NUMB_BITS;
   189        if (maxdbits > nbits)
   190  	dbits = random_word (rands) % nbits + 1;
   191        else
   192  	dbits = random_word (rands) % maxdbits + 1;
   193  
   194  #if RAND_UNIFORM
   195  #define RANDFUNC mpz_urandomb
   196  #else
   197  #define RANDFUNC mpz_rrandomb
   198  #endif
   199  
   200        do
   201  	{
   202  	  RANDFUNC (n, rands, nbits);
   203  	  do
   204  	    {
   205  	      RANDFUNC (d, rands, dbits);
   206  	    }
   207  	  while (mpz_sgn (d) == 0);
   208  
   209  	  np = PTR (n);
   210  	  dp = PTR (d);
   211  	  nn = SIZ (n);
   212  	  dn = SIZ (d);
   213  	}
   214        while (nn < dn);
   215  
   216        dp[0] |= 1;
   217  
   218        mpz_urandomb (tz, rands, 32);
   219        t = mpz_get_ui (tz);
   220  
   221        if (t % 17 == 0)
   222  	dp[0] = GMP_NUMB_MAX;
   223  
   224        switch ((int) t % 16)
   225  	{
   226  	case 0:
   227  	  clearn = random_word (rands) % nn;
   228  	  for (i = 0; i <= clearn; i++)
   229  	    np[i] = 0;
   230  	  break;
   231  	case 1:
   232  	  mpn_sub_1 (np + nn - dn, dp, dn, random_word (rands));
   233  	  break;
   234  	case 2:
   235  	  mpn_add_1 (np + nn - dn, dp, dn, random_word (rands));
   236  	  break;
   237  	}
   238  
   239        test++;
   240  
   241        binvert_limb (dinv, dp[0]);
   242  
   243        rran0 = random_word (rands);
   244        rran1 = random_word (rands);
   245        qran0 = random_word (rands);
   246        qran1 = random_word (rands);
   247  
   248        qp[-1] = qran0;
   249        qp[nn - dn + 1] = qran1;
   250        rp[-1] = rran0;
   251  
   252        ran = random_word (rands);
   253  
   254        if ((double) (nn - dn) * dn < 1e5)
   255  	{
   256  	  if (nn > dn)
   257  	    {
   258  	      /* Test mpn_sbpi1_bdiv_qr */
   259  	      MPN_ZERO (qp, nn - dn);
   260  	      MPN_ZERO (rp, dn);
   261  	      MPN_COPY (rp, np, nn);
   262  	      rh = mpn_sbpi1_bdiv_qr (qp, rp, nn, dp, dn, -dinv);
   263  	      ASSERT_ALWAYS (qp[-1] == qran0);  ASSERT_ALWAYS (qp[nn - dn + 1] == qran1);
   264  	      ASSERT_ALWAYS (rp[-1] == rran0);
   265  	      check_one (qp, rp + nn - dn, rh, np, nn, dp, dn, "mpn_sbpi1_bdiv_qr");
   266  	    }
   267  
   268  	  if (nn > dn)
   269  	    {
   270  	      /* Test mpn_sbpi1_bdiv_q */
   271  	      MPN_COPY (rp, np, nn);
   272  	      MPN_ZERO (qp, nn - dn);
   273  	      mpn_sbpi1_bdiv_q (qp, rp, nn - dn, dp, MIN(dn,nn-dn), -dinv);
   274  	      ASSERT_ALWAYS (qp[-1] == qran0);  ASSERT_ALWAYS (qp[nn - dn + 1] == qran1);
   275  	      ASSERT_ALWAYS (rp[-1] == rran0);
   276  	      check_one (qp, NULL, 0, np, nn, dp, dn, "mpn_sbpi1_bdiv_q");
   277  	    }
   278  	}
   279  
   280        if (dn >= 4 && nn - dn >= 2)
   281  	{
   282  	  /* Test mpn_dcpi1_bdiv_qr */
   283  	  MPN_COPY (rp, np, nn);
   284  	  MPN_ZERO (qp, nn - dn);
   285  	  rh = mpn_dcpi1_bdiv_qr (qp, rp, nn, dp, dn, -dinv);
   286  	  ASSERT_ALWAYS (qp[-1] == qran0);  ASSERT_ALWAYS (qp[nn - dn + 1] == qran1);
   287  	  ASSERT_ALWAYS (rp[-1] == rran0);
   288  	  check_one (qp, rp + nn - dn, rh, np, nn, dp, dn, "mpn_dcpi1_bdiv_qr");
   289  	}
   290  
   291        if (dn >= 4 && nn - dn >= 2)
   292  	{
   293  	  /* Test mpn_dcpi1_bdiv_q */
   294  	  MPN_COPY (rp, np, nn);
   295  	  MPN_ZERO (qp, nn - dn);
   296  	  mpn_dcpi1_bdiv_q (qp, rp, nn - dn, dp, MIN(dn,nn-dn), -dinv);
   297  	  ASSERT_ALWAYS (qp[-1] == qran0);  ASSERT_ALWAYS (qp[nn - dn + 1] == qran1);
   298  	  ASSERT_ALWAYS (rp[-1] == rran0);
   299  	  check_one (qp, NULL, 0, np, nn, dp, dn, "mpn_dcpi1_bdiv_q");
   300  	}
   301  
   302        if (nn > dn)
   303  	{
   304  	  /* Test mpn_bdiv_qr */
   305  	  itch = mpn_bdiv_qr_itch (nn, dn);
   306  	  if (itch + 1 > alloc)
   307  	    {
   308  	      scratch = __GMP_REALLOCATE_FUNC_LIMBS (scratch, alloc, itch + 1);
   309  	      alloc = itch + 1;
   310  	    }
   311  	  scratch[itch] = ran;
   312  	  MPN_ZERO (qp, nn - dn);
   313  	  MPN_ZERO (rp, dn);
   314  	  rp[dn] = rran1;
   315  	  rh = mpn_bdiv_qr (qp, rp, np, nn, dp, dn, scratch);
   316  	  ASSERT_ALWAYS (ran == scratch[itch]);
   317  	  ASSERT_ALWAYS (qp[-1] == qran0);  ASSERT_ALWAYS (qp[nn - dn + 1] == qran1);
   318  	  ASSERT_ALWAYS (rp[-1] == rran0);  ASSERT_ALWAYS (rp[dn] == rran1);
   319  
   320  	  check_one (qp, rp, rh, np, nn, dp, dn, "mpn_bdiv_qr");
   321  	}
   322  
   323        if (nn - dn < 2 || dn < 2)
   324  	continue;
   325  
   326        /* Test mpn_mu_bdiv_qr */
   327        itch = mpn_mu_bdiv_qr_itch (nn, dn);
   328        if (itch + 1 > alloc)
   329  	{
   330  	  scratch = __GMP_REALLOCATE_FUNC_LIMBS (scratch, alloc, itch + 1);
   331  	  alloc = itch + 1;
   332  	}
   333        scratch[itch] = ran;
   334        MPN_ZERO (qp, nn - dn);
   335        MPN_ZERO (rp, dn);
   336        rp[dn] = rran1;
   337        rh = mpn_mu_bdiv_qr (qp, rp, np, nn, dp, dn, scratch);
   338        ASSERT_ALWAYS (ran == scratch[itch]);
   339        ASSERT_ALWAYS (qp[-1] == qran0);  ASSERT_ALWAYS (qp[nn - dn + 1] == qran1);
   340        ASSERT_ALWAYS (rp[-1] == rran0);  ASSERT_ALWAYS (rp[dn] == rran1);
   341        check_one (qp, rp, rh, np, nn, dp, dn, "mpn_mu_bdiv_qr");
   342  
   343        /* Test mpn_mu_bdiv_q */
   344        itch = mpn_mu_bdiv_q_itch (nn, dn);
   345        if (itch + 1 > alloc)
   346  	{
   347  	  scratch = __GMP_REALLOCATE_FUNC_LIMBS (scratch, alloc, itch + 1);
   348  	  alloc = itch + 1;
   349  	}
   350        scratch[itch] = ran;
   351        MPN_ZERO (qp, nn - dn + 1);
   352        mpn_mu_bdiv_q (qp, np, nn - dn, dp, dn, scratch);
   353        ASSERT_ALWAYS (ran == scratch[itch]);
   354        ASSERT_ALWAYS (qp[-1] == qran0);  ASSERT_ALWAYS (qp[nn - dn + 1] == qran1);
   355        check_one (qp, NULL, 0, np, nn, dp, dn, "mpn_mu_bdiv_q");
   356      }
   357  
   358    __GMP_FREE_FUNC_LIMBS (scratch, alloc);
   359  
   360    TMP_FREE;
   361  
   362    mpz_clear (n);
   363    mpz_clear (d);
   364    mpz_clear (tz);
   365  
   366    tests_end ();
   367    return 0;
   368  }