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

     1  /*
     2  
     3  Copyright 2012, Free Software Foundation, Inc.
     4  
     5  This file is part of the GNU MP Library test suite.
     6  
     7  The GNU MP Library test suite is free software; you can redistribute it
     8  and/or modify it under the terms of the GNU General Public License as
     9  published by the Free Software Foundation; either version 3 of the License,
    10  or (at your option) any later version.
    11  
    12  The GNU MP Library test suite is distributed in the hope that it will be
    13  useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
    14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General
    15  Public License for more details.
    16  
    17  You should have received a copy of the GNU General Public License along with
    18  the GNU MP Library test suite.  If not, see https://www.gnu.org/licenses/.  */
    19  
    20  #include <limits.h>
    21  #include <stdlib.h>
    22  #include <stdio.h>
    23  
    24  #include "testutils.h"
    25  
    26  #define MAXBITS 400
    27  #define COUNT 10000
    28  
    29  /* Called when g is supposed to be gcd(a,b), and g = s a + t b. */
    30  static int
    31  gcdext_valid_p (const mpz_t a, const mpz_t b,
    32  		const mpz_t g, const mpz_t s, const mpz_t t)
    33  {
    34    mpz_t ta, tb, r;
    35  
    36    /* It's not clear that gcd(0,0) is well defined, but we allow it and
    37       require that gcd(0,0) = 0. */
    38    if (mpz_sgn (g) < 0)
    39      return 0;
    40  
    41    if (mpz_sgn (a) == 0)
    42      {
    43        /* Must have g == abs (b). Any value for s is in some sense "correct",
    44  	 but it makes sense to require that s == 0. */
    45        return mpz_cmpabs (g, b) == 0 && mpz_sgn (s) == 0;
    46      }
    47    else if (mpz_sgn (b) == 0)
    48      {
    49        /* Must have g == abs (a), s == sign (a) */
    50        return mpz_cmpabs (g, a) == 0 && mpz_cmp_si (s, mpz_sgn (a)) == 0;
    51      }
    52  
    53    if (mpz_sgn (g) <= 0)
    54      return 0;
    55  
    56    mpz_init (ta);
    57    mpz_init (tb);
    58    mpz_init (r);
    59  
    60    mpz_mul (ta, s, a);
    61    mpz_mul (tb, t, b);
    62    mpz_add (ta, ta, tb);
    63  
    64    if (mpz_cmp (ta, g) != 0)
    65      {
    66      fail:
    67        mpz_clear (ta);
    68        mpz_clear (tb);
    69        mpz_clear (r);
    70        return 0;
    71      }
    72    mpz_tdiv_qr (ta, r, a, g);
    73    if (mpz_sgn (r) != 0)
    74      goto fail;
    75  
    76    mpz_tdiv_qr (tb, r, b, g);
    77    if (mpz_sgn (r) != 0)
    78      goto fail;
    79  
    80    /* Require that 2 |s| < |b/g|, or |s| == 1. */
    81    if (mpz_cmpabs_ui (s, 1) > 0)
    82      {
    83        mpz_mul_2exp (r, s, 1);
    84        if (mpz_cmpabs (r, tb) > 0)
    85  	goto fail;
    86      }
    87  
    88    /* Require that 2 |t| < |a/g| or |t| == 1*/
    89    if (mpz_cmpabs_ui (t, 1) > 0)
    90      {
    91        mpz_mul_2exp (r, t, 1);
    92        if (mpz_cmpabs (r, ta) > 0)
    93  	return 0;
    94      }
    95  
    96    mpz_clear (ta);
    97    mpz_clear (tb);
    98    mpz_clear (r);
    99  
   100    return 1;
   101  }
   102  
   103  void
   104  testmain (int argc, char **argv)
   105  {
   106    unsigned i;
   107    mpz_t a, b, g, s, t;
   108  
   109    mpz_init (a);
   110    mpz_init (b);
   111    mpz_init (g);
   112    mpz_init (s);
   113    mpz_init (t);
   114  
   115    for (i = 0; i < COUNT; i++)
   116      {
   117        mini_random_op3 (OP_GCD, MAXBITS, a, b, s);
   118        mpz_gcd (g, a, b);
   119        if (mpz_cmp (g, s))
   120  	{
   121  	  fprintf (stderr, "mpz_gcd failed:\n");
   122  	  dump ("a", a);
   123  	  dump ("b", b);
   124  	  dump ("r", g);
   125  	  dump ("ref", s);
   126  	  abort ();
   127  	}
   128      }
   129  
   130    for (i = 0; i < COUNT; i++)
   131      {
   132        unsigned flags;
   133        mini_urandomb (a, 32);
   134        flags = mpz_get_ui (a);
   135        mini_rrandomb (a, MAXBITS);
   136        mini_rrandomb (b, MAXBITS);
   137  
   138        if (flags % 37 == 0)
   139  	mpz_mul (a, a, b);
   140        if (flags % 37 == 1)
   141  	mpz_mul (b, a, b);
   142  
   143        if (flags & 1)
   144  	mpz_neg (a, a);
   145        if (flags & 2)
   146  	mpz_neg (b, b);
   147  
   148        mpz_gcdext (g, s, t, a, b);
   149        if (!gcdext_valid_p (a, b, g, s, t))
   150  	{
   151  	  fprintf (stderr, "mpz_gcdext failed:\n");
   152  	  dump ("a", a);
   153  	  dump ("b", b);
   154  	  dump ("g", g);
   155  	  dump ("s", s);
   156  	  dump ("t", t);
   157  	  abort ();
   158  	}
   159  
   160        mpz_gcd (s, a, b);
   161        if (mpz_cmp (g, s))
   162  	{
   163  	  fprintf (stderr, "mpz_gcd failed:\n");
   164  	  dump ("a", a);
   165  	  dump ("b", b);
   166  	  dump ("r", g);
   167  	  dump ("ref", s);
   168  	  abort ();
   169  	}
   170      }
   171    mpz_clear (a);
   172    mpz_clear (b);
   173    mpz_clear (g);
   174    mpz_clear (s);
   175    mpz_clear (t);
   176  }