Tracy,
n refers to the number of bits in the fraction computed. You posted the routine. f refers to the result of that routine. when n = 11, f must be shifted 5 bits left, not 1 bit. The case with n = 16 doesn't require any shifting of f.
Fred,
As soon as I clean out some other stuff and spiff it up a bit I'll post it. Maybe later today.
Hi Tom, I thought so, but that is my confusion. Maybe it has to do with the way that value is used in a followup subroutine that expects it left justified in 16 bits, and that will become clear when you post the cleaned up code. The idea of making N=11 is to create the 11-bit address for table lookup directly, and then shift one more left +$D000 to turn it into a word address. If N=16, then there are 5 extra lsbs for interpolation.
Tracy,
That seems to be correct. You can use 16 bits without interpolating if you want as the last 5 bits aren't used when accessing the table. It will cost you five loops through your routine though.
It'll be a few more hours before I can post the code.
Well, here it is. I would like to hear from the experts in assembly coding about this object. I started one month ago on assembly programming so I positively know it can be done better.
Some warnings: The denominator and numerator are both unsigned integers. The numerator must be larger than the denominator. My application is such that the values I plug in are near ideal so no error or range checking is done. The ratio of numerator over denominator must be less than 32 or overflow will occur (2^32 is 32 bits). For low ratios the number of usable bits may be too low to be useful, ie 2^(1/2) would yield 1 for an answer. Not particularly informative. My range of use is approximately 2^18 to 2^24, yielding at a minimum 6 decimal digits to a max of 8 decimal digits. Accuracy in this range seems to be on the order of 15 bits, quite good actually. Also better than I thought could be achieved. Many thanks to all who jumped in with advice.
One more thing. I copied and pasted code from my application so I haven't tested the standalone version. I would of course like to hear of any errors.
Tom, If you won't mind I should take your code as an example to comment a little bit about "style" (But I also can do it in a PM..) On the whole it looks quite excellent, but there are always quibbles...
But there is one remark to logic, however:
Your code is a "dataflow server", computing the result of two cells into a third contninously.
However the client does not know to what arguments the result belongs!
When he reads the result to fast (but maybe this can be excluded?), the result can belong to
- the previous set of arguments
- the first new and the second old argument
depending on the position of the loop the server COG is in when you change the arguments
This will stabilize after a few microseconds of course
In physical dataflow applications this not a problem, as the arguments change slowly (continously with limited growth) only.
In "non-analog" applications one generally uses a command flag, saying " Hello! New set of arguments!", and the server clears this flag after completion. When the results come from limited set of values, one can also use one of the impossible results (-1, 0, MAXP, MINP..) as command flag, which then is reset implicitely when storing the result...
Thanks for the reply. In my application the denominator is a constant so only the numerator changes avoiding the problem you stated. I can see where in a real time system it could be a problem and would have to be addressed.
It looks fine to me and the structure is easy to follow. There are probably little tweaks that could be done to speed it up a little bit if that is an issue. Nice job on the interpolation, too!
It occurs to me that it would be possible to improve the accuracy of the interpolation by a couple of bits by using quadratic instead of linear interpolation. The linear interpolation on the b^x curve (concave upward) always overestimates the true value by a small amount that is maximum near the center point of the interpolation.
The speed of the asm routine is more than fast enough, though spin would be too slow.
OK, how do you do quadratic interpolation? I'll be glad to add in in. Any further accuracy improvements would be icing on the cake.
BTW, the posted code is not in the form I use it. My working code was modified so that others could easily test it. I must be lucky because the range of the numerator just happens to fall in a very nice spot. Like I said earlier, my demoninator is a constant so I don't have to wait on it.
Would attaching a hardware floating point math co-processor to the propeller solve this problem? Parallax sells one called the uM-FPU V3 for $20 USD, see Parallax part-number 604-00050. There's also a V2 part for $15 USD, see part-number 604-00030; I don't know what the difference is between V2 and V3.
I think these co-processors on offer from Parallax are from MicroMega Corp. See their site at: www.micromegacorp.com
MicroMega now seems to offer a V3.1 part that does not appear on the Parallax site; at least not as of today.
Attention! This is not a Co-Processor in the strict sense, but a cleverly programmed PIC (I think); the performance is somewhat below one COG (i.e. A Prop running FLOAT32 in 4 COGs will outperform this Chip considerably!!
However there are some goodies:
- 2 (medium speed!) 12 bit AD-Converters
- Very extensive set of mathematical routines, including matrix multiplication and a FFT (!), of course, as Badgers put it in a different context: "... memory permitting it..."
Nice toy, but an additional Prop gives you 5x the power for half the price.
However a very god companion for the BasicStamp, or when you just need the algorihms and not much speed
Carefully read the good manual before you buy (I didn't )
Thanks but no thanks. I am quite happy with the performance of my routine. However the uM-FPU looks like an interesting chip if you need its capabilities.
DeSilva,
Thanks for the comments. "Hau weg die Scheisser!"
The attached graph shows the error curve for linear interpolation, The difference between values interpolated between steps of x=1/(2^11) and the actual value of b^x in between those steps. The possible quadratic fit for each segment is evident. Like the second term of the taylor series at each point.
This is b^x with x near 0.5. The peak error at x=1 is twice as much as it is at x=0. I don't know if it is really worthwhile to work it out in order to apply a second order correction, if what you have now is adequate. The additional code might be more than would be required for direct computation with the shift/add methods.
Tracy,
My accuracy is more than adequate but I might try it for fun. I'm away from home right now so it will have to wait. If I come up with any code I'll post it in this thread.
Franz, DeSilva,
It's a quote from "Werner". My son taught it to me.
Fred Hawkins said...
Alas, new laptop is even lacks MSWorks. I was hoping for a tale of the best baked n fresh download.
Fred,
I've used google docs in the past, you can get to it here: docs.google.com, it's quite nice for making graphs, although not as advanced as Excel. Another office suite you could use is OpenOffice, it's a free, Open Source office "clone". It works quite nicely, I am sure the Excel-like app (calc) does most of everything you need.
Gr,
Mightor
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
| To know recursion, you must first know recursion.
In looking at it more, I realized that the interpolation itself was giving errors of only 0.00000002, which is around 24 bits of precision. A much larger error then is due to quantization. That's good news, because the result can be improved by finer linear interpolation. Instead of 2^5=32 steps, go in whole hog with 2^13=8192 steps. That amounts to 11 bits from the HUB table + 13 bits of linear interpolation, equals 24 total. Only at around 24 bits would a quadratic correction to the linear result make any difference.
Of course, if there are no more than 16 bits of precision going in, it does not make any sense to expect more than 16 out. But with, say, 24 significant bits in a floating point mantissa, it could make a big improvement.
Seeing as how I get the remainder through 16 bit division, it wouldn't make sense to go any further. Actually I'm quite happy with the results I'm getting thanks to you and the others.
Yes, you're right. My excel "what-if" was assuming that the interpolation points are perfect (within the accuracy of excel), but the points in the HUB table are only accurate to 16 bits, so that will in fact limit the resulting interpolation. I want to thank you too, for bringing up the topic for discussion.
One comment on your Excel graphs. It sure looks like the absolute value of a sine wave. Is there any way to mathematically relate the two? Even if it isn't a sine wave it looks like a good approximation to one.
It is a good parabola. Here is a graph of one cycle (red) and the least squares fit to a second order polynomial. R^2 =1.
I don't think there is any direct connection with sine, but who knows? For example, there is a deep connection in terms of periodicity, something that people don't usually think of when looking at the exponential and log funtions, which seem on the surface to be monotonic over their whole range. When calculating the trig functions using CORDIC, the reference table is the arctangent at binary arguments, but for CORDIC exponentials, the reference table is the hyperbolic arctangent.
The graph is from Excel, with the interpolation points calculated to the accuracy that Excel allows (15 significant digits, about 50 binary bits, I think). Wiith the HUB table, the interpolation points are accurate to 16 bits, so there is corresponding quantization error in the end points of the interpolation. So what I said in an earlier post about interpolating to 8192 steps and using a quadratic interpolation is BS. 16 bits with 32 linear interpolation points is about as good as it gets.
Comments
n refers to the number of bits in the fraction computed. You posted the routine. f refers to the result of that routine. when n = 11, f must be shifted 5 bits left, not 1 bit. The case with n = 16 doesn't require any shifting of f.
Fred,
As soon as I clean out some other stuff and spiff it up a bit I'll post it. Maybe later today.
Tom
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
That seems to be correct. You can use 16 bits without interpolating if you want as the last 5 bits aren't used when accessing the table. It will cost you five loops through your routine though.
It'll be a few more hours before I can post the code.
Tom
Some warnings: The denominator and numerator are both unsigned integers. The numerator must be larger than the denominator. My application is such that the values I plug in are near ideal so no error or range checking is done. The ratio of numerator over denominator must be less than 32 or overflow will occur (2^32 is 32 bits). For low ratios the number of usable bits may be too low to be useful, ie 2^(1/2) would yield 1 for an answer. Not particularly informative. My range of use is approximately 2^18 to 2^24, yielding at a minimum 6 decimal digits to a max of 8 decimal digits. Accuracy in this range seems to be on the order of 15 bits, quite good actually. Also better than I thought could be achieved. Many thanks to all who jumped in with advice.
One more thing. I copied and pasted code from my application so I haven't tested the standalone version. I would of course like to hear of any errors.
Feedback wanted!
Tom
Post Edited (TomS) : 8/26/2007 10:33:50 PM GMT
Tom, If you won't mind I should take your code as an example to comment a little bit about "style" (But I also can do it in a PM..) On the whole it looks quite excellent, but there are always quibbles...
But there is one remark to logic, however:
Your code is a "dataflow server", computing the result of two cells into a third contninously.
However the client does not know to what arguments the result belongs!
When he reads the result to fast (but maybe this can be excluded?), the result can belong to
- the previous set of arguments
- the first new and the second old argument
depending on the position of the loop the server COG is in when you change the arguments
This will stabilize after a few microseconds of course
In physical dataflow applications this not a problem, as the arguments change slowly (continously with limited growth) only.
In "non-analog" applications one generally uses a command flag, saying " Hello! New set of arguments!", and the server clears this flag after completion. When the results come from limited set of values, one can also use one of the impossible results (-1, 0, MAXP, MINP..) as command flag, which then is reset implicitely when storing the result...
Post Edited (deSilva) : 8/27/2007 5:56:30 AM GMT
Thanks for the reply. In my application the denominator is a constant so only the numerator changes avoiding the problem you stated. I can see where in a real time system it could be a problem and would have to be addressed.
Tom
It looks fine to me and the structure is easy to follow. There are probably little tweaks that could be done to speed it up a little bit if that is an issue. Nice job on the interpolation, too!
It occurs to me that it would be possible to improve the accuracy of the interpolation by a couple of bits by using quadratic instead of linear interpolation. The linear interpolation on the b^x curve (concave upward) always overestimates the true value by a small amount that is maximum near the center point of the interpolation.
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
The speed of the asm routine is more than fast enough, though spin would be too slow.
OK, how do you do quadratic interpolation? I'll be glad to add in in. Any further accuracy improvements would be icing on the cake.
BTW, the posted code is not in the form I use it. My working code was modified so that others could easily test it. I must be lucky because the range of the numerator just happens to fall in a very nice spot. Like I said earlier, my demoninator is a constant so I don't have to wait on it.
Tom
Would attaching a hardware floating point math co-processor to the propeller solve this problem? Parallax sells one called the uM-FPU V3 for $20 USD, see Parallax part-number 604-00050. There's also a V2 part for $15 USD, see part-number 604-00030; I don't know what the difference is between V2 and V3.
I think these co-processors on offer from Parallax are from MicroMega Corp. See their site at: www.micromegacorp.com
MicroMega now seems to offer a V3.1 part that does not appear on the Parallax site; at least not as of today.
Best Regards,
David
However there are some goodies:
- 2 (medium speed!) 12 bit AD-Converters
- Very extensive set of mathematical routines, including matrix multiplication and a FFT (!), of course, as Badgers put it in a different context: "... memory permitting it..."
Nice toy, but an additional Prop gives you 5x the power for half the price.
However a very god companion for the BasicStamp, or when you just need the algorihms and not much speed
Carefully read the good manual before you buy (I didn't
Post Edited (deSilva) : 8/28/2007 4:57:48 PM GMT
Thanks but no thanks. I am quite happy with the performance of my routine. However the uM-FPU looks like an interesting chip if you need its capabilities.
DeSilva,
Thanks for the comments. "Hau weg die Scheisser!"
Tom
This is b^x with x near 0.5. The peak error at x=1 is twice as much as it is at x=0. I don't know if it is really worthwhile to work it out in order to apply a second order correction, if what you have now is adequate. The additional code might be more than would be required for direct computation with the shift/add methods.
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
Post Edited (Tracy Allen) : 8/28/2007 5:10:12 PM GMT
not a teacher from the southern part of germany because then it's: "Hau weg die Sau!"
I think its more up, maybe Hamburg
Franz
I have been admiring your graphs. (Very nice) What sort of program produces them?
Fred
Hehe. Good ol' Excel charts.
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
My accuracy is more than adequate but I might try it for fun. I'm away from home right now so it will have to wait. If I come up with any code I'll post it in this thread.
Franz, DeSilva,
It's a quote from "Werner". My son taught it to me.
Tom
I've used google docs in the past, you can get to it here: docs.google.com, it's quite nice for making graphs, although not as advanced as Excel. Another office suite you could use is OpenOffice, it's a free, Open Source office "clone". It works quite nicely, I am sure the Excel-like app (calc) does most of everything you need.
Gr,
Mightor
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
| To know recursion, you must first know recursion.
In looking at it more, I realized that the interpolation itself was giving errors of only 0.00000002, which is around 24 bits of precision. A much larger error then is due to quantization. That's good news, because the result can be improved by finer linear interpolation. Instead of 2^5=32 steps, go in whole hog with 2^13=8192 steps. That amounts to 11 bits from the HUB table + 13 bits of linear interpolation, equals 24 total. Only at around 24 bits would a quadratic correction to the linear result make any difference.
Of course, if there are no more than 16 bits of precision going in, it does not make any sense to expect more than 16 out. But with, say, 24 significant bits in a floating point mantissa, it could make a big improvement.
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
Seeing as how I get the remainder through 16 bit division, it wouldn't make sense to go any further. Actually I'm quite happy with the results I'm getting thanks to you and the others.
Tom
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
One comment on your Excel graphs. It sure looks like the absolute value of a sine wave. Is there any way to mathematically relate the two? Even if it isn't a sine wave it looks like a good approximation to one.
Tom
It is a good parabola. Here is a graph of one cycle (red) and the least squares fit to a second order polynomial. R^2 =1.
I don't think there is any direct connection with sine, but who knows? For example, there is a deep connection in terms of periodicity, something that people don't usually think of when looking at the exponential and log funtions, which seem on the surface to be monotonic over their whole range. When calculating the trig functions using CORDIC, the reference table is the arctangent at binary arguments, but for CORDIC exponentials, the reference table is the hyperbolic arctangent.
The graph is from Excel, with the interpolation points calculated to the accuracy that Excel allows (15 significant digits, about 50 binary bits, I think). Wiith the HUB table, the interpolation points are accurate to 16 bits, so there is corresponding quantization error in the end points of the interpolation. So what I said in an earlier post about interpolating to 8192 steps and using a quadratic interpolation is BS. 16 bits with 32 linear interpolation points is about as good as it gets.
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com