Subject: Re: 16 bit fixed point SQRT for PIC16
Date: Tue, 19 Feb 2002 20:12:56 +0000 (UTC)
From: "Alexandr A. Redchuck"
Organization: ReAl at home
Newsgroups: fido7.ru.embedded19-Feb-02 18:05 Dmitry Orlov wrote to All:
DO> А не завалялось ли у кого асмовской реализации (желательно или вообще
DO> без
DO> использования стека или с минимальным) subj? Я пока что использую
DO> где-то откопанный такой код:
[^Y]
DO> for(count = 0; count < 8; count++) {
DO> src &= (mask >> count*2);
DO> tmp_src <<= 2;
DO> tmp_src |= src >> (7 - count)*2;
DO> tmp_dst = dst;
DO> dst <<= 1;
DO> if((tmp_src != 0 )||(dst!=0)) {
DO> tmp_dst <<= 2;
DO> tmp_dst += 1;
DO> if(tmp_dst <= tmp_src) {
DO> dst += 1;
DO> tmp_src -= tmp_dst;
DO> }
DO> }
DO> }
DO> tmp_src <<= 2;
DO> tmp_dst = dst;
DO> tmp_dst <<= 2;
DO> tmp_dst += 1;
DO> if (tmp_dst <= tmp_src) dst += 1;
DO> В принципе он меня устраивает, но может быть есть что-то существенно
DO> более эффективное по скорости, или размеру?
Думаю, что и по скорости, и по размеру. Здесь совсем не оптимизированы
сдвиги, разве можно делать столько многобайтных сдвигов на переменную
величину. Мог бы хоть их почистить...
В ru.algorithms когда-то этот вопрос поднимался, Винокуров Андрей приводил
объяснение "школьного" метода "цифра за цифрой" ("корень в столбик") для
десятичной системы и его асм386 реализацию в двоичном виде.
Ниже выдержка из моего письма туда же с несколько другим методом (регистр
последовательных приближений).
В асмовой реализации для 386 он дал сравнимое быстродействие.
Тут я приведу только С-шный кусок письма.
From: "Alexandr A. Redchuck"
Newsgroups: fido7.ru.algorithms
Subject: [NEWS] isqrt
Date: 12 Apr 1998 16:51:28 +0400> РПП:
> поднимаем в регистре старший бит. Если функция содержимого регистра
> меньше цели, оставляем этот бит, если больше - сбрасываем.
> Повторяем вправо со всеми битами до младшего.
>
> В тупом варианте для корня:
>
> unsigned short sar_sqrt(unsigned long from) {
> unsigned short sqr,mask=0x8000;
> do {
> sqr |= mask;
> if( sqr*sqr > from ) sqr &= ~mask;
> } while( mask >>= 1 );
> return sqr;
> }
>
> Дальше для ухода от возведения в квадрат переходим к постепенному
> приближению. Используем
> (sqr+mask)*(sqr+mask) - sqr*sqr = 2*sqr*mask + mask*mask
> и учитываем, что в mask поднята только одна единичка, т.е. умножение
> на нее это сдвиг (в десятичном варианте таки надо множить -- стоИт
> (10*2*R+t)*t и t приходится искать "на глаз" среди 10 весов разряда).
> Можно гнать цикл по j от 15 до 0 и превратить
> 2*sqr*mask + mask*mask
> в
> (sqr << (j+1)) + ( 1 << (j+j) )
> но при этом много сдвигов на переменную величину, которые не везде
> легко делаются.
>
> После оптимизации сдвигов
>
> unsigned short isqrt( unsigned long from) {
> unsigned long mask = 0x40000000, sqr = 0, temp;
> do {
> temp = sqr | mask;
> sqr >>= 1;
> if( temp < = from ) {
> sqr |= mask;
> from -= temp;
> }
> } while( mask >>= 2 );
> // можно дать еще округление if( sqr < from ) ++sqr;
> return (unsigned short)sqr;
> }
Тут для корня из длинного, для корня из 16-битного надо соответственно
уменьшить разрядность переменных и начать с маски 0x4000.
wbr,
--
/* Alexandr Redchuck, Kyiv, Ukraine */
/* real на real тчк kiev тчк ua */
