algorithm - 浮点除法器硬件实现细节

Question

我正在尝试在硬件中实现一个 32 位浮点硬件除法器，我想知道是否可以就不同算法之间的一些权衡获得任何建议？

我的浮点单元目前支持乘法和加法/减法，但我不打算将其切换到融合乘法加法 (FMA) 浮点架构，因为这是一个嵌入式平台，我试图在其中尽量减少面积使用。

Answer 1

很久以前，我遇到了那个时期军用 FPU 中使用的这种简洁且易于实现的浮点/定点除法算法：

1. 输入必须是无符号和移位的，所以x < y两者都在范围内< 0.5 ; 1 >

   不要忘记存储班次sh = shx - shy和原始标志的差异

2. 找到f（通过迭代）所以y*f -> 1....之后 x*f -> x/y是除法结果

3. x*f向后移动sh并恢复结果符号(sig=sigx*sigy)

   可以x*f像这样轻松计算：

   z=1-y
   (x*f)=(x/y)=x*(1+z)*(1+z^2)*(1+z^4)*(1+z^8)*(1+z^16)...(1+z^2n)
   

   在哪里

   n = log2(num of fractional bits for fixed point, or mantisa bit size for floating point)
   

   您也可以z^2n在固定位宽数据类型为零时停止。

[Edit2] 有一点时间和心情，所以这里是 32 位 IEEE 754 C++ 实现

我删除了旧的（bignum）示例以避免将来的读者混淆（如果需要，它们仍然可以在编辑历史记录中访问）

//---------------------------------------------------------------------------
// IEEE 754 single masks
const DWORD _f32_sig    =0x80000000;    // sign
const DWORD _f32_exp    =0x7F800000;    // exponent
const DWORD _f32_exp_sig=0x40000000;    // exponent sign
const DWORD _f32_exp_bia=0x3F800000;    // exponent bias
const DWORD _f32_exp_lsb=0x00800000;    // exponent LSB
const DWORD _f32_exp_pos=        23;    // exponent LSB bit position
const DWORD _f32_man    =0x007FFFFF;    // mantisa
const DWORD _f32_man_msb=0x00400000;    // mantisa MSB
const DWORD _f32_man_bits=       23;    // mantisa bits
//---------------------------------------------------------------------------
float f32_div(float x,float y)
    {
    union _f32          // float bits access
        {
        float f;        // 32bit floating point
        DWORD u;        // 32 bit uint
        };
    _f32 xx,yy,zz; int sh; DWORD zsig; float z;
    //      result signum        abs value
    xx.f=x; zsig =xx.u&_f32_sig; xx.u&=(0xFFFFFFFF^_f32_sig);
    yy.f=y; zsig^=yy.u&_f32_sig; yy.u&=(0xFFFFFFFF^_f32_sig);
    // initial exponent difference sh and normalize exponents to speed up shift in range
    sh =0;
    sh-=((xx.u&_f32_exp)>>_f32_exp_pos)-(_f32_exp_bia>>_f32_exp_pos); xx.u&=(0xFFFFFFFF^_f32_exp); xx.u|=_f32_exp_bia;
    sh+=((yy.u&_f32_exp)>>_f32_exp_pos)-(_f32_exp_bia>>_f32_exp_pos); yy.u&=(0xFFFFFFFF^_f32_exp); yy.u|=_f32_exp_bia;
    // shift input in range
    while (xx.f> 1.0f) { xx.f*=0.5f; sh--; }
    while (xx.f< 0.5f) { xx.f*=2.0f; sh++; }
    while (yy.f> 1.0f) { yy.f*=0.5f; sh++; }
    while (yy.f< 0.5f) { yy.f*=2.0f; sh--; }
    while (xx.f<=yy.f) { yy.f*=0.5f; sh++; }
    // divider block
    z=(1.0f-yy.f);
    zz.f=xx.f*(1.0f+z);
    for (;;)
        {
        z*=z; if (z==0.0f) break;
        zz.f*=(1.0f+z);
        }
    // shift result back
    for (;sh>0;) { sh--; zz.f*=0.5f; }
    for (;sh<0;) { sh++; zz.f*=2.0f; }
    // set signum
    zz.u&=(0xFFFFFFFF^_f32_sig);
    zz.u|=zsig;
    return zz.f;
    }
//---------------------------------------------------------------------------

我想保持简单，所以它还没有优化。例如，您可以用指数替换所有*=0.5和...如果您与运算符上的 FPU 结果进行比较，这将不太精确，因为大多数 FPU 以 80 位内部格式计算，而此实现仅在 32 位上。*=2.0inc/decfloat/

如您所见，我只是从 FPU 使用+,-,*. 这些东西可以通过使用快速的 sqr 算法来加速，比如

• 快速大数平方计算

特别是如果你想使用大位宽......

不要忘记实施规范化和/或上溢/下溢校正。