1217309Snwhitehorn/* Test mpn_mul function for all sizes up to a selected limit.
2217309Snwhitehorn
3217309SnwhitehornCopyright 2011, 2012 Free Software Foundation, Inc.
4217309Snwhitehorn
5217309SnwhitehornThis file is part of the GNU MP Library test suite.
6217309Snwhitehorn
7217309SnwhitehornThe GNU MP Library test suite is free software; you can redistribute it
8217309Snwhitehornand/or modify it under the terms of the GNU General Public License as
9217309Snwhitehornpublished by the Free Software Foundation; either version 3 of the License,
10217309Snwhitehornor (at your option) any later version.
11217309Snwhitehorn
12217309SnwhitehornThe GNU MP Library test suite is distributed in the hope that it will be
13217309Snwhitehornuseful, but WITHOUT ANY WARRANTY; without even the implied warranty of
14217309SnwhitehornMERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General
15217309SnwhitehornPublic License for more details.
16217309Snwhitehorn
17217309SnwhitehornYou should have received a copy of the GNU General Public License along with
18217309Snwhitehornthe GNU MP Library test suite.  If not, see https://www.gnu.org/licenses/.  */
19217309Snwhitehorn
20217309Snwhitehorn
21217309Snwhitehorn#include <stdlib.h>
22217309Snwhitehorn#include <stdio.h>
23217309Snwhitehorn
24217309Snwhitehorn#include "gmp-impl.h"
25217309Snwhitehorn#include "tests.h"
26217309Snwhitehorn
27217309Snwhitehornstatic unsigned
28217309Snwhitehornisqrt (unsigned t)
29{
30  unsigned s, b;
31
32  for (b = 0, s = t;  b++, s >>= 1; )
33    ;
34
35  s = 1 << (b >> 1);
36  if (b & 1)
37    s += s >> 1;
38
39  do
40    {
41      b = t / s;
42      s = (s + b) >> 1;
43    }
44  while (b < s);
45
46  return s;
47}
48
49int
50main (int argc, char **argv)
51{
52  mp_ptr ap, bp, rp, refp;
53  mp_size_t max_n, an, bn, rn;
54  int reps;
55  TMP_DECL;
56  TMP_MARK;
57
58  reps = 1;
59
60  tests_start ();
61  TESTS_REPS (reps, argv, argc);
62
63  /* Re-interpret reps argument as a size argument.  */
64  max_n = isqrt (reps * 25000);
65
66  ap = TMP_ALLOC_LIMBS (max_n + 1);
67  bp = TMP_ALLOC_LIMBS (max_n + 1);
68  rp = TMP_ALLOC_LIMBS (2 * max_n);
69  refp = TMP_ALLOC_LIMBS (2 * max_n);
70
71  for (an = 1; an <= max_n; an += 1)
72    {
73      for (bn = 1; bn <= an; bn += 1)
74	{
75	  mpn_random2 (ap, an + 1);
76	  mpn_random2 (bp, bn + 1);
77
78	  refmpn_mul (refp, ap, an, bp, bn);
79	  mpn_mul (rp, ap, an, bp, bn);
80
81	  rn = an + bn;
82	  if (mpn_cmp (refp, rp, rn))
83	    {
84	      printf ("ERROR, an = %d, bn = %d, rn = %d\n",
85		      (int) an, (int) bn, (int) rn);
86	      printf ("a: "); mpn_dump (ap, an);
87	      printf ("b: "); mpn_dump (bp, bn);
88	      printf ("r:   "); mpn_dump (rp, rn);
89	      printf ("ref: "); mpn_dump (refp, rn);
90	      abort();
91	    }
92	}
93    }
94  TMP_FREE;
95  tests_end ();
96  return 0;
97}
98