I'm trying to write some CORDIC routines for log2 and exp2 for single-width floating point input.
The code I see online is usually for trigonometric functions or hypot. I saw an example for log and exp, but I don't want to do a re-scaling to get the result in base 2. Basically it's for a low powered processor and the input x is normalized to just the mantissa 1 <= x < 2 with the exponent handled separately.
In the log CORDIC routine I found, it uses a table of pre-computed values and it says elsewhere that "the algorithm is the same for all bases". For log these pre-computed values are exp((1/2)^n) for each n = 1 ... iteration-limit, so I thought for log2 I would use exp2((1/2)^n), but that gives very wrong answers.
Here's an example I thought up in PHP just to test, it uses the identity log(xy) = log(x) + log(y) and I picked values 1.5, 1.25, 1.125, etc. It's reasonably accurate, kind of looks like a CORDIC routine, and seems to be doing the same number of divisions as the CORDIC log I have. Is this about what one would be doing in a CORDIC log2 function?
<?php
$C = array(
array(1.50000000000000E+00, 5.84962500721156E-01),
array(1.25000000000000E+00, 3.21928094887362E-01),
array(1.12500000000000E+00, 1.69925001442312E-01),
array(1.06250000000000E+00, 8.74628412503394E-02),
array(1.03125000000000E+00, 4.43941193584534E-02),
array(1.01562500000000E+00, 2.23678130284545E-02),
array(1.00781250000000E+00, 1.12272554232541E-02),
array(1.00390625000000E+00, 5.62454919387811E-03),
array(1.00195312500000E+00, 2.81501560705404E-03),
array(1.00097656250000E+00, 1.40819439280839E-03),
array(1.00048828125000E+00, 7.04269011246643E-04),
array(1.00024414062500E+00, 3.52177480301027E-04),
array(1.00012207031250E+00, 1.76099486442506E-04),
array(1.00006103515625E+00, 8.80524301221769E-05),
array(1.00003051757813E+00, 4.40268868273167E-05),
array(1.00001525878906E+00, 2.20136113603405E-05),
array(1.00000762939453E+00, 1.10068476674814E-05),
array(1.00000381469727E+00, 5.50343433064860E-06),
array(1.00000190734863E+00, 2.75171978956128E-06),
array(1.00000095367432E+00, 1.37586055084114E-06),
array(1.00000047683716E+00, 6.87930439435850E-07),
array(1.00000023841858E+00, 3.43965260721765E-07),
array(1.00000011920929E+00, 1.71982640611845E-07),
array(1.00000005960464E+00, 8.59913228686632E-08),
array(1.00000002980232E+00, 4.29956620750169E-08),
array(1.00000001490116E+00, 2.14978311976798E-08),
array(1.00000000745058E+00, 1.07489156388827E-08),
array(1.00000000372529E+00, 5.37445782945206E-09),
array(1.00000000186265E+00, 2.68722891722871E-09),
array(1.00000000093132E+00, 1.34361445924002E-09),
array(1.00000000046566E+00, 6.71807229776429E-10),
array(1.00000000023283E+00, 3.35903614927319E-10),
array(1.00000000011642E+00, 1.67951807473435E-10),
array(1.00000000005821E+00, 8.39759037391618E-11),
array(1.00000000002910E+00, 4.19879518701919E-11),
array(1.00000000001455E+00, 2.09939759352487E-11),
array(1.00000000000728E+00, 1.04969879676625E-11),
array(1.00000000000364E+00, 5.24849398384081E-12),
array(1.00000000000182E+00, 2.62424699192279E-12),
array(1.00000000000091E+00, 1.31212349596199E-12),
array(1.00000000000045E+00, 6.56061747981146E-13),
array(1.00000000000023E+00, 3.28030873990610E-13),
array(1.00000000000011E+00, 1.64015436995314E-13),
array(1.00000000000006E+00, 8.20077184976596E-14),
array(1.00000000000003E+00, 4.10038592488304E-14),
array(1.00000000000001E+00, 2.05019296244153E-14)
);
$y = $x = 1.999;
$r = 0.0;
foreach ($C as $c) {
if ($x >= $c[0]) {
echo '!';
$x /= $c[0];
$r += $c[1];
}
}
echo ' '.$y.' '.$x,' '.$r.' '.log($y)/log(2.0).PHP_EOL;
Output (! for every division):
!!!!!!!!!!!!!!!!!!!!!!!! 1.999 1 0.99927847208251 0.99927847208254
And this is what I came up with for exp2. The first column is equivalent to checking a bit in the mantissa.
<?php
$C = array(
array(1.00000000000000E+00, 2.00000000000000E+00),
array(5.00000000000000E-01, 1.41421356237310E+00),
array(2.50000000000000E-01, 1.18920711500272E+00),
array(1.25000000000000E-01, 1.09050773266526E+00),
array(6.25000000000000E-02, 1.04427378242741E+00),
array(3.12500000000000E-02, 1.02189714865412E+00),
array(1.56250000000000E-02, 1.01088928605170E+00),
array(7.81250000000000E-03, 1.00542990111280E+00),
array(3.90625000000000E-03, 1.00271127505020E+00),
array(1.95312500000000E-03, 1.00135471989211E+00),
array(9.76562500000000E-04, 1.00067713069307E+00),
array(4.88281250000000E-04, 1.00033850805268E+00),
array(2.44140625000000E-04, 1.00016923970530E+00),
array(1.22070312500000E-04, 1.00008461627269E+00),
array(6.10351562500000E-05, 1.00004230724140E+00),
array(3.05175781250000E-05, 1.00002115339696E+00),
array(1.52587890625000E-05, 1.00001057664255E+00),
array(7.62939453125000E-06, 1.00000528830729E+00),
array(3.81469726562500E-06, 1.00000264415015E+00),
array(1.90734863281250E-06, 1.00000132207420E+00),
array(9.53674316406250E-07, 1.00000066103688E+00),
array(4.76837158203125E-07, 1.00000033051839E+00),
array(2.38418579101563E-07, 1.00000016525918E+00),
array(1.19209289550781E-07, 1.00000008262959E+00),
array(5.96046447753906E-08, 1.00000004131479E+00),
array(2.98023223876953E-08, 1.00000002065740E+00),
array(1.49011611938477E-08, 1.00000001032870E+00),
array(7.45058059692383E-09, 1.00000000516435E+00),
array(3.72529029846191E-09, 1.00000000258217E+00),
array(1.86264514923096E-09, 1.00000000129109E+00),
array(9.31322574615479E-10, 1.00000000064554E+00),
array(4.65661287307739E-10, 1.00000000032277E+00),
array(2.32830643653870E-10, 1.00000000016139E+00),
array(1.16415321826935E-10, 1.00000000008069E+00),
array(5.82076609134674E-11, 1.00000000004035E+00),
array(2.91038304567337E-11, 1.00000000002017E+00),
array(1.45519152283669E-11, 1.00000000001009E+00),
array(7.27595761418343E-12, 1.00000000000504E+00),
array(3.63797880709171E-12, 1.00000000000252E+00),
array(1.81898940354586E-12, 1.00000000000126E+00),
array(9.09494701772928E-13, 1.00000000000063E+00),
array(4.54747350886464E-13, 1.00000000000032E+00),
array(2.27373675443232E-13, 1.00000000000016E+00),
array(1.13686837721616E-13, 1.00000000000008E+00),
array(5.68434188608080E-14, 1.00000000000004E+00),
array(2.84217094304040E-14, 1.00000000000002E+00),
array(1.42108547152020E-14, 1.00000000000001E+00)
);
$y = $x = 1.999;
$r = 1.0;
foreach ($C as $c) {
if ($x >= $c[0]) {
echo '!';
$x -= $c[0];
$r *= $c[1];
}
}
echo ' '.$y.' '.$x,' '.$r.' '.pow(2.0, $y).PHP_EOL;
Output (! for every multiplication):
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1.999 4.8849810504606E-15 3.9972283719618 3.9972283719618
Update
So far, have had success using a CORDIC routine for log2 which is accurate to 7.2 decimal digits (i.e. the maximum for IEEE 754 32-bit floating point).
For exp2, no matter what I do: CORDIC, Taylor series about various points (including averaging independent approximations), and polynomial approximations, the best I can get is 6.6 decimal digits of accuracy.