Fix mkFixedPointSquareRooter and add safe-guard to mkSquareRooter

[thotypous]

This patch:

• Adds provisos as a safe-guard to mkSquareRooter, forcing the designer to choose an even number of bits for the integer type (UInt#(m)). The implementation works in m/2 iterations and doesn't seem to work for odd values of m.

• Fixes the logic which computes the amount of shift in mkFixedPointSquareRooter, which currently gives results far off from what would be expected (e.g. sqrt(0.5)==0.0 and sqrt(0.25)==0.0).

• Fixes the test bench, which currently displays the runtime computed values as "expected", and compile-time computed values as the operation result. Also, adds two new test values for mkFixedPointSquareRooter.

I don't know if my solution is optimal. I have written it some years ago and submitted it to the forums when Bluespec was proprietary. Suggestions are welcome!

[quark17]

Thank you for this. I apologize that your submission on the forums fell through the cracks. I have some test cases that confirm that mkSquareRooter only works for even sizes and that mkFixedSizeSquareRooter is broken, and that your version does work. I'm not familiar with this design or algorithm, so I can't say if this is the best way to fix it; but it seems reasonable to accept it (with a few cleanups). Since it was broken, I expect that no one was using it, so there's no one to mind.

[quark17]

I showed your changes to someone who is more familiar with the module. He welcomes this fix but thinks it is written in an unnecessarily complicated way. He proposes the following patch (to the original code):

@@ -118,19 +121,23 @@ module mkFixedPointSquareRooter#(Integer n)(Server#(FixedPoint#(isize,fsize),Tup
       UInt#(m) value = unpack({op.i,op.f});
 
       let zeros = countZerosMSB(pack(value));
-      Int#(TAdd#(TLog#(m),1)) shift;
 
-      // align input
-      value = value << (zeros - 1);
-
-      // compute shift for output
-      shift = (fromInteger(valueOf(isize)) - cExtend(zeros)) >> 1;
-      shift = shift + 1;
-      if ((shift & 1) == 0) begin
-        value = value >> 1;
-      end
+      // The decimal place in 'value' is at 'fsize' bit place.
+      // After alignment, it is 'zeros + fsize' bit place.
+      // After square root, it is '(zeros + fsize) / 2' bit place.
+      //
+      // The shift restores the decimal to 'fsize' bit place, thus a
+      // right (left, if negative) shift by '(zeros + fsize) / 2 -
+      // fsize' or '(zeros - fsize) / 2'.
+      //
+      // Since it is not possible to shift by half a bit, adjust
+      // zeros, if necessary, so that shift is an integer.
+      Int#(TAdd#(TLog#(m),2)) shift2 = cExtend(zeros) - fromInteger(valueOf(fsize));
+      Int#(TAdd#(TLog#(m),1)) shift = cExtendLSB(shift2);
+      UInt#(1) lsb = cExtend(shift2);
 
-      shift = fromInteger(valueOf(isize)) - shift;
+      // align input
+      value = value << (zeros - extend(lsb));
 
       sqrt.request.put(value);
       fShift.enq(shift);
@@ -141,7 +148,10 @@ module mkFixedPointSquareRooter#(Integer n)(Server#(FixedPoint#(isize,fsize),Tup
       let shift <- toGet(fShift).get;
 
       // shift result as necessary
-      result = result >> shift;
+      if (shift >= 0)
+         result = result >> shift;
+      else  // invert direction if shift is negative
+         result = result << (-shift);
 
       FixedPoint#(isize,fsize) fx;
       fx.i = cExtendLSB(result);

In place of this proviso that your version needs:

Log#(TAdd#(1, m), TLog#(TAdd#(m, 1)))

this version only needs this proviso:

Add#(a__, 1, TLog#(TAdd#(1, m)))

which I think is just a roundabout way to require that m is non-zero.

What do you think of this version? I have not had a chance to digest these new changes, yet, to understand the difference.

[thotypous]

What do you think of this version? I have not had a chance to digest these new changes, yet, to understand the difference.

Looks like it implements the same logic as mine but in a much less convoluted way. Please say thanks to this person.

These cExtend / cExtendLSB functions are very handful. I should have paid more attention to the BUtils package.

I have force-pushed a new commit integrating the suggested changes to my other fixes.

[quark17]

Thanks again for this. I hacked together some tests for these modules, which I'll add to the testsuite (in bsc.lib/SquareRoot/).

The tests include random inputs on a few instances of each module and a check that each result is within expected error. For mkSquareRooter, that check seems simple. For mkFixedPointSquareRooter, I'm less clear on what the expected error range is -- presumably it's a function of the size parameters of the FixedPoint. For FixedPoint#(48,16), for example, the error isn't just in the fractional part but extends into the whole part -- so my naive attempt at a self-check thinks that the result is too far away. I may still check in my tests, which the current output as the expected output, since it seems reasonable (and a lot better than what the library was computing before).