Перейти к содержанию
    

Очень нужна помощь в оптимизации алгоритма извлечения квадратного корня

Гость TSerg

Топик-стартер слушает только себя в попытке изобрести велосипед.:)

 

Вам уже приводились идеи и работающие варианты.

Науке известны разные способы ивлечения целочисленного корня на разные вкусы.

 

Еще раз приведу (для 16р сетки), число иттераций всегда 9 и никаких умножений.

 

function QuickSqrt_(Val : word) : word;

var bitSqr : word;

begin

bitSqr := $10000;

Result := 0;

 

While bitSqr <> 0 do begin

if (Val >= bitSqr + Result) then begin

Dec(Val, bitSqr + Result);

Result := (Result shr 1) or bitSqr;

end

else Result := Result shr 1;

bitSqr := bitSqr shr 2;

end;

end;

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

Топик-стартер слушает только себя в попытке изобрести велосипед.:)

 

Вам уже приводились идеи и работающие варианты.

Науке известны разные способы извлечения целочисленного корня на разные вкусы.

 

Уважаемый TSerg, ну тё Вы ругаетесь? Себя я, конечно же, слушаю))) Велосипед? Да, изобретаю! Не мешайте мне побыть гениальным:biggrin: :biggrin: :biggrin: :biggrin: :biggrin:

 

Ваш вариант действительно очень эффективен. Да, выполняется за всего за 9 циклов! Я бы добавил его в начало темы, но не могу. Первое сообщение не доступно уже для изменения.

 

Более того, в последнем Вашем коде есть ошибка: bitSqr : word. Это говорит о том, что Вы лично его не проверяли в деле.

 

Судя по тому, что маске присваивается значение 0x010000, маска bitSqr должна быть 32-х разрядной (как и было в первом Вашем варианте). Далее в коде все операции с ней происходят как с 32-х разрядным числом (компилятор не будет динамически менять разрядность), это требует дополнительной памяти данных и инструкций. Это первое.

 

Второе. Объем кода в Вашем варианте мне показался больше, чем в моем (к тому же, учитывая операции с 32-х разрядным числом в микроконтроллере с 8-ми разрядной архитектурой (даже с аппаратным умножителем 8х8), например). Но это на вскидку. Если все совсем делать серьезно, то, по идее, оба алгоритма нужно перевести на ассемблер микроконтроллера PIC16, например, и сравнить.

 

Но, Ваш код прост. Мой имеет один (?) недостаток - операция умножения. В микроконтроллере с аппаратным умножителем это не будет проблемой, в остальных да.

 

Вариантов, действительно, было множество. Но ищется оптимальный (имеющий малый объем кода, простоту вычислений, быстрый). Не все они отвечают этим требованиям. И вкус, часто определяется архитектурой вычислительного ядра. Ищу золотую середину)))

 

ЗЫ. Кстати, очень интересно какие все же идеи есть у beaRTS)))

Изменено пользователем Rostislav

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

TSerg, Ваш вариант описывается тут: http://en.wikipedia.org/wiki/Methods_of_co...ng_square_roots

 

Ребят, кто-нибудь понимает предложенные там еще два метода:

 

          1  2. 3  4 
        ___________
       /
     \/  01 52.27 56

         01                   1*1 <= 1 < 2*2                 x = 1
         01                     y = x*x = 1*1 = 1
         00 52                22*2 <= 52 < 23*3              x = 2
         00 44                  y = (20+x)*x = 22*2 = 44
            08 27             243*3 <= 827 < 244*4           x = 3
            07 29               y = (240+x)*x = 243*3 = 729
               98 56          2464*4 <= 9856 < 2465*5        x = 4
               98 56            y = (2460+x)*x = 2464*4 = 9856
               00 00          Algorithm terminates: Answer is 12.34

          1. 4  1  4  2
          _____________
       /
     \/  02.00 00 00 00

         02                  1*1 <= 2 < 2*2                 x = 1
         01                    y = x*x = 1*1 = 1
         01 00               24*4 <= 100 < 25*5             x = 4
         00 96                 y = (20+x)*x = 24*4 = 96
            04 00            281*1 <= 400 < 282*2           x = 1
            02 81              y = (280+x)*x = 281*1 = 281
            01 19 00         2824*4 <= 11900 < 2825*5       x = 4
            01 12 96           y = (2820+x)*x = 2824*4 = 11296
               06 04 00      28282*2 <= 60400 < 28283*3     x = 2
                             The desired precision is achieved:
                             The square root of 2 is about 1.4142

?

Изменено пользователем Rostislav

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

ЗЫ. Кстати, очень интересно какие все же идеи есть у beaRTS)))

да какие там идеи... =) тож начитался много чего (по диагонали)

мне корень нужен для нахождения модуля вектора, т.е.: имею синфазную и квадратурную составляющую сигнала (I и Q), т.е. как бы координаты вектора, ну и надо модуль найти, чтобы два канала I и Q соединить в один. Так вот. У меня на проце TMS320F28335 есть FPU (в общем, floating-point понимает и перемножитель даже есть аппаратный), работаю с числами float (вещественными одинарной точности).

 

Ну и показался метод Ньютона (он же вроде вырождается в метод Вавилонский) более подходящий из всех, т.к. ему пофигу что считать: либо целые либо вещественны числа. Стал думать как делать для него первое приближение g0 (в книге Уоррен Генри. "Алгоритмические трюки для программистов такие обозначения" ). Ну а тут прям в масть ознакомился с неким алгоритмом "максимум альфа плюс минимум бета", который линейно аппроксимирует выражение sqrt( I^2 + Q^2) с приемлемой точностью и быстротой (точность зависит от коэффициентов альфа и бета). ну и решил первое приближение g0 делать при помощи него.

 

итого, программка работает. только мне еще вкурить надо как должным образом проверить ее производительность во времени (что-то типа профилирования), сколько итераций расчета по методу ньютона производит с таким начальным/первым приближением. и сравнить по времени выполнения со стандартными функциями С++. Ну и величину ошибки тоже вычислить.

 

Что успел бегло подметить: что количество итераций для разных чисел - разное, что многие значения высчитываются точно в сравнении со стандартными функциями, а некоторые с некоторой ошибкой...

 

В общем статистику еще не проводил и не оформлял. но стопудово сделаю какую-нить статейку.

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

ЗЫ. Кстати, очень интересно какие все же идеи есть у beaRTS)))

да какие там идеи... =) тож начитался много чего (по диагонали)

мне корень нужен для нахождения модуля вектора, т.е.: имею синфазную и квадратурную составляющую сигнала (I и Q), т.е. как бы координаты вектора, ну и надо модуль найти, чтобы два канала I и Q соединить в один. Так вот. У меня на проце TMS320F28335 есть FPU (в общем, floating-point понимает и перемножитель даже есть аппаратный), работаю с числами float (вещественными одинарной точности).

 

Ну и показался метод Ньютона (он же вроде вырождается в метод Вавилонский) более подходящий из всех, т.к. ему пофигу что считать: либо целые либо вещественны числа. Стал думать как делать для него первое приближение g0 (в книге Уоррен Генри. "Алгоритмические трюки для программистов такие обозначения" ). Ну а тут прям в масть ознакомился с неким алгоритмом "максимум альфа плюс минимум бета", который линейно аппроксимирует выражение sqrt( I^2 + Q^2) с приемлемой точностью и быстротой (точность зависит от коэффициентов альфа и бета). ну и решил первое приближение g0 делать при помощи него.

 

итого, программка работает. только мне еще вкурить надо как должным образом проверить ее производительность во времени (что-то типа профилирования), сколько итераций расчета по методу ньютона производит с таким начальным/первым приближением. и сравнить по времени выполнения со стандартными функциями С++. Ну и величину ошибки тоже вычислить.

 

Что успел бегло подметить: что количество итераций для разных чисел - разное, что многие значения высчитываются точно в сравнении со стандартными функциями, а некоторые с некоторой ошибкой...

 

В общем статистику еще не проводил и не оформлял. но стопудово сделаю какую-нить статейку.

 

ну вот сырые данные в картинке во вложении.

Под итерациями я здесь подразумеваю сколько раз алгоритму требуется пройтись по циклу while() :

    while( firstEstimation - nextEstimation != 0 ) {    
        firstEstimation = nextEstimation;

        nextEstimation = ( firstEstimation + quotient )/2;

        //for next iteration check in while()
        quotient = radicand/firstEstimation;

        //for debug
        ++iterationCnt;
    } //while

 

 

а в main() прогоняю такие числа:

    for( int k = 0, g = 50; k < 50; ++k, --g) {
        temp = k*k + g*g;
        i = ( float )k;
        q = ( float )g;
        std::cout << "<cmath> sqrt(): " << sqrt(temp) << "     my sqrt(): " << sqrtIQ(i, q) << "   perform in " << iterationCnt << " iterations"<< std::endl;
    } // for

post-69111-1356166744_thumb.png

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

по поводу своей идеи начал новую тему. вот здесь кому интересно - зайдите , посоветуйте

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

Гость TSerg
Ваш вариант действительно очень эффективен. Да, выполняется за всего за 9 циклов! Я бы добавил его в начало темы, но не могу. Первое сообщение не доступно уже для изменения.

 

Более того, в последнем Вашем коде есть ошибка: bitSqr : word. Это говорит о том, что Вы лично его не проверяли в деле.

 

Это не мой вариант, но давно известный.

 

bitSqr := $FFFF;

 

Да, чисто механически оттранслировал на 16р:)

Delphi исправляет мою погрешность, но в микроконтроллере надо быть внимательнее, конечно же.

 

 

Поделиться сообщением


Ссылка на сообщение
Поделиться на другие сайты

Присоединяйтесь к обсуждению

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

Гость
Ответить в этой теме...

×   Вставлено с форматированием.   Вставить как обычный текст

  Разрешено использовать не более 75 эмодзи.

×   Ваша ссылка была автоматически встроена.   Отображать как обычную ссылку

×   Ваш предыдущий контент был восстановлен.   Очистить редактор

×   Вы не можете вставлять изображения напрямую. Загружайте или вставляйте изображения по ссылке.

×
×
  • Создать...