Digit recurrence square root

Question

I need to implement a digit recurrence square root for generic floating point format such that exp_size + mant_size + 1 <= 64.

I basically followed the implementation suggested here handbook of floating point arithmetic in the software implementation of floating point operator.

I've tried to test my implementation (not an exhaustive test) and basically for format like 32 bit it looks like to work fine, while for format like mantissa = 10, exponent = 5 for the input x = 0.25 instead to give me 0.5 it gives me apparently 0.707031.

So i was wandering if for small format maybe the digit recurrence approach has some limits or not or... simply my implementation is bad...

I hope you can help me... it's a pain to implement this stuff from 0...

Answer 1

it is extremly hard to look at your code but you should:

1. test all the operand combinations

   • if it works for single example does not mean it works for all of them

2. check bit masks

   • you wrote when you use 32bit then result is fine
   • when use 10 then not
   • that is hinting overflow somewhere
   • are you sure you have the right bit counts reserved/masked for R?
   • R should be 2 bits more then Q (+1 bit for accuracy and +1 bit for sign)
   • and also you should handle R as twos complement
   • Q is half of the D bits and unsigned

Could not find your algorithm (that book you linked does not allow me further then page 265 where SQRT starts may be some incompatibility I Use good old Opera) but this is The closest one I found in Google (Non-Restoring-SQRT) in some PDF research and HW implementation on FPGA and after clearing the bugs and testing this is what I code in C++ and tested:

DWORD u32_sqrt(DWORD D) // 32bit
    {
    const int   _bits  =32;
    const DWORD _R_mask=(4<<(_bits>>1))-1;
    const DWORD _R_sign= 2<<(_bits>>1);
    DWORD Q=0; // u(_bits/2    ) result (quotient)
    DWORD R=0; // i(_bits/2 + 2) 2os complement (remainder) R=D-Q*Q
    for (int i=_bits-2;i>=0;i-=2)
        {
        if (DWORD(R&_R_sign)){ R=(R<<2)|((D>>i)&3); R+=(Q<<2)|3; }  // R< 0
         else                { R=(R<<2)|((D>>i)&3); R-=(Q<<2)|1; }  // R>=0
        R&=_R_mask; Q<<=1; if (!DWORD(R&_R_sign)) Q|=1;             // R>=0
        }
    return Q;
    }