Le moyen le plus rapide d’obtenir la partie entière de sqrt (n)?

Question

Comme nous le savons si n n'est pas un carré parfait, alors sqrt(n) ne serait pas un nombre entier. Car j'ai besoin d'uniquement la partie entière, j'ai l'impression que l'appelant sqrt(n) ne serait pas rapide, il faut du temps pour calculer la fraction de la partie également.

Donc ma question est :

Pouvons-nous obtenir uniquement la partie entière de sqrt(n) sans calculer la valeur réelle de sqrt(n)? L'algorithme doit être plus rapide qu' sqrt(n) (définie en <math.h> ou <cmath>)?

Si possible, vous pouvez écrire le code en asm bloc également.

Answer 1

Je voudrais essayer le Rapide Inverse de la Racine Carrée de truc.

C'est un moyen d'obtenir une très bonne approximation de la 1/sqrt(n) sans aucune direction, basée sur quelques bits de se tourner afin de ne pas portable (notamment entre les 32-bits et 64-bits pour les plates-formes).

Une fois que vous l'obtenez, vous avez juste besoin d'inverser le résultat, et prend la partie entière.

Il y a peut-être plus rapide des astuces, des cours, puisque celui-ci est un peu ronde.

EDIT: let's do it!

D'abord une petite aide:

// benchmark.h
#include <sys/time.h>

template <typename Func>
double benchmark(Func f, size_t iterations)
{
  f();

  timeval a, b;
  gettimeofday(&a, 0);
  for (; iterations --> 0;)
  {
    f();
  }
  gettimeofday(&b, 0);
  return (b.tv_sec * (unsigned int)1e6 + b.tv_usec) -
         (a.tv_sec * (unsigned int)1e6 + a.tv_usec);
}

Puis le corps principal:

#include <iostream>

#include <cmath>

#include "benchmark.h"

class Sqrt
{
public:
  Sqrt(int n): _number(n) {}

  int operator()() const
  {
    double d = _number;
    return static_cast<int>(std::sqrt(d) + 0.5);
  }

private:
  int _number;
};

// http://www.codecodex.com/wiki/Calculate_an_integer_square_root
class IntSqrt
{
public:
  IntSqrt(int n): _number(n) {}

  int operator()() const 
  {
    int remainder = _number;
    if (remainder < 0) { return 0; }

    int place = 1 <<(sizeof(int)*8 -2);

    while (place > remainder) { place /= 4; }

    int root = 0;
    while (place)
    {
      if (remainder >= root + place)
      {
        remainder -= root + place;
        root += place*2;
      }
      root /= 2;
      place /= 4;
    }
    return root;
  }

private:
  int _number;
};

// http://en.wikipedia.org/wiki/Fast_inverse_square_root
class FastSqrt
{
public:
  FastSqrt(int n): _number(n) {}

  int operator()() const
  {
    float number = _number;

    float x2 = number * 0.5F;
    float y = number;
    long i = *(long*)&y;
    //i = (long)0x5fe6ec85e7de30da - (i >> 1);
    i = 0x5f3759df - (i >> 1);
    y = *(float*)&i;

    y = y * (1.5F - (x2*y*y));
    y = y * (1.5F - (x2*y*y)); // let's be precise

    return static_cast<int>(1/y + 0.5f);
  }

private:
  int _number;
};


int main(int argc, char* argv[])
{
  if (argc != 3) {
    std::cerr << "Usage: %prog integer iterations\n";
    return 1;
  }

  int n = atoi(argv[1]);
  int it = atoi(argv[2]);

  assert(Sqrt(n)() == IntSqrt(n)() &&
          Sqrt(n)() == FastSqrt(n)() && "Different Roots!");
  std::cout << "sqrt(" << n << ") = " << Sqrt(n)() << "\n";

  double time = benchmark(Sqrt(n), it);
  double intTime = benchmark(IntSqrt(n), it);
  double fastTime = benchmark(FastSqrt(n), it);

  std::cout << "Number iterations: " << it << "\n"
               "Sqrt computation : " << time << "\n"
               "Int computation  : " << intTime << "\n"
               "Fast computation : " << fastTime << "\n";

  return 0;
}

Et les résultats:

sqrt(82) = 9
Number iterations: 4096
Sqrt computation : 56
Int computation  : 217
Fast computation : 119

// Note had to tweak the program here as Int here returns -1 :/
sqrt(2147483647) = 46341 // real answer sqrt(2 147 483 647) = 46 340.95
Number iterations: 4096
Sqrt computation : 57
Int computation  : 313
Fast computation : 119

Où comme prévu, le Rapide calcul effectue beaucoup mieux que l' Int calcul.

Oh, et par la manière, sqrt est plus rapide :)

Answer 2

Edit: cette réponse est stupide utiliser (int) sqrt(i)

Après profilage avec les bon paramètres (-march=native -m64 -O3), le plus haut était un beaucoup plus rapide.

---

Bon, un peu vieille question, mais le "le plus rapide" réponse n'a pas encore été donné. La manière la plus rapide (je pense) est le Binaire de la Racine Carrée de l'algorithme, expliqué en détail dans ce Embedded.com l'article.

Il fondamentalement revient à ceci:

unsigned short isqrt(unsigned long a) {
    unsigned long rem = 0;
    int root = 0;
    int i;

    for (i = 0; i < 16; i++) {
        root <<= 1;
        rem <<= 2;
        rem += a >> 30;
        a <<= 2;

        if (root < rem) {
            root++;
            rem -= root;
            root++;
        }
    }

    return (unsigned short) (root >> 1);
}

Sur ma machine (Q6600, Ubuntu 10.10) je profilé en prenant la racine carrée des nombres 1-100000000. À l'aide de iqsrt(i) a pris 2750 mme. À l'aide de (unsigned short) sqrt((float) i) a pris 3600ms. Cela a été fait à l'aide de g++ -O3. À l'aide de l' -ffast-math option de compilation les temps étaient 2100ms et 3100ms respectivement. Remarque c'est sans même en utilisant une seule ligne d'assembleur, de sorte qu'il peut être encore plus rapide.

Le code ci-dessus fonctionne en C et en C++ et avec de légères modifications de syntaxe aussi pour Java.

Ce qui fonctionne encore mieux pour une gamme limitée est une binaire de recherche. Sur ma machine ce coups la version ci-dessus de l'eau par un facteur 4. C'est malheureusement très limitée dans la gamme:

#include <stdint.h>

const uint16_t squares[] = {
    0, 1, 4, 9,
    16, 25, 36, 49,
    64, 81, 100, 121,
    144, 169, 196, 225,
    256, 289, 324, 361,
    400, 441, 484, 529,
    576, 625, 676, 729,
    784, 841, 900, 961,
    1024, 1089, 1156, 1225,
    1296, 1369, 1444, 1521,
    1600, 1681, 1764, 1849,
    1936, 2025, 2116, 2209,
    2304, 2401, 2500, 2601,
    2704, 2809, 2916, 3025,
    3136, 3249, 3364, 3481,
    3600, 3721, 3844, 3969,
    4096, 4225, 4356, 4489,
    4624, 4761, 4900, 5041,
    5184, 5329, 5476, 5625,
    5776, 5929, 6084, 6241,
    6400, 6561, 6724, 6889,
    7056, 7225, 7396, 7569,
    7744, 7921, 8100, 8281,
    8464, 8649, 8836, 9025,
    9216, 9409, 9604, 9801,
    10000, 10201, 10404, 10609,
    10816, 11025, 11236, 11449,
    11664, 11881, 12100, 12321,
    12544, 12769, 12996, 13225,
    13456, 13689, 13924, 14161,
    14400, 14641, 14884, 15129,
    15376, 15625, 15876, 16129,
    16384, 16641, 16900, 17161,
    17424, 17689, 17956, 18225,
    18496, 18769, 19044, 19321,
    19600, 19881, 20164, 20449,
    20736, 21025, 21316, 21609,
    21904, 22201, 22500, 22801,
    23104, 23409, 23716, 24025,
    24336, 24649, 24964, 25281,
    25600, 25921, 26244, 26569,
    26896, 27225, 27556, 27889,
    28224, 28561, 28900, 29241,
    29584, 29929, 30276, 30625,
    30976, 31329, 31684, 32041,
    32400, 32761, 33124, 33489,
    33856, 34225, 34596, 34969,
    35344, 35721, 36100, 36481,
    36864, 37249, 37636, 38025,
    38416, 38809, 39204, 39601,
    40000, 40401, 40804, 41209,
    41616, 42025, 42436, 42849,
    43264, 43681, 44100, 44521,
    44944, 45369, 45796, 46225,
    46656, 47089, 47524, 47961,
    48400, 48841, 49284, 49729,
    50176, 50625, 51076, 51529,
    51984, 52441, 52900, 53361,
    53824, 54289, 54756, 55225,
    55696, 56169, 56644, 57121,
    57600, 58081, 58564, 59049,
    59536, 60025, 60516, 61009,
    61504, 62001, 62500, 63001,
    63504, 64009, 64516, 65025
};

inline int isqrt(uint16_t x) {
    const uint16_t *p = squares;

    if (p[128] <= x) p += 128;
    if (p[ 64] <= x) p +=  64;
    if (p[ 32] <= x) p +=  32;
    if (p[ 16] <= x) p +=  16;
    if (p[  8] <= x) p +=   8;
    if (p[  4] <= x) p +=   4;
    if (p[  2] <= x) p +=   2;
    if (p[  1] <= x) p +=   1;

    return p - squares;
}

Une version 32 bits peut être téléchargé ici: https://gist.github.com/3481770

Answer 3

Je pense que Google search fournit de bons articles comme Calculate an integer square root qui ont discuté aussi de nombreuses façons de calcul rapide et il y a de bons articles de référence, je pense qu'ici, personne ne peut fournir de meilleurs qu'eux (et si quelqu'un peut première sera de produire du papier à ce sujet), mais si vous les lisez et il y a ambiguïté avec eux, alors peut-être nous pouvons vous aider.

Answer 4

Bien que je soupçonne que vous pouvez trouver un beaucoup de les options de la recherche pour "fast entier racine carrée", voici quelques potentiellement de nouvelles idées qui pourraient bien fonctionner (chaque indépendant, ou peut-être vous pouvez combiner entre eux):

1. Faire un static const tableau de tous les carrés parfaits dans le nom de domaine que vous souhaitez soutenir, et d'effectuer une rapide sans branches binaire de recherche. L'indice dans le tableau est la racine carrée.
2. Convertir le nombre à virgule flottante et de le briser en mantisse et l'exposant. Réduire de moitié l'exposant et de multiplier la mantisse par quelque magie facteur (à vous de trouver). Ce devrait être en mesure de vous donner une très bonne approximation. Inclure une étape finale pour l'ajuster si ce n'est pas exact (ou l'utiliser comme un point de départ pour le binaire de recherche ci-dessus).

Answer 5

Pourquoi personne ne suggère la méthode la plus rapide?

Si:

1. la gamme de nombre est limité
2. la consommation de mémoire n'est pas crucial
3. lancement de l'application le temps n'est pas critique

puis créez int[MAX_X] rempli (au lancement) avec sqrt(x) (vous n'avez pas besoin d'utiliser la fonction sqrt() pour elle).

L'ensemble de ces conditions s'adapter à mon programme assez bien. En particulier, une int[10000000] tableau va consommer 40MB.

Quelles sont vos pensées sur cette question?