1/* Test mpf_sqrt_ui. 2 3Copyright 2004, 2015 Free Software Foundation, Inc. 4 5This file is part of the GNU MP Library test suite. 6 7The GNU MP Library test suite is free software; you can redistribute it 8and/or modify it under the terms of the GNU General Public License as 9published by the Free Software Foundation; either version 3 of the License, 10or (at your option) any later version. 11 12The GNU MP Library test suite is distributed in the hope that it will be 13useful, but WITHOUT ANY WARRANTY; without even the implied warranty of 14MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 15Public License for more details. 16 17You should have received a copy of the GNU General Public License along with 18the GNU MP Library test suite. If not, see https://www.gnu.org/licenses/. */ 19 20#include <stdio.h> 21#include <stdlib.h> 22 23#include "gmp-impl.h" 24#include "tests.h" 25 26 27void 28check_rand (void) 29{ 30 unsigned long max_prec = 15; 31 unsigned long min_prec = __GMPF_BITS_TO_PREC (1); 32 gmp_randstate_ptr rands = RANDS; 33 unsigned long x, prec; 34 mpf_t r, s; 35 int i; 36 37 mpf_init (r); 38 mpf_init (s); 39 refmpf_set_prec_limbs (s, 2*max_prec+10); 40 41 for (x = 0; x < 2; x++) 42 { 43 mpf_sqrt_ui (r, x); 44 MPF_CHECK_FORMAT (r); 45 if (mpf_cmp_ui (r, x) != 0) 46 { 47 printf ("mpf_sqrt_ui wrong for special case:\n"); 48 printf (" x=%lu\n", x); 49 mpf_trace (" r", r); 50 abort (); 51 } 52 } 53 54 for (i = 0; i < 50; i++) 55 { 56 /* input, a random non-zero ulong, exponentially distributed */ 57 do { 58 x = gmp_urandomb_ui (rands, 59 gmp_urandomm_ui (rands, BITS_PER_ULONG) + 1); 60 } while (x <= 1); 61 62 /* result precision */ 63 prec = gmp_urandomm_ui (rands, max_prec-min_prec) + min_prec; 64 refmpf_set_prec_limbs (r, prec); 65 66 mpf_sqrt_ui (r, x); 67 MPF_CHECK_FORMAT (r); 68 69 /* Expect to prec limbs of result. 70 In the current implementation there's no stripping of low zero 71 limbs in mpf_sqrt_ui, not even on perfect squares, so size should 72 be exactly prec. */ 73 if (SIZ(r) != prec) 74 { 75 printf ("mpf_sqrt_ui result not enough result limbs\n"); 76 printf (" x=%lu\n", x); 77 printf (" want prec=%lu\n", prec); 78 mpf_trace (" r", r); 79 printf (" r size %ld\n", (long) SIZ(r)); 80 printf (" r prec %ld\n", (long) PREC(r)); 81 abort (); 82 } 83 84 /* Must have r^2 <= x, since r has been truncated. */ 85 mpf_mul (s, r, r); 86 if (! (mpf_cmp_ui (s, x) <= 0)) 87 { 88 printf ("mpf_sqrt_ui result too big\n"); 89 printf (" x=%lu\n", x); 90 printf (" want prec=%lu\n", prec); 91 mpf_trace (" r", r); 92 mpf_trace (" s", s); 93 abort (); 94 } 95 96 /* Must have (r+ulp)^2 > x. 97 No overflow from refmpf_add_ulp since r is only prec limbs. */ 98 refmpf_add_ulp (r); 99 mpf_mul (s, r, r); 100 if (! (mpf_cmp_ui (s, x) > 0)) 101 { 102 printf ("mpf_sqrt_ui result too small\n"); 103 printf (" x=%lu\n", x); 104 printf (" want prec=%lu\n", prec); 105 mpf_trace (" r+ulp", r); 106 mpf_trace (" s", s); 107 abort (); 108 } 109 } 110 111 mpf_clear (r); 112 mpf_clear (s); 113} 114 115int 116main (int argc, char **argv) 117{ 118 tests_start (); 119 mp_trace_base = -16; 120 121 check_rand (); 122 123 tests_end (); 124 exit (0); 125} 126