Fix mkFixedPointSquareRooter and add safe-guard to mkSquareRooter

Author: thotypous

src/Libraries/Base3-Math/SquareRoot.bsv

@@ -17,6 +17,9 @@ export mkNonPipelinedSquareRooter;
 
 module mkSquareRooter#(Integer n)(Server#(UInt#(m),Tuple2#(UInt#(m),Bool)))
    provisos(
+      Div#(m, 2, m2),
+      // m must be even for this implementation to work
+      Mul#(m2, 2, m),
       // per request of bsc
       Add#(a__, 2, m),
       Log#(TAdd#(1, m), TLog#(TAdd#(m, 1)))
@@ -48,7 +51,7 @@ module mkSquareRooter#(Integer n)(Server#(UInt#(m),Tuple2#(UInt#(m),Bool)))
    FIFO#(Tuple4#(Maybe#(Bit#(m)),Bit#(m),Bit#(m),Bit#(m))) fThis = fFirst;
    FIFO#(Tuple4#(Maybe#(Bit#(m)),Bit#(m),Bit#(m),Bit#(m))) fNext;
 
-   for (Integer i = 0; i < (valueOf(m)/2)/n+1; i = i + 1) begin
+   for (Integer i = 0; i < valueOf(m2)/n+1; i = i + 1) begin
       fNext <- mkLFIFO;
       rule work;
 	 //match {.res, .s, .r, .b} <- toGet(fThis).get;
@@ -59,7 +62,7 @@ module mkSquareRooter#(Integer n)(Server#(UInt#(m),Tuple2#(UInt#(m),Bool)))
 	 Bit#(m) b = tpl_4(x);
 
 	 for (Integer j = 0; j < n; j = j + 1) begin
-	    if ((i + j) <= (valueOf(m)/2)) begin
+	    if ((i + j) <= valueOf(m2)) begin
 	       if (res matches tagged Invalid) begin
 		  if (b == 0) begin
 		     res = tagged Valid r;
@@ -101,6 +104,7 @@ module mkFixedPointSquareRooter#(Integer n)(Server#(FixedPoint#(isize,fsize),Tup
    provisos(
       Add#(isize,fsize,m),
       // per request of bsc
+      Mul#(TDiv#(m, 2), 2, m),
       Add#(a__, 2, m),
       Add#(b__, TLog#(TAdd#(1, m)), TAdd#(TLog#(m), 1)),
       Log#(TAdd#(1, m), TLog#(TAdd#(m, 1))),
@@ -120,17 +124,24 @@ module mkFixedPointSquareRooter#(Integer n)(Server#(FixedPoint#(isize,fsize),Tup
       let zeros = countZerosMSB(pack(value));
       Int#(TAdd#(TLog#(m),1)) shift;
 
+      // if (fsize + zeros) is not even, round zeros down
+      if (valueOf(fsize) % 2 == 0)
+         zeros = zeros & (~1);
+      else
+         zeros = (zeros & 1) == 0 ? zeros - 1 : zeros;
+
       // align input
-      value = value << (zeros - 1);
+      value = value << zeros;
 
-      // compute shift for output
-      shift = (fromInteger(valueOf(isize)) - cExtend(zeros)) >> 1;
-      shift = shift + 1;
-      if ((shift & 1) == 0) begin
-	 value = value >> 1;
-      end
+      // A conversion to UInt behaves like the value was multiplied by 2**fsize.
+      // Then, the value is shifted by zeros, which multiplies it by 2**zeros.
+      // Therefore, after computing the square root of the resulting value,
+      // it must be divided by sqrt(2**(fsize + zeros)) = 2**((fsize + zeros) / 2)
+      shift = extend(unpack(pack( (zeros + fromInteger(valueOf(fsize))) >> 1 )));
 
-      shift = fromInteger(valueOf(isize)) - shift;
+      // Converting the value back from UInt to FixedPoint will already divide it
+      // by 2**fsize, therefore this is excluded from the shift
+      shift = shift - fromInteger(valueOf(fsize));
 
       sqrt.request.put(value);
       fShift.enq(shift);
@@ -141,11 +152,14 @@ 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);
-      fx.f = cExtend(result);
+      fx.i = truncateLSB(pack(result));
+      fx.f = truncate(pack(result));
 
       fResponse.enq(tuple2(fx,inexact));
    endrule
@@ -287,6 +301,8 @@ module mkTb(Empty);
       testSqrtPipe('hffff0000_00000000);
       testSqrtPipe('hfffe0000);
 
+      testSqrtPipe(536870912.9375);
+      testSqrtPipe(1073741824.96875);
       testSqrtPipe(0.5);
       testSqrtPipe(0.25);
 
@@ -307,10 +323,10 @@ module mkTb(Empty);
       end
 
       if (rfx != rrfx) begin
-	 $display("sqrtfx(", fshow(nfx), ") = ", fshow(rfx), " (expected ", fshow(rrfx) ,")");
+	 $display("sqrtfx(", fshow(nfx), ") = ", fshow(rrfx), " (expected ", fshow(rfx) ,")");
       end
       else begin
-	 $display("sqrtfx(", fshow(nfx), ") = ", fshow(rfx));
+	 $display("sqrtfx(", fshow(nfx), ") = ", fshow(rrfx));
       end
    endrule