On 6 February 2018 at 16:48, Alex Bennée <alex.ben...@linaro.org> wrote:
> This is a little bit of a departure from softfloat's original approach
> as we skip the estimate step in favour of a straight iteration.
>
> Suggested-by: Richard Henderson <richard.hender...@linaro.org>
> Signed-off-by: Alex Bennée <alex.ben...@linaro.org>
>
> ---
> v3
>   - added to series
>   - fixed renames of structs
> v4
>   - fix up comments
>   - use is_nan
>   - use return_nan instead of pick_nan(a,a)
> ---
>  fpu/softfloat.c         | 201 
> ++++++++++++++++++++++--------------------------
>  include/fpu/softfloat.h |   1 +
>  2 files changed, 91 insertions(+), 111 deletions(-)
>
> diff --git a/fpu/softfloat.c b/fpu/softfloat.c
> index 8fc1c2a8d9..80301d8e04 100644
> --- a/fpu/softfloat.c
> +++ b/fpu/softfloat.c
> @@ -1897,6 +1897,96 @@ float64 float64_scalbn(float64 a, int n, float_status 
> *status)
>      return float64_round_pack_canonical(pr, status);
>  }
>
> +/*
> + * Square Root
> + *
> + * The old softfloat code did an approximation step before zeroing in
> + * on the final result. However for simpleness we just compute the
> + * square root by iterating down from the implicit bit to enough extra
> + * bits to ensure we get a correctly rounded result.

So is this new approach slower?

> + */
> +
> +static FloatParts sqrt_float(FloatParts a, float_status *s, const FloatFmt 
> *p)
> +{
> +    uint64_t a_frac, r_frac, s_frac;
> +    int bit, last_bit;
> +
> +    if (is_nan(a.cls)) {
> +        return return_nan(a, s);
> +    }
> +    if (a.cls == float_class_zero) {
> +        return a;  /* sqrt(+-0) = +-0 */
> +    }
> +    if (a.sign) {
> +        s->float_exception_flags |= float_flag_invalid;
> +        a.cls = float_class_dnan;
> +        return a;
> +    }
> +    if (a.cls == float_class_inf) {
> +        return a;  /* sqrt(+inf) = +inf */
> +    }
> +
> +    assert(a.cls == float_class_normal);
> +
> +    /* We need two overflow bits at the top.  Adding room for that is
> +       a right shift.  If the exponent is odd, we can discard the low
> +       bit by multiplying the fraction by 2; that's a left shift.
> +       Combine those and we shift right if the exponent is even.  */
> +    a_frac = a.frac;
> +    if (!(a.exp & 1)) {
> +        a_frac >>= 1;
> +    }
> +    a.exp >>= 1;

Comment says "shift right if the exponent is even", but code
says "shift right by 1 if exponent is odd, by 2 if exponent is even".

> +
> +    /* Bit-by-bit computation of sqrt.  */
> +    r_frac = 0;
> +    s_frac = 0;
> +
> +    /* Iterate from implicit bit down to the 3 extra bits to compute a
> +     * properly rounded result.  Remember we've inserted one more bit
> +     * at the top, so these positions are one less.  */

Some consistency in multiline comment formatting would be nice.
The top-of-function header and these two in the body of
the function have 3 different styles between them...

Otherwise
Reviewed-by: Peter Maydell <peter.mayd...@linaro.org>

thanks
-- PMM

Reply via email to