Searching \ for '[PIC]: [MATH]: Power function, part 2' in subject line. ()
Help us get a faster server
FAQ page: www.piclist.com/techref/microchip/math/index.htm?key=math
Search entire site for: '[MATH]: Power function, part 2'.

Exact match. Not showing close matches.
'[PIC]: [MATH]: Power function, part 2'
2000\09\29@175856 by

Hi List / Nikolai,

Thanks again for the log2() function Nikolai.

I've implemented it in my program and it works perfectly now (took a
little time to get it work because of a really stupid mistake of me).
I have also implemented the multiplication by y (compared to the log2
function, this is a piece of cake of course), so the the
power-function is now finished till :

y * log2(x)

where the result is in an unsigned Q6.10 (implicit minus) notation.

Now I only have to take care of the exp() function. But the problem is
that I'm not a math-wizard like Nikolai. So I would really appreciate
it if you could help me out with that routine also. The output of the
exp() function has to be in the unsigned Q0.16 format again.

Again thanks a lot.

Daniel...

--

part 1 6460 bytes content-type:text/plain; charset=ISO-8859-1 (decoded quoted-printable)

Hi Daniel,

Well, let's calculate y=2^x! x is Q6.10 notation (always
negative with implicit sign). Output y has Q0.16 notation
(always less then 1). Therefore x should not be zero (2^0 =
1).

y is in the range from 0.9993 (x = -1/1024) to 5e-20 (x
= -64 + 1/1024). Okay....

y=2^x

y=2^-(32*x.15+16*x.14+...1*x.10 + 0.5*x.9+0.25*x.8+...)
\========================/   \===================/
x_int                         x_frac

x_int = 0:1:63
x_frac = 0:1/1024:(1-1/1024)

y=2^-(x_int+x_frac)

y=2^(-x_int) * 2^(-x_frac)

y= (2^(-x_frac)) >> x_int
=========================

So we first find 2^-x_frac using a table and then shift the
result right x_int times.

First we should convert the fraction part through a table to
a new fraction part. 2^(-0..-1) = 0.5..1

; y = 2^x2
;Input: x2hi:x2lo
;Output: yhi:ylo

;Form table index from 4 higher bits of x2 fraction part
rlf x2lo, f
rlf x2ho, w
movwf index
rlf x2lo, w
rlf index, f
;here middle 6 bits of x2lo contain multiplier for linear
;interpolation, and index has 4 lower bits of table pointer
rlf index, w    ;multiply index by 2 for 16 bit
andlw 0x1E      ;and clear lsb and higher bits
movwf index
call exp2_table

;clear result and subtract the first point from it
clrf yhi
clrf ylo
movf templo, w
subwf ylo, f
movf temphi, w
skpc
incfsz temphi, w
subwf yhi, f
incf index, f
incf index, f
call exp2_table
;find difference with previous point to temphi:templo
movf ylo, w
movf yhi, w
skpnc
incfsz yhi, w
;multiply the difference by the next 6 fraction bits in
;x2lo and divide by 64. Then subtract it from the current
;result in yhi:ylo
movlw 6
movwf index     ;use index as a loop counter
rrf x2lo, f     ;align the 6 bits to the lsb

clrf temp       ;we will use another temp
;for mul purposes

rrf yhi, f      ;rotate y right 2 times
rrf ylo, f      ;to make easy multiplication.
rrf temp, f     ;lower bits will go to temp
rrf yhi, f
rrf ylo, f
rrf temp, f
exp2loop
rrf x2lo, f      ;get next multiplier bit
bnc exp2loopnext ;skip subtraction if carry is clear

movf templo, w
subwf temp, f
movf temphi, w
skpc
incfsz temphi, w
subwf ylo, f
skpc
decf yhi, f

exp2loopnext

rrf yhi, f
rrf ylo, f
rrf temp, f
decfsz index, f
goto exp2loop
;shift the result 8 bits left
movf ylo, w
movwf yhi
movf temp, w
movwf ylo

;Now we will shift y right the number of times
;equal to the integer part of x2
clrc            ;align integer bits
rrf x2hi, f     ;to lsb
clrc
rrf x2hi, f

test x2hi
skpnz
return         ;done
exp2loop2
clrc
rrf y2hi, f
rrf y2lo, f
decfsz x2hi, f
goto exp2loop2
return          ;done

Now we need to calculate the table...
1-2.^(0:-1/16:-1)
ans =
Columns 1 through 4                   0   0.04239671930143   0.08299595679533   0.12187391981335
Columns 5 through 8    0.15910358474629   0.19475483402537   0.22889458729603   0.26158692703025
Columns 9 through 12    0.29289321881345   0.32287222653155   0.35158022267450   0.37907109396326
Columns 13 through 16    0.40539644249864   0.43060568262165   0.45474613366737   0.47786310878629
Column 17    0.50000000000000

Multiplied by 65536 and rounded
round(65536*(1-2.^(0:-1/16:-1)))
ans =
Columns 1 through 6            0        2779        5439        7987       10427       12763
Columns 7 through 12        15001       17143       19195       21160       23041       24843
Columns 13 through 17        26568       28220       29802       31317       32768
»
So, the table would be:

exp2_table
movf index, w
call exp2tableStart
movwf templo
incf index, w
call exp2tableStart
movwf temphi
return

exp2tableStart
DT     0 & 0xFF,     0 >> 8,  2779 & 0xFF,  2779 >> 8
DT  5439 & 0xFF,  5439 >> 8,  7987 & 0xFF,  7987 >> 8
DT 10427 & 0xFF, 10427 >> 8, 12763 & 0xFF, 12763 >> 8
DT 15001 & 0xFF, 15001 >> 8, 17143 & 0xFF, 17143 >> 8
DT 19195 & 0xFF, 19195 >> 8, 21160 & 0xFF, 21160 >> 8
DT 23041 & 0xFF, 23041 >> 8, 24843 & 0xFF, 24843 >> 8
DT 26568 & 0xFF, 26568 >> 8, 28220 & 0xFF, 28220 >> 8
DT 29802 & 0xFF, 29802 >> 8, 31317 & 0xFF, 31317 >> 8
DT 32768 & 0xFF, 32768 >> 8
IF ((\$ - 1) >> 8) - ((exp2tableStart + 1) >> 8) != 0
error 'exp2 table crossed 8-bit boundary'
ENDIF

And the cleaned up version is attached... Yeah, this one was
harder and not really the reverse of log!

Don't forget to check x and y for zero in x^y.

Hope it works for you :)
Nikolai

---- Original Message ----
From: D. Schouten <danielsxs4all.nl>
Sent: Saturday, September 30, 2000 0:56:03
To: Nikolai Golovchenko
Subj: [PIC]: [MATH]: Power function, part 2

{Quote hidden}

part 2 8615 bytes content-type:application/octet-stream; name="log2exp2.asm" (decode)

part 3 144 bytes
--
http://www.piclist.com hint: The list server can filter out subtopics
(like ads or off topics) for you. See http://www.piclist.com/#topics

'[PIC]: [MATH]: Power function, part 2'
2000\10\03@131101 by
Hi List / Nikolai,

Thanks for the Exp2 routine Nikolai!

I've also implemented this last piece of the power function. It works
perfectly now but it took a little time to debug. I found two little
bugs in the Exp2 routine in the log2exp2.asm attachment:

1.) the code :

tstf x2hi
skpnz
return         ;done

has to be placed UNDER :

clrc            ;align integer bits
rrf x2hi, f     ;to lsb
clrc
rrf x2hi, f

instead of above. This is done right in the explanation posted on the
list but wrong in the attached log2exp2.asm.

2.) in the piece of code :

movlw 17     ;limit shifts number to 16 maximum
subwf x2hi, w
movlw 16
skpc
movwf x2hi

the instruction :

skpc

has to be changed in :

skpnc

to limit shifts number to 16 like Nikolai says in the comment.

Again, thank you very much Nikolai. The routine is not only very
accurate, it's also efficient. The whole power function takes less
then 1000 cycles to calculate and uses less then 300 words of
codespace. WOW!

Again thanks!

Bye,

Daniel...

--
http://www.piclist.com hint: To leave the PICList
piclist-unsubscribe-requestmitvma.mit.edu

More... (looser matching)
- Last day of these posts
- In 2000 , 2001 only
- Today
- New search...