busybox/libbb/isqrt.c
<<
>>
Prefs
   1/*
   2 * Copyright (C) 2017 Denys Vlasenko <vda.linux@googlemail.com>
   3 *
   4 * Licensed under GPLv2, see file LICENSE in this source tree.
   5 */
   6
   7//kbuild:lib-y += isqrt.o
   8
   9#ifndef ISQRT_TEST
  10# include "libbb.h"
  11#else
  12/* gcc -DISQRT_TEST -Wall -O2 isqrt.c -oisqrt && ./isqrt $((RANDOM*RANDOM)) */
  13# include <stdlib.h>
  14# include <stdio.h>
  15# include <time.h>
  16# define FAST_FUNC /* nothing */
  17#endif
  18
  19/* Returns such x that x+1 > sqrt(N) */
  20unsigned long FAST_FUNC isqrt(unsigned long long N)
  21{
  22        unsigned long x;
  23        unsigned shift;
  24#define LL_WIDTH_BITS (unsigned)(sizeof(N)*8)
  25
  26        shift = LL_WIDTH_BITS - 2;
  27        x = 0;
  28        do {
  29                x = (x << 1) + 1;
  30                if ((unsigned long long)x * x > (N >> shift))
  31                        x--; /* whoops, that +1 was too much */
  32                shift -= 2;
  33        } while ((int)shift >= 0);
  34        return x;
  35}
  36
  37#ifdef ISQRT_TEST
  38int main(int argc, char **argv)
  39{
  40        unsigned long long n = argv[1] ? strtoull(argv[1], NULL, 0) : time(NULL);
  41        for (;;) {
  42                unsigned long h;
  43                n--;
  44                h = isqrt(n);
  45                if (!(n & 0xffff))
  46                        printf("isqrt(%llx)=%lx\n", n, h);
  47                if ((unsigned long long)h * h > n) {
  48                        printf("BAD1: isqrt(%llx)=%lx\n", n, h);
  49                        return 1;
  50                }
  51                h++;
  52                if ((unsigned long long)h * h != 0 /* this can overflow to 0 - not a bug */
  53                 && (unsigned long long)h * h <= n)
  54                {
  55                        printf("BAD2: isqrt(%llx)=%lx\n", n, h);
  56                        return 1;
  57                }
  58        }
  59}
  60#endif
  61