1  /* Test mpf_sqrt_ui.
       2  
       3  Copyright 2004, 2015 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 <stdio.h>
      21  #include <stdlib.h>
      22  
      23  #include "gmp-impl.h"
      24  #include "tests.h"
      25  
      26  
      27  void
      28  check_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  
     115  int
     116  main (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  }