Как я могу улучшить этот метод квадратного корня?

Можете попробовать r = x >> 1;

вместо / 2 (также в другом месте вы используете 2). Это может дать вам небольшое преимущество. Я бы также переместил 100 в петлю. Вероятно, ничего, но мы говорим здесь о тиках.

просто проверяю сейчас.

РЕДАКТИРОВАТЬ: Исправлен> в >>, но он не работает для двойников, так что не обращайте внимания. встраивание 100 не дало мне увеличения скорости.

Ну, собственная функция Sqrt (), вероятно, не реализована в C #, скорее всего, она будет реализована на низкоуровневом языке и, безусловно, будет использовать более эффективный алгоритм. Так что пытаться соответствовать его скорости, вероятно, бесполезно.

Однако, что касается простой попытки оптимизировать вашу функцию для heckuvit, страница Википедии, на которую вы ссылаетесь, рекомендует "начальное предположение" равным 2 этажу (D / 2), где D представляет количество двоичных цифр в номере. Вы можете попробовать это, я не вижу ничего другого, что можно было бы значительно оптимизировать в вашем коде.

Поскольку вы сказали, что приведенный ниже код был недостаточно быстрым, попробуйте следующее:

    static double guess(double n)
    {
        return Math.Pow(10, Math.Log10(n) / 2);
    }

Он должен быть очень точным и, надеюсь, быстрым.

Вот код для описанной начальной оценки здесь . Вроде неплохо. Используйте этот код, а затем вы также должны выполнять итерацию, пока значения не сойдутся в пределах эпсилона разницы.

    public static double digits(double x)
    {
        double n = Math.Floor(x);
        double d;

        if (d >= 1.0)
        {
            for (d = 1; n >= 1.0; ++d)
            {
                n = n / 10;
            }
        }
        else
        {
            for (d = 1; n < 1.0; ++d)
            {
                n = n * 10;
            }
        }


        return d;
    }

    public static double guess(double x)
    {
        double output;
        double d = Program.digits(x);

        if (d % 2 == 0)
        {
            output = 6*Math.Pow(10, (d - 2) / 2);
        }
        else
        {
            output = 2*Math.Pow(10, (d - 1) / 2);
        }

        return output;
    }

Вместо того, чтобы прервать цикл и затем вернуть r, вы можете просто вернуть r. Не может обеспечить заметного увеличения производительности.

Замена деления на 2 битовым сдвигом вряд ли приведет к такой большой разнице; учитывая, что деление происходит на константу, я надеюсь, что компилятор достаточно умен, чтобы сделать это за вас, но вы также можете попробовать его, чтобы увидеть.

У вас гораздо больше шансов получить улучшение, выйдя из цикл раньше, поэтому либо сохраните новое r в переменной и сравните со старым r, либо сохраните x / r в переменной и сравните это с r перед выполнением сложения и деления.

Определение допуска и возврат раньше, когда последующие итерации попадают в пределы этого допуска.

With your method, each iteration doubles the number of correct bits.

Using a table to obtain the initial 4 bits (for example), you will have 8 bits after the 1st iteration, then 16 bits after the second, and all the bits you need after the fourth iteration (since a double stores 52+1 bits of mantissa).

For a table lookup, you can extract the mantissa in [0.5,1[ and exponent from the input (using a function like frexp), then normalize the mantissa in [64,256[ using multiplication by a suitable power of 2.

mantissa *= 2^K
exponent -= K

After this, your input number is still mantissa*2^exponent. K must be 7 or 8, to obtain an even exponent. You can obtain the initial value for the iterations from a table containing all the square roots of the integral part of mantissa. Perform 4 iterations to get the square root r of mantissa. The result is r*2^(exponent/2), constructed using a function like ldexp.

EDIT. I put some C++ code below to illustrate this. The OP's function sr1 with improved test takes 2.78s to compute 2^24 square roots; my function sr2 takes 1.42s, and the hardware sqrt takes 0.12s.

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

double sr1(double x)
{
  double last = 0;
  double r = x * 0.5;
  int maxIters = 100;
  for (int i = 0; i < maxIters; i++) {
    r = (r + x / r) / 2;
    if ( fabs(r - last) < 1.0e-10 )
      break;
    last = r;
  }
  return r;
}

double sr2(double x)
{
  // Square roots of values in 0..256 (rounded to nearest integer)
  static const int ROOTS256[] = {
    0,1,1,2,2,2,2,3,3,3,3,3,3,4,4,4,4,4,4,4,4,5,5,5,5,5,5,5,5,5,5,6,6,6,6,6,6,6,6,6,6,6,6,
    7,7,7,7,7,7,7,7,7,7,7,7,7,7,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,9,9,9,9,9,9,9,9,9,9,9,9,9,
    9,9,9,9,9,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,11,11,11,11,11,
    11,11,11,11,11,11,11,11,11,11,11,11,11,11,11,11,11,12,12,12,12,12,12,12,12,12,12,12,12,
    12,12,12,12,12,12,12,12,12,12,12,12,13,13,13,13,13,13,13,13,13,13,13,13,13,13,13,13,13,
    13,13,13,13,13,13,13,13,13,14,14,14,14,14,14,14,14,14,14,14,14,14,14,14,14,14,14,14,14,
    14,14,14,14,14,14,14,14,15,15,15,15,15,15,15,15,15,15,15,15,15,15,15,15,15,15,15,15,15,
    15,15,15,15,15,15,15,15,15,16,16,16,16,16,16,16,16,16,16,16,16,16,16,16,16 };

  // Normalize input
  int exponent;
  double mantissa = frexp(x,&exponent); // MANTISSA in [0.5,1[ unless X is 0
  if (mantissa == 0) return 0; // X is 0
  if (exponent & 1) { mantissa *= 128; exponent -= 7; } // odd exponent
  else { mantissa *= 256; exponent -= 8; } // even exponent
  // Here MANTISSA is in [64,256[

  // Initial value on 4 bits
  double root = ROOTS256[(int)floor(mantissa)];

  // Iterate
  for (int it=0;it<4;it++)
    {
      root = 0.5 * (root + mantissa / root);
    }

  // Restore exponent in result
  return ldexp(root,exponent>>1);
}

int main()
{
  // Used to generate the table
  // for (int i=0;i<=256;i++) printf(",%.0f",sqrt(i));

  double s = 0;
  int mx = 1<<24;
  // for (int i=0;i<mx;i++) s += sqrt(i); // 0.120s
  // for (int i=0;i<mx;i++) s += sr1(i);  // 2.780s
  for (int i=0;i<mx;i++) s += sr2(i);  // 1.420s
}

Здесь вы выполняете метод Ньютона для поиска корня . Так что вы можете просто использовать более эффективный алгоритм поиска корня. Вы можете начать поиск здесь .

Я тоже смотрел на это в учебных целях. Вас могут заинтересовать две модификации, которые я пробовал.

Первый заключался в использовании приближения ряда Тейлора первого порядка по x0:

    Func<double, double> fNewton = (b) =>
    {
        // Use first order taylor expansion for initial guess
        // http://www27.wolframalpha.com/input/?i=series+expansion+x^.5
        double x0 = 1 + (b - 1) / 2;
        double xn = x0;
        do
        {
            x0 = xn;
            xn = (x0 + b / x0) / 2;
        } while (Math.Abs(xn - x0) > Double.Epsilon);
        return xn;
    };

Второй - попробовать третий порядок (более дорогой), итерация

    Func<double, double> fNewtonThird = (b) =>
    {
        double x0 = b/2;
        double xn = x0;
        do
        {
            x0 = xn;
            xn = (x0*(x0*x0+3*b))/(3*x0*x0+b);
        } while (Math.Abs(xn - x0) > Double.Epsilon);
        return xn;
    };

Я создал вспомогательный метод для измерения времени функций

public static class Helper
{
    public static long Time(
        this Func<double, double> f,
        double testValue)
    {
        int imax = 120000;
        double avg = 0.0;
        Stopwatch st = new Stopwatch();
        for (int i = 0; i < imax; i++)
        {
            // note the timing is strictly on the function
            st.Start();
            var t = f(testValue);
            st.Stop();
            avg = (avg * i + t) / (i + 1);
        }
        Console.WriteLine("Average Val: {0}",avg);
        return st.ElapsedTicks/imax;
    }
}

оригинальный метод был быстрее, но опять же, может быть интересно :)

float InvSqrt (float x){
  float xhalf = 0.5f*x;
  int i = *(int*)&x;
  i = 0x5f3759df - (i>>1);
  x = *(float*)&i;
  x = x*(1.5f - xhalf*x*x);
  return x;}

Это мой любимый быстрый квадратный корень. На самом деле это обратный квадратный корень, но вы можете перевернуть его после, если хотите....I не может сказать, быстрее ли он, если вы хотите квадратный корень, а не обратный квадратный корень, но это чертовски круто, как раз то же самое.
http://www.beyond3d.com/content/articles/8/