log2, C++ (gcc) - rextester

Language: Layout:

#include <cstdio>
#include <cmath>
#include <cstring>
#include <cstdint>
#include <ctime>

float log2_v10(float x)
{
    float y=-127.0f;
    uint32_t bits,bias=0x00400000u;
    memcpy(&bits,&x,sizeof(float));
    if((bits>>23)-1u>253u)
    {
        if(!(x>0.0f)) y=(x==0.0f?-INFINITY:nanf(""));
        else if(x==+INFINITY) y=+INFINITY;
        else {x*=4294967296.0f;memcpy(&bits,&x,sizeof(float));y-=32.0f;}
    }
    bits+=bias;
    y+=float(int32_t(bits>>23));
    bits=((bits&0x007FFFFFu)|0x3F800000u)-bias;
    memcpy(&x,&bits,sizeof(float));
    x=(x-1.0f)/(x+1.0f);
    float x2=x*x;
    return y+x*(2.8853900432586669921875f
           +x2*(0.961808860301971435546875f
           +x2*(0.57604229450225830078125f
           +x2*(0.44434440135955810546875f))));
}

float ulp(float x)
{
    float y;
    uint32_t bits;
    memcpy(&bits,&x,sizeof(float));
    bits^=1;
    memcpy(&y,&bits,sizeof(float));
    return fabsf(x-y);
}

int main()
{
    float constants[]={
        nanf(""),
       -INFINITY, 
      -10.0f,
       -0.0f,
        0.0f,
        1.0e-40f,
        0.07f,
        0.5f,
        1.0f,
        nextafterf(2.0f,1.0f),
        2.0f,
        nextafterf(2.0f,4.0f),
        4.0f,
        5.0f,
        +INFINITY
    };
    for(size_t i=0;i<sizeof(constants)/sizeof(constants[0]);++i)
        printf("log2(%.9g)=%.9g [%.9g]\n",constants[i],log2_v10(constants[i]),log2(constants[i]));
    for(float low=1.0f;low>0.0f;low-=0.5f)
    {
        int n=0;
        double ra=0.0,rr=0.0,ru=0.0;
        double la=0.0,lr=0.0,lu=0.0;
        printf("[%.1f; 2.0):\n",low);
        for(float x=low;x<2.0f;x=nextafterf(x,INFINITY))
        {
            float y=log2_v10(x);
            double z=log2(x);
            double da=(y-z);
            double dr=(y-z)/z;
            if(dr!=dr) dr=0.0f;
            double du=da/ulp(float(z));
            la=fmax(la,fabs(da));
            lr=fmax(lr,fabs(dr));
            lu=fmax(lu,fabs(du));
            ra=ra+da*da;
            rr=rr+dr*dr;
            ru=ru+du*du;
            ++n;
        }
        ra=sqrtf(ra/double(n));
        rr=sqrtf(rr/double(n));
        ru=sqrtf(ru/double(n));
        printf("  Linf: abs=%8.3g, rel=%8.3g, %.3f ulp\n",la,lr,lu);
        printf("  RMS : abs=%8.3g, rel=%8.3g, %.3f ulp\n",ra,rr,ru);
    }
    int N=100000;
    volatile float src[1024]={0.0f};
    volatile float dst[1024]={0.0f};
    for(int i=0;i<1024;++i) src[i]=0.5f+i/1024.0f;
    double t=clock();
    for(int i=0;i<N;++i)
        dst[i&1023]=log2_v10(src[i&1023]);
    t=clock()-t;
    printf("Speed: %.1f ns/call\n",1.0e+9*t/CLOCKS_PER_SEC/N);
    return 0;
}

Show compiler warnings [ + ] Compiler args [ + ] Show input
Compilation time: 0.22 sec, absolute running time: 1.15 sec, cpu time: 1.55 sec, memory peak: 5 Mb, absolute service time: 1,5 sec