Welcome, spiff! I think it's difficult to compare these two cases, because a 128 bit result demands at least 128 precision in each of the steps, and on an 8/16 bit machine each step is going to take quite a few cycles.

The usual approach is polynomials: rational polynomials (or continued fractions) in the faster implementations, and Chebyshev or some other tweaked Taylor-like series in the slower implementations. A pure Taylor series is, I think, not the wisest tactic - it just happens to be the simplest approach which most people come across first. The derivation of the 'best' coefficients is rather deep, as far as I can tell. We had a bit of a look into this over here.

I'm working on a calculator project and just got done with BCD CORDIC. This site shows how the HP-35 did it, and some guys on the HP forum showed how to implement CORDIC for sine and cosine without having to derive it from an identity involving atan like the first page shows. I ended up saving a lot of cycles this way since division by 10 is easy with BCD. Shifting one bit to divide by two is slow with BCD since you still have to correct an 8 shifted into the MSB of every nibble to a 5. Also, because you reuse every table value ten times, you end up with a smaller table than division by two.

Cool! Please keep us updated.

True. Code size was probably a more important priority than execution speed back then.


Thanks for sharing your insights. I know you have spent some time looking at different implementations.

I will try to look further into this and do some experiments of my own...

Hi!



For derivation of coefficients of rational/polynomial interpolations, I found that http://sollya.gforge.inria.fr/ is very good. It allows to specify the target floating point format in number of bits, so you can use it to generate polynomials appropriate for your specific implementation.

Using it I derived best polynomials for my SIN / COS implementation in FastBasic, see https://github.com/dmsc/fastbasic/blob/ ... os.asm#L38 . As an optimization of code space, I also fixed the last coefficient to exactly PI/2, this works as I'm calculating SIN( PI/2 * y ). It is very fast as it only uses 6 multiplications and produces relative error < 1.23*10^-8

A few years ago I built a hobby robot that was programmed in Forth. My Forth implementation was 16 bit integer only, and I needed trig functions for inverse kinematics. I used brads (binary scaling) to implement pseudo floating point for both CORDIC and lookup tables. CORDIC requires multiple iterations (e.g. 16 or so) to converge to a reasonable approximation, while a lookup is an order one operation. So I ended up using lookup tables with binary search of the lookup tables for inverse trig functions, as it was much much faster.

You only need one quadrant of sine table samples to build all the trig functions, and linear interpolation between samples lets you further compress tables size. I've been meaning to get this code back in shape again for another project as I am using a different version of Forth.