🌱 log2(x)

we have sin(x) and exp2(x) now. only log2(x) is missing. well, our sin(x) is stolen from wikipedia for now. maybe we'll eventually learn how to make our own.

just like with exp2(x), there are binary-level hacks that gets us pretty far. maybe it makes more than a little sense because

x = log2(exp2(x))
x = exp2(log2(x))

last time we got decent results for exp2(x), almost for free, by lerping x into large integers and re-interpreting them as float32 values. this time we're doing it the other way around, by converting the number directly to binary, and then returning the value, almost as-is, after "unlerping" it.

we can also do a stairstep version (the stairstep version was easiest last time?)

if you have trouble understanding what's going on with the bit-shifts, binary-level stuff and such, we went into more detail in the previous article.

in any case, look at the stairstep version and convince yourself that the "missing correction" this time is not a multiplier, but an "adder". that's because the steps are all 1 unit tall today. with exp2(x), each step was twice as tall as the previous one. so last time our solution was:

our_exp2(x) = stairstep_exp2(x) * taylorseries_exp2(x-floor(x))

our_log2(x) = stairstep_log2(x) + taylorseries_log2(...?)

• exp2(x) - steps are 1 unit wide, and twice as tall as the previous one
• log2(x) - steps are 1 unit tall, and twice as wide as the previous one

however, we almost have the answer already for how wide a step is: look at the step that's 4 wide. its y-value is 2. so the width is 1<<2=4. so let's cheat and make a "corrective adder" from log2(x) itself, just so we get a feel for what we have to approximate next.

so we also need an exp2(x), huh (for 2**±stairstep). but the input is always an integer, so we can copy it into the exponent of a float32 and apply bias (so exp2(floor(x)) = u32_to_f32((x+127)<<23)). (read ahead for a much simpler solution to calculate xc. hint: it also involves the other bit-hackery function, u32_to_f32())

on the train without internet i couldn't figure out if maxima has a specific function for log2 (or logn: log in any base), but instead we can just use the math equivalence that log2(x)=log(x)/log(2) (with log(x) being the "natural logarithm", or the base-e logarithm):

(%i1) diff(log(x)/log(2),x,0);
                                    log(x)
(%o1)                               ──────
                                    log(2)
(%i2) diff(log(x)/log(2),x,1);
                                      1
(%o2)                              ────────
                                   log(2) x
(%i3) diff(log(x)/log(2),x,2);
                                        1
(%o3)                             - ─────────
                                            2
                                    log(2) x
(%i4) diff(log(x)/log(2),x,3);
                                       2
(%o4)                              ─────────
                                           3
                                   log(2) x
(%i5) diff(log(x)/log(2),x,4);
                                        6
(%o5)                             - ─────────
                                            4
                                    log(2) x
      

• -(1*1)=-1
• -(-1*2)=2
• -(2*3)=-6
• -(-6*4)=24
• -(24*5)=-120
• ...

note: maxima doesn't give you an error for log2(x). it just assumes it's a function of some kind and gives you the differentiation syntax (in Leibniz notation), like saying "yeah the 3rd derivative of log2(x) is the 3rd derivative of log2(x)":

(%i9) diff(log2(x),x,3);
                                  3
                                 d
(%o9)                            ─── (log2(x))
                                   3
                                 dx
      

less talk more taylor order=3
a=1.500
auto-calibrate

not bad! i get the best results with auto-calibrate on, but it's interesting that it usually has a sweet spot that's not quite a=1.5. for order=3 it's around a=1.473, improving maxerr=0.0016 to maxerr=0.0013. it's also wild how much worse the taylor series performs for a=1 where it hardly converges at all (try turning auto-calibrate off to see how bad it is), and even a=2 is much better. before i started plotting and making tools i couldn't see these things, and just kinda assumed taylor series were terrible for log2(x), and maybe they are a little bit, but not as much as i thought.

the one dough.c "steals" at the time of writing uses division, and taylor series don't (being polynomials). so apparently there are other ways of approximating functions :) i also tried correcting the smoother non-stairstep one from the top, but the error is the same, so i prefer the stairstep one (slightly simpler code)

let's have a look at the "stolen one" using division.

roughly as good as a 6th order above, but for a lot less work. also, the line

const f0 = u32_to_f32((u0 & ((1<<23)-1)) | (126<<23));

const stairstep = (f32_to_u32(x) >> 23) - 127;
const stairexp = u32_to_f32((stairstep+127)<<23);
const invstairexp = u32_to_f32((-stairstep+127)<<23);
const xc = 1+(x-stairexp)*invstairexp;

const xc = u32_to_f32((f32_to_u32(x) & ((1<<23)-1)) | (127<<23));