cmicrocontrollernatural-logarithm

Issue with implementation of natural logarithm (ln) and exponentiation


I try to follow the topic "Efficient implementation of natural logarithm (ln) and exponentiation" to be able to implement a log function without math.h. The described algorithm works well for values between 1 and 2 (normalized values). However, if the values are not normalized and I follow the normalization instructions, then I get wrong values.

Link: Click here

If I follow the code with an exemplary integer value of 12510, I get the following results:

y = 12510 (0x30DE), log2 = 13, divisor = 26, x = 481,1538

float ln(float y) {
    int log2;
    float divisor, x, result;

    log2 = msb((int)y); // See: https://stackoverflow.com/a/4970859/6630230
    divisor = (float)(1 << log2);
    x = y / divisor;    // normalized value between [1.0, 2.0]

    result = -1.7417939 + (2.8212026 + (-1.4699568 + (0.44717955 - 0.056570851 * x) * x) * x) * x;
    result += ((float)log2) * 0.69314718; // ln(2) = 0.69314718

    return result;
}

The expected result for x should be a normalized value of 1 < x < 2. However, I fail in this calculation, because the received result is 481,1538.

Thanks in advance for any help


Solution

  • Out of curiosity, I tried to reproduce:

    #include <stdio.h>
    
    int msb(unsigned int v) {
      unsigned int r = 0;
      while (v >>= 1) r++;
      return r;
    }
    
    float ln(float y)
    {
        int log2;
        float divisor, x, result;
    
        log2 = msb((int)y); // See: https://stackoverflow.com/a/4970859/6630230
        printf("log2: %d\n", log2);
        divisor = (float)(1 << log2);
        printf("divisor: %f\n", divisor);
        x = y / divisor;    // normalized value between [1.0, 2.0]
        printf("x: %f\n", x);
        result = -1.7417939 + (2.8212026 + (-1.4699568 + (0.44717955 - 0.056570851 * x) * x) * x) * x;
        result += ((float)log2) * 0.69314718; // ln(2) = 0.69314718
        return result;
    }
    
    int main()
    {
      printf("ln(12510): %f\n", ln(12510));
    }
    

    Output:

    log2: 13
    divisor: 8192.000000
    x: 1.527100
    ln(12510): 9.434252
    

    Live Demo on coliru

    I just tried this in my Windows 7 Pocket calculator and got:

    9.434283603460956823997266847405
    

    9,434283603460956823997266847405

    The first 5 digits are identical. – The rest I would consider as rounding issues whithout knowing which one got closer.

    However, there is a typo (or mistake) in the question:

    y = 12510 (0x30DE), log2 = 13, divisor = 26, x = 481,1538

    divisor = (float)(1 << log2); with log2 = 13 yields 8192.

    log2 << 1 would result in 26.

    Just for fun, I changed the line into divisor = (float)(log2 << 1); and got the following output:

    log2: 13
    divisor: 26.000000
    x: 481.153839
    ln(12510): -2982522368.000000
    

    So, this leaves me a bit puzzled:

    The exposed code seems to be correct but the OP seems to have it interpreted (or resembled) wrong.