cheng/mpir
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

added more tests for roots , corrected existing rootrem error

Browse Source
This commit is contained in:
jasonmoxham 2009-09-11 13:49:34 +00:00
parent 8063a5b84d
commit fff313847b
2 changed files with 124 additions and 11 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
mpz/rootrem.c
View File

@ -81,4 +81,5 @@ mpz_rootrem (mpz_ptr root, mpz_ptr rem, mpz_srcptr u, unsigned long int nth)
__GMP_FREE_FUNC_LIMBS (rootp, rootn);
SIZ(rem) = remn;
if(us<0)SIZ(rem)=-SIZ(rem);
}

134
tests/mpz/t-root.c
View File

@ -28,6 +28,19 @@ MA 02110-1301, USA. */
void debug_mp _PROTO ((mpz_t, int));
void
mpz_add_si (mpz_ptr x, mpz_srcptr y, signed long z)
{
if (z >= 0)
{
mpz_add_ui (x, y, z);
return;
}
mpz_sub_ui (x, y, -z);
return;
}
int
main (int argc, char **argv)
{
@ -35,7 +48,7 @@ main (int argc, char **argv)
mpz_t root1, root2, root3, rem2;
mpz_t temp, temp2;
mp_size_t x2_size;
int i;
int i, r1, j, k;
int reps = 1000;
unsigned long nth;
gmp_randstate_ptr rands;
@ -48,7 +61,7 @@ main (int argc, char **argv)
mpz_init (bs);
if (argc == 2)
reps = atoi (argv[1]);
reps = atoi (argv[1]);
mpz_init (x2);
mpz_init (root1);
@ -71,9 +84,13 @@ main (int argc, char **argv)
nth = mpz_getlimbn (bs, 0) % mpz_sizeinbase (x2, 2) + 2;
mpz_root (root1, x2, nth);
mpz_nthroot(root3,x2,nth);
if(mpz_cmp(root1,root3)!=0){printf("nthroot to root mismatch\n");abort();}
mpz_nthroot (root3, x2, nth);
if (mpz_cmp (root1, root3) != 0)
{
printf ("nthroot to root mismatch\n");
abort ();
}
mpz_urandomb (bs, rands, 4);
bsi = mpz_get_ui (bs);
if ((bsi & 1) != 0)
@ -88,8 +105,12 @@ main (int argc, char **argv)
else
mpz_add_ui (x2, x2, bsi >> 2);
mpz_root (root1, x2, nth);
mpz_nthroot(root3,x2,nth);
if(mpz_cmp(root1,root3)!=0){printf("nthroot to root mismatch\n");abort();}
mpz_nthroot (root3, x2, nth);
if (mpz_cmp (root1, root3) != 0)
{
printf ("nthroot to root mismatch\n");
abort ();
}
}
/* printf ("%ld %lu\n", SIZ (x2), nth); */
@ -99,7 +120,8 @@ main (int argc, char **argv)
mpz_add (temp2, temp, rem2);
/* Is power of result > argument? */
if (mpz_cmp (root1, root2) != 0 || mpz_cmp (x2, temp2) != 0 || mpz_cmp (temp, x2) > 0)
if (mpz_cmp (root1, root2) != 0 || mpz_cmp (x2, temp2) != 0
|| mpz_cmp (temp, x2) > 0)
{
fprintf (stderr, "ERROR after test %d\n", i);
debug_mp (x2, 10);
@ -109,7 +131,7 @@ main (int argc, char **argv)
abort ();
}
if (nth > 1 && mpz_cmp_ui (temp, 1L) > 0 && ! mpz_perfect_power_p (temp))
if (nth > 1 && mpz_cmp_ui (temp, 1L) > 0 && !mpz_perfect_power_p (temp))
{
fprintf (stderr, "ERROR in mpz_perfect_power_p after test %d\n", i);
debug_mp (temp, 10);
@ -134,7 +156,96 @@ main (int argc, char **argv)
abort ();
}
}
for (j = 0; j < reps; j++)
{ // pick a number bs and a power nth
mpz_urandomb (bs, rands, GMP_LIMB_BITS);
nth = mpz_get_ui (bs) % 30 + 1;
mpz_urandomb (bs, rands, GMP_LIMB_BITS);
x2_size = mpz_get_ui (bs) % 500;
mpz_urandomb (bs, rands, x2_size);
mpz_pow_ui (temp, bs, nth);
for (k = 0; k < 2; k++)
{
if (k != 0)
mpz_neg (temp, temp);
for (i = -10; i < 10; i++)
{
mpz_add_si (temp2, temp, i);
if (nth % 2 == 0 && mpz_cmp_ui (temp2, 0) < 0)
continue;
//printf("i is %d nth is %d\n",i,nth);
mpz_nthroot (root1, temp2, nth);
r1 = mpz_root (root2, temp2, nth);
mpz_rootrem (root3, rem2, temp2, nth);
if (mpz_cmp (root1, root2) != 0)
{
printf ("root12 mismatch\n");
abort ();
}
if (mpz_cmp (root2, root3) != 0)
{
printf ("root23 mismatch\n");
abort ();
}
mpz_pow_ui (root3, root3, nth);
mpz_add (root3, root3, rem2);
//printf("bs is ");mpz_out_str(stdout,10,bs);printf("\n");
//printf("temp2 is ");mpz_out_str(stdout,10,temp2);printf("\n");
//printf("root is ");mpz_out_str(stdout,10,root1);printf("\n");
//printf("rem is ");mpz_out_str(stdout,10,rem2);printf("\n");
if (mpz_cmp (root3, temp2) != 0)
{
printf ("rootrem mismatch\n");
abort ();
}
if (r1 && mpz_cmp_ui (rem2, 0) != 0)
{
printf ("rootexact1 mismatch %d\n", r1);
abort ();
}
if (!r1 && mpz_cmp_ui (rem2, 0) == 0)
{
printf ("rootexact2 mismatch\n");
abort ();
}
mpz_add_ui (root2, root2, 1);
mpz_pow_ui (root2, root2, nth);
if (mpz_cmp (root2, temp2) <= 0)
{
printf ("root wrong\n");
abort ();
}
if (!r1 && i == 0)
{
printf ("rootexact2 wrong\n");
abort ();
}
}
}
}
// check small values and stupid stuff
for (k = 0; k < 2; k++)
{
for (nth = 1; nth < 3 * reps; nth++)
{
if (k % 2 != 0 && nth % 2 == 0)
continue;
for (i = 0; i < 2 * reps; i++)
{
mpz_set_ui (temp, i);
if (k % 2 != 0)
mpz_neg (temp, temp);
mpz_rootrem (root1, rem2, temp, nth);
mpz_pow_ui (temp2, root1, nth);
mpz_add (temp2, temp2, rem2);
if (mpz_cmp (temp2, temp) != 0)
{
printf ("root small mismatch\n");
abort ();
}
}
}
}
mpz_clear (bs);
mpz_clear (x2);
mpz_clear (root1);
@ -151,5 +262,6 @@ main (int argc, char **argv)
void
debug_mp (mpz_t x, int base)
{
mpz_out_str (stderr, base, x); fputc ('\n', stderr);
mpz_out_str (stderr, base, x);
fputc ('\n', stderr);
}