C++ C 语言中 CORDIC 对数的问题

为了开始使用CORDIC for log10，我实现了在this PDF, pp6中派生的算法:

#include <stdio.h>
#include <math.h>

// https://www.mikrocontroller.net/attachment/31117/cordic1.pdf
float logf_cordic (float x, int B)
{
    float z = 0.0f;

    for (int i = 1; i <= B; i++)
    {
        if (x > 1.0f)
        {
            printf ("-");
            x = x - ldexpf (x, -i);
            z = z - log10f (1.0f - ldexpf (1.0f, -i));
        }
        else
        {
            printf ("+");
            x = x + ldexpf (x, -i);
            z = z - log10f (1.0f + ldexpf (1.0f, -i));
        }
    }
    return z;
}

这是本文中代码的直接C99实现，其中我使用ldexpf(x,n)来计算x·2ⁿ.文中称，该方法对0.4194 < x < 3.4627收敛.该程序使用1 ≤ x ≤ 2.打印以下完整的C99代码:

+--+--+--+-+-+++--++ x = 1.00000: y = -0.000000, dy = -1.487272e-07
-+++++++++++++++++++ x = 1.12500: y = +0.099773, dy = +4.862081e-02
-+++++++++++++++++++ x = 1.25000: y = +0.099773, dy = +2.863325e-03
-+++-+--+-+--+-++--+ x = 1.37500: y = +0.138302, dy = -4.023314e-07
-++-+--++----+++--+- x = 1.50000: y = +0.176091, dy = -2.831221e-07
-+-+++++-++-++-++++- x = 1.62500: y = +0.210854, dy = +2.831221e-07
-+-+-+-+++++--+---++ x = 1.75000: y = +0.243038, dy = +2.235174e-07
-+--++-+--+---+---+- x = 1.87500: y = +0.273001, dy = +0.000000e+00
-+---+--+++--------- x = 2.00000: y = +0.301030, dy = -5.960464e-08

因此，它的工作与预期的一样，除了x = 1.125, 1.25，在那里误差很大，并且在计算更多迭代时不会减少.我现在盯着那个代码看了几个小时，但找不到我错过的东西……

• Code in Compiler Explorer (gcc v11)个

完整的C99代码

#include <stdio.h>
#include <math.h>

float logf_cordic (float x, int B)
{
    float z = 0.0f;

    for (int i = 1; i <= B; i++)
    {
        if (x > 1.0f)
        {
            printf ("-");
            x = x - ldexpf (x, -i);
            z = z - log10f (1.0f - ldexpf (1.0f, -i));
        }
        else
        {
            printf ("+");
            x = x + ldexpf (x, -i);
            z = z - log10f (1.0f + ldexpf (1.0f, -i));
        }
    }
    return z;
}

int main (int argc, char *argv[])
{
    int ex = 3;
    int B = 20;

    if (argc >= 2)
        sscanf (argv[1], "%i", &ex);

    if (argc >= 3)
        sscanf (argv[2], "%i", &B);

    if (ex < 0) ex = 0;
    if (ex > 16) ex = 16;
    if (B > 100) B = 100;

    int n = 1 << ex;
    float dx = 1.0f / n;
    for (int i = 0; i <= n; ++i)
    {
        float x = 1.0f + i * dx;
        float y = logf_cordic (x, B);
        float dy = y - log10f (x);
        printf (" x = %.5f: y = %+f, dy = %+e\n",
                (double) x, (double) y, (double) dy);
    }
    return 0;
}

以下是本文介绍的算法，以供参考:

log10(x){
    z = 0;
    for ( i=1;i=<B;i++ ){
        if (x > 1)
            x = x - x*2^(-i);
            z = z - log10(1-2^(-i));
        else
            x = x + x*2^(-i);
            z = z - log10(1+2^(-i));
    }
return(z)
}