sqrt af long long's?

Det kommer an på hvad du mener med bedre.
På moderne maskiner er det sansynligvis det hurtigste. Du skal dog være opmærksom på at der kan blive tabt præsision ved konvertering fra long long til double og tilbage igen.

Hvad mener du med "på moderne maskiner"? Jeg ville da tro at det tager flere clock cykels eller hvordan det nu staves at convertere det end hvis man havde en funtion (måske noget Asambler) til at lave det direkte på variabel?

Med moderne maskiner mener jeg f.ex en PC der er mindre end 5 år gammel.
Jeg ved ikke nok om assembler på PC'er til at vide om det kan gøres hurtige ved at operere direkte på long long en ved at konvertere til og fra double, men jeg tvivler.
Jeg er ret sikker på at man ikke kan lave en funktion der udregner sqrt på long long's i C eller C++ der er hurtigere end compilerens sqrt.

Nu kan man altid diskutere hvad bedre er.

Taler vi performance eller kode vedligeholdelses omkostninger.

Hvis du spørger om du selv kan kode det, så er svaret ja.

Her er et eksempel:

#include <stdio.h>
#include <math.h>

long long llsqrt(long long v)
{
  long long tmp,res;
  res = v / 2;
  do
  {
      tmp = res;
      res = (v / tmp + tmp) / 2;
  } while((res > tmp + 1) || (tmp > res + 1));
  tmp = res;
  res = (v / tmp + tmp) / 2;
  return (res < tmp ) ? res : tmp;
}

int main()
{
  int i;
  long long v;
  for(i = 2; i <= 20; i++) {
      v = i;
      printf("%d %ld\n",i,(long int)llsqrt(v));
  }
  v = 1234567890;
  printf("%ld %ld %ld\n",1234567890,(long int)sqrt((double)v),(long)llsqrt(v));
  return 0;
}

Er de nye typer med int_64 ikke direkte understøttet af matematik library'et?

segmose >> it_64 og long long's (som man burde bruge efter som de er de "orginale" 64 bit varibler) kan ikke bruges sammen med math.h funtionerne... men nogle gang convertere compileren dem til en double. -og tilbage igen.
Arne_v >> llsqrt() ville jeg sige er ret ineffektiv mod sqrt() og en conventering... Men det er godt at se lidt rigtig hacker ånd :-)

Du har ret med hensyn til effektiviteten.

Jeg troede faktisk at det ville være rimeligt hurtigt, fordi integer
kvadratrod principielt burde være hurtigere end floating point
kvadratrod (fordi der gennemsnitligt skal regnes på færre bit).

Men x86 har jo kvadratrods instruktion og macro kode kan nu engang
ikke konkurrere med micro kode.

Og jeg har heller ikke kunnet optimere det nok.

Test kode:

#include <stdio.h>
#include <math.h>
#include <time.h>

long long llsqrt1(long long v)
{
  long long tmp,res;
  res = v / 2;
  do
  {
      tmp = res;
      res = (v / tmp + tmp) / 2;
  } while((res > tmp + 1) || (tmp > res + 1));
  tmp = res;
  res = (v / tmp + tmp) / 2;
  return (res < tmp ) ? res : tmp;
}

long long llsqrt2(long long v)
{
  long long tmp,res;
  if(v < 10000)
  {
    if(v < 100)
    {
        if(v < 10)
        {
          res = 2;
        }
        else
        {
          res = 7;
        }
    }
    else
    {
        if(v < 1000)
        {
          res = 22;
        }
        else
        {
          res = 71;
        }
    }
  }
  else
  {
    if(v < 1000000)
    {
        if(v < 100000)
        {
          res = 316;
        }
        else
        {
          res = 707;
        }
    }
    else
    {
        if(v < 10000000)
        {
          res = 2236;
        }
        else
        {
          res = 7071;
        }
    }
  }
  do
  {
      tmp = res;
      res = (v / tmp + tmp) >> 1;
  } while((res > tmp + 1) || (tmp > res + 1));
  tmp = res;
  res = (v / tmp + tmp) >> 1;
  return (res < tmp ) ? res : tmp;
}

int main()
{
  int i;
  long long v,res;
  time_t t1,t2;
  t1 = time(NULL);
  for(i = 2; i <= 100000000; i++) {
      v = i;
      res = llsqrt1(v);
  }
  t2 = time(NULL);
  printf("homemade base %ld %d sec\n",(long)res,t2-t1);
  t1 = time(NULL);
  for(i = 2; i <= 100000000; i++) {
      v = i;
      res = llsqrt2(v);
  }
  t2 = time(NULL);
  printf("homemade optmized %ld %d sec\n",(long)res,t2-t1);
  t1 = time(NULL);
  for(i = 2; i <= 100000000; i++) {
      v = i;
      res = (long long)sqrt((double)v);
  }
  t2 = time(NULL);
  printf("builtin %ld %d sec\n",(long)res,t2-t1);
  return 0;
}

llsqrt kan dog stadigvæk regne op til 63 bits integers, mens konvertering
til double og sqrt vil give fejl på mere end 53 bit.

Man kan måske lave en mere optimal løsning ved at bruge følgende formel til at udregne sqrt: pow2(log2(n)/2), det er specielt anvendeligt hvis man ikke behøver en helt præsis løsning.
Man kan finde log2 ved at finde MSB (mange processorer kan gøre det i en instruktion) og bruge et tabel opslag til de resternde bits, størrelsen af tabellen afhænger af hvilken præsision der ønskes.
pow2 er også ret let at lave.
Jeg roder lidt med en løsning...

Er der ikke for stor usikkerhed på log2(n)==MSB ?

int(log2(n))===MSB, men så gælder formlen ikke og ellers er det jo
kun MSB <= log2(n) < MSB+1.

Man kunne selvfølgelig bruge det som et udgangspunkt for Newton.

int(log2(n))===MSB-1
MSB-1 <= log2(n) < MSB

Man kan ikke nøjes med at finde MSB, man skal også bruge de efterfølgende bits.
Princippet er at man skifter værdien op intil MSB er sat, derefter tager man de efterfølgende bits i værdien og bruger disse i et tabel opslag:

int64 log2(int64 val)
{
  int64 n, m;

  for(n = 63; !(*(((unsigned char *)&val) + 7)); n -= 8, val <<= 8 );
  for(; !(*(((unsigned char *)&val) + 7) & 0x80) && n; n--, val <<= 1);
  /* Nu har man skiftet val op så MSB er sat, n er 63 - antal pladser der er skiftet */

  m = (unsigned int )((val >> 55) & (256 - 1));
  /* m er nu de 8 bit efter MSB */

  return (n << 8) + nn_256[m]; /* nn_256 er en tabel med 256 værdier */
}

Mere info følger.

Her er et program der udregner kvadrat rødder. Test programmet udregner roden af 1*1, 2*2 ... 1023*1023 og beregner hvor mange der fejl der er.
Der er fejl på 598 af disse, ønskes der større presision skal man lave tabellerne nn_256 og pow_tab større. Der kan også opnåes lidt færre fejl ved at håndtrimme tabellerne.
Det er mere end 5 år siden jeg arbejdede med det sidst og jeg må indrømme at jeg har glemt hvor tallet 850 kommer fra og hvordan jeg beregnede værdierne i pow_tab.

#include <windows.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define int64 __int64
#define TABLEN 256

#ifdef __GNUC__
#define LLARG "%lld"
#else
#define LLARG "%I64d"
#endif

int64 nn_256[TABLEN];
static int64 pow_tab[256] =
{/*          0    1    2    3    4    5    6    7    8    9 */
/* 000 */  512,  513,  515,  516,  518,  519,  520,  522,  523,  525,
/* 010 */  526,  527,  529,  530,  532,  533,  535,  536,  538,  539,
/* 020 */  540,  540,  543,  545,  546,  548,  549,  551,  552,  554,
/* 030 */  555,  557,  558,  560,  561,  563,  564,  566,  567,  569,
/* 040 */  571,  572,  573,  575,  577,  578,  580,  581,  583,  585,
/* 050 */  586,  588,  589,  591,  593,  595,  596,  597,  599,  601,
/* 060 */  602,  604,  605,  607,  609,  611,  612,  614,  616,  617,
/* 070 */  619,  621,  622,  624,  626,  627,  629,  631,  632,  634,
/* 080 */  636,  638,  639,  641,  643,  644,  646,  648,  650,  652,
/* 090 */  653,  655,  657,  659,  660,  663,  664,  666,  668,  669,
/* 100 */  672,  673,  675,  677,  679,  680,  682,  684,  686,  688,
/* 110 */  690,  692,  693,  695,  697,  699,  700,  703,  705,  707,
/* 120 */  709,  711,  712,  714,  716,  719,  720,  722,  724,  727,
/* 130 */  728,  730,  732,  734,  736,  738,  740,  742,  744,  746,
/* 140 */  748,  750,  752,  754,  756,  758,  760,  762,  764,  766,
/* 150 */  769,  771,  773,  775,  777,  779,  781,  783,  785,  787,
/* 160 */  790,  792,  794,  796,  799,  800,  803,  805,  807,  809,
/* 170 */  811,  813,  816,  818,  820,  823,  825,  827,  828,  831,
/* 180 */  834,  836,  839,  840,  843,  845,  847,  849,  852,  855,
/* 190 */  856,  859,  861,  863,  866,  868,  875,  873,  875,  878,
/* 200 */  880,  882,  885,  887,  890,  892,  895,  897,  899,  902,
/* 210 */  904,  907,  909,  911,  914,  916,  919,  921,  924,  927,
/* 220 */  929,  931,  935,  936,  939,  942,  944,  947,  949,  952,
/* 230 */  954,  956,  960,  962,  965,  967,  970,  973,  975,  978,
/* 240 */  981,  983,  986,  989,  991,  994,  997,  999,  1002, 1005,
/* 250 */  1007, 1010, 1013, 1016, 1018, 1020
};

int64 pow2(int64 val)
{
  int64 a;
  a = val >> 8;

  if(a >= 10)
    return (pow_tab[val & 0x00FF]) << (a - 9);

  return ((int64 )pow_tab[val & 0x00FF] + (9 - a)/2) >> (9 - a);
}

int64 log2(int64 val)
{
  int64 n, m;

  for(n = 63; !(*(((unsigned char *)&val) + 7)); n -= 8, val <<= 8 );
  for(; !(*(((unsigned char *)&val) + 7) & 0x80) && n; n--, val <<= 1);

  m = (unsigned int )((val >> 55) & (TABLEN - 1));

  return (n << 8) + nn_256[m];
}

int64 MySqrt(int64 val)
{
  return pow2(log2(val)/2);
}

int main(void)
{
  int64 n;
  int64 r1;
  int64 c = 0;
  int i;
  double d;

  for(i = 0, d = 1.001953125; i < 256; i++, d += 1.0/256.0)
    nn_256[i] = 850*log10(d);

  for(n = 1; n < 1024; n++)
  {
    r1 = MySqrt(n*n);

    if(r1 != n)
    {
      fprintf(stderr, "N: " LLARG  " " LLARG "\n", n, r1);
      c++;
    }
  }
  printf("Errors: " LLARG " of " LLARG "\nHit enter key to exit\n", c, n);
  getchar();
  return 0;
}