Patchwork [06/12] RISC-V: Generic <math.h> and soft-fp Routines

login
register
mail settings
Submitter Palmer Dabbelt
Date June 14, 2017, 6:30 p.m.
Message ID <20170614183048.11040-7-palmer@dabbelt.com>
Download mbox | patch
Permalink /patch/21017/
State New
Headers show

Comments

Palmer Dabbelt - June 14, 2017, 6:30 p.m.
This patch contains the miscellaneous math routines and headers we have
implemented for RISC-V.  This includes things from <math.h> that aren't
completely ISA-generic, floating-point bit manipulation, and soft-fp
hooks.
---
 sysdeps/riscv/bits/fenv.h           |  65 ++++++++
 sysdeps/riscv/bits/mathdef.h        |  21 +++
 sysdeps/riscv/bits/mathinline.h     | 142 ++++++++++++++++
 sysdeps/riscv/fpu_control.h         |  77 +++++++++
 sysdeps/riscv/ieee754.h             | 325 ++++++++++++++++++++++++++++++++++++
 sysdeps/riscv/soft-fp/e_sqrtl.c     |   1 +
 sysdeps/riscv/soft-fp/sfp-machine.h |  93 +++++++++++
 7 files changed, 724 insertions(+)
 create mode 100644 sysdeps/riscv/bits/fenv.h
 create mode 100644 sysdeps/riscv/bits/mathdef.h
 create mode 100644 sysdeps/riscv/bits/mathinline.h
 create mode 100644 sysdeps/riscv/fpu_control.h
 create mode 100644 sysdeps/riscv/ieee754.h
 create mode 100644 sysdeps/riscv/soft-fp/e_sqrtl.c
 create mode 100644 sysdeps/riscv/soft-fp/sfp-machine.h
Joseph Myers - June 14, 2017, 8:59 p.m.
On Wed, 14 Jun 2017, Palmer Dabbelt wrote:

> diff --git a/sysdeps/riscv/bits/fenv.h b/sysdeps/riscv/bits/fenv.h

> +enum
> +  {
> +    FE_INEXACT = 0x01,
> +#define FE_INEXACT	FE_INEXACT
> +    FE_UNDERFLOW = 0x02,
> +#define FE_UNDERFLOW	FE_UNDERFLOW
> +    FE_OVERFLOW = 0x04,
> +#define FE_OVERFLOW	FE_OVERFLOW
> +    FE_DIVBYZERO = 0x08,
> +#define FE_DIVBYZERO	FE_DIVBYZERO
> +    FE_INVALID = 0x10,
> +#define FE_INVALID	FE_INVALID
> +  };

This clearly hasn't passed the glibc testsuite.  See how other bits/fenv.h 
headers now define macros to make them usable in #if (which the conform/ 
tests test for).

> diff --git a/sysdeps/riscv/bits/mathdef.h b/sysdeps/riscv/bits/mathdef.h

You shouldn't need this file any more.

> diff --git a/sysdeps/riscv/bits/mathinline.h b/sysdeps/riscv/bits/mathinline.h

> +/* Test for negative number.  Used in the signbit() macro.  */
> +__MATH_INLINE int
> +__NTH (__signbitf (float __x))
> +{
> +  return __builtin_signbitf (__x);
> +}
> +__MATH_INLINE int
> +__NTH (__signbit (double __x))
> +{
> +  return __builtin_signbit (__x);
> +}

Is this really of any use now, given that signbit uses __builtin_signbit 
directly and that the signbit macro, rather than __signbit etc. functions, 
is used for internal calls within glibc?

> +/* Leave it to the compiler to optimize these if __NO_MATH_ERRNO__.  */
> +# if defined __riscv_flen && !defined __NO_MATH_ERRNO__
> +
> +#  if __riscv_xlen >= 64
> +
> +__MATH_INLINE long long int
> +__NTH (llrintf (float __x))
> +{
> +  long long int __res;
> +  __asm__ __volatile__ ("fcvt.l.s %0, %1" : "=r" (__res) : "f" (__x));
> +  return __res;

This all looks rather like making inline functions replicate the bugs of 
out-of-line functions (bug 6797, 6798), i.e. a bad idea.

> diff --git a/sysdeps/riscv/ieee754.h b/sysdeps/riscv/ieee754.h

My reading of the ABI is that long double is 128-bit unconditionally for 
RISC-V.  Thus, you shouldn't need this header; the default 
sysdeps/ieee754/ldbl-128/ieee754.h should suffice.

To pass the math/ tests for soft-float, I'd expect you to need to have 
s_fma.c and s_fmaf.c that include the soft-fp versions, as on various 
other architectures.  (You also need a libm-test-ulps file or files, which 
appear to be missing from this patch series, to reflect the actual ulps 
seen on RISC-V.)
Palmer Dabbelt - June 20, 2017, 9:14 p.m.
On Wed, 14 Jun 2017 13:59:48 PDT (-0700), joseph@codesourcery.com wrote:
> On Wed, 14 Jun 2017, Palmer Dabbelt wrote:
>
>> diff --git a/sysdeps/riscv/bits/fenv.h b/sysdeps/riscv/bits/fenv.h
>
>> +enum
>> +  {
>> +    FE_INEXACT = 0x01,
>> +#define FE_INEXACT	FE_INEXACT
>> +    FE_UNDERFLOW = 0x02,
>> +#define FE_UNDERFLOW	FE_UNDERFLOW
>> +    FE_OVERFLOW = 0x04,
>> +#define FE_OVERFLOW	FE_OVERFLOW
>> +    FE_DIVBYZERO = 0x08,
>> +#define FE_DIVBYZERO	FE_DIVBYZERO
>> +    FE_INVALID = 0x10,
>> +#define FE_INVALID	FE_INVALID
>> +  };
>
> This clearly hasn't passed the glibc testsuite.  See how other bits/fenv.h
> headers now define macros to make them usable in #if (which the conform/
> tests test for).

Yep -- I went ahead and fixed this failure while cleaning up the test suite..
Thanks for the heads up!

>> diff --git a/sysdeps/riscv/bits/mathdef.h b/sysdeps/riscv/bits/mathdef.h
>
> You shouldn't need this file any more.
>
>> diff --git a/sysdeps/riscv/bits/mathinline.h b/sysdeps/riscv/bits/mathinline.h
>
>> +/* Test for negative number.  Used in the signbit() macro.  */
>> +__MATH_INLINE int
>> +__NTH (__signbitf (float __x))
>> +{
>> +  return __builtin_signbitf (__x);
>> +}
>> +__MATH_INLINE int
>> +__NTH (__signbit (double __x))
>> +{
>> +  return __builtin_signbit (__x);
>> +}
>
> Is this really of any use now, given that signbit uses __builtin_signbit
> directly and that the signbit macro, rather than __signbit etc. functions,
> is used for internal calls within glibc?
>
>> +/* Leave it to the compiler to optimize these if __NO_MATH_ERRNO__.  */
>> +# if defined __riscv_flen && !defined __NO_MATH_ERRNO__
>> +
>> +#  if __riscv_xlen >= 64
>> +
>> +__MATH_INLINE long long int
>> +__NTH (llrintf (float __x))
>> +{
>> +  long long int __res;
>> +  __asm__ __volatile__ ("fcvt.l.s %0, %1" : "=r" (__res) : "f" (__x));
>> +  return __res;
>
> This all looks rather like making inline functions replicate the bugs of
> out-of-line functions (bug 6797, 6798), i.e. a bad idea.
>
>> diff --git a/sysdeps/riscv/ieee754.h b/sysdeps/riscv/ieee754.h
>
> My reading of the ABI is that long double is 128-bit unconditionally for
> RISC-V.  Thus, you shouldn't need this header; the default
> sysdeps/ieee754/ldbl-128/ieee754.h should suffice.
>
> To pass the math/ tests for soft-float, I'd expect you to need to have
> s_fma.c and s_fmaf.c that include the soft-fp versions, as on various
> other architectures.  (You also need a libm-test-ulps file or files, which
> appear to be missing from this patch series, to reflect the actual ulps
> seen on RISC-V.)

Thanks for the feedback.  Andrew went through and fixed a bunch of these up.
I'm triaging the test suite to see where we stand, and I'll make sure we
generate the ULPs files.  I believe we need three (soft, single, and double).
Joseph Myers - June 20, 2017, 9:21 p.m.
On Tue, 20 Jun 2017, Palmer Dabbelt wrote:

> Thanks for the feedback.  Andrew went through and fixed a bunch of these up.
> I'm triaging the test suite to see where we stand, and I'll make sure we
> generate the ULPs files.  I believe we need three (soft, single, and double).

Whatever libm-test-ulps files you have, you should have corresponding 
libm-test-ulps-name files with the short description that goes in the 
table of ulps in the manual.
Palmer Dabbelt - June 20, 2017, 9:22 p.m.
On Tue, 20 Jun 2017 14:21:29 PDT (-0700), joseph@codesourcery.com wrote:
> On Tue, 20 Jun 2017, Palmer Dabbelt wrote:
>
>> Thanks for the feedback.  Andrew went through and fixed a bunch of these up.
>> I'm triaging the test suite to see where we stand, and I'll make sure we
>> generate the ULPs files.  I believe we need three (soft, single, and double).
>
> Whatever libm-test-ulps files you have, you should have corresponding
> libm-test-ulps-name files with the short description that goes in the
> table of ulps in the manual.

I saw those too -- thanks for checking, though, as I didn't know about the
manual part.

Patch

diff --git a/sysdeps/riscv/bits/fenv.h b/sysdeps/riscv/bits/fenv.h
new file mode 100644
index 0000000000..a7a8c604ca
--- /dev/null
+++ b/sysdeps/riscv/bits/fenv.h
@@ -0,0 +1,65 @@ 
+/* Copyright (C) 1998, 1999, 2000, 2017 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Lesser General Public
+   License as published by the Free Software Foundation; either
+   version 2.1 of the License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Lesser General Public License for more details.
+
+   You should have received a copy of the GNU Lesser General Public
+   License along with the GNU C Library; if not, write to the Free
+   Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
+   02111-1307 USA.  */
+
+#ifndef _FENV_H
+# error "Never use <bits/fenv.h> directly; include <fenv.h> instead."
+#endif
+
+enum
+  {
+    FE_INEXACT = 0x01,
+#define FE_INEXACT	FE_INEXACT
+    FE_UNDERFLOW = 0x02,
+#define FE_UNDERFLOW	FE_UNDERFLOW
+    FE_OVERFLOW = 0x04,
+#define FE_OVERFLOW	FE_OVERFLOW
+    FE_DIVBYZERO = 0x08,
+#define FE_DIVBYZERO	FE_DIVBYZERO
+    FE_INVALID = 0x10,
+#define FE_INVALID	FE_INVALID
+  };
+
+#define FE_ALL_EXCEPT \
+	(FE_INEXACT | FE_DIVBYZERO | FE_UNDERFLOW | FE_OVERFLOW | FE_INVALID)
+
+enum
+  {
+    FE_TONEAREST = 0x0,
+#define FE_TONEAREST	FE_TONEAREST
+    FE_TOWARDZERO = 0x1,
+#define FE_TOWARDZERO	FE_TOWARDZERO
+    FE_DOWNWARD = 0x2,
+#define FE_DOWNWARD	FE_DOWNWARD
+    FE_UPWARD = 0x3
+#define FE_UPWARD	FE_UPWARD
+  };
+
+
+typedef unsigned int fexcept_t;
+typedef unsigned int fenv_t;
+
+/* If the default argument is used we use this value.  */
+#define FE_DFL_ENV	((__const fenv_t *) -1)
+
+#if __GLIBC_USE (IEC_60559_BFP_EXT)
+/* Type representing floating-point control modes.  */
+typedef unsigned int femode_t;
+
+/* Default floating-point control modes.  */
+# define FE_DFL_MODE	((const femode_t *) -1L)
+#endif
diff --git a/sysdeps/riscv/bits/mathdef.h b/sysdeps/riscv/bits/mathdef.h
new file mode 100644
index 0000000000..c0b9be2bd5
--- /dev/null
+++ b/sysdeps/riscv/bits/mathdef.h
@@ -0,0 +1,21 @@ 
+/* Copyright (C) 2017 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Lesser General Public
+   License as published by the Free Software Foundation; either
+   version 2.1 of the License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Lesser General Public License for more details.
+
+   You should have received a copy of the GNU Lesser General Public
+   License along with the GNU C Library; if not, write to the Free
+   Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
+   02111-1307 USA.  */
+
+#if !defined _MATH_H && !defined _COMPLEX_H
+# error "Never use <bits/mathdef.h> directly; include <math.h> instead"
+#endif
diff --git a/sysdeps/riscv/bits/mathinline.h b/sysdeps/riscv/bits/mathinline.h
new file mode 100644
index 0000000000..9fd009e362
--- /dev/null
+++ b/sysdeps/riscv/bits/mathinline.h
@@ -0,0 +1,142 @@ 
+/* Copyright (C) 2011-2015, 2017 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Chris Metcalf <cmetcalf@tilera.com>, 2011.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Lesser General Public
+   License as published by the Free Software Foundation; either
+   version 2.1 of the License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Lesser General Public License for more details.
+
+   You should have received a copy of the GNU Lesser General Public
+   License along with the GNU C Library.  If not, see
+   <http://www.gnu.org/licenses/>.  */
+
+#ifndef _MATH_H
+# error "Never use <bits/mathinline.h> directly; include <math.h> instead."
+#endif
+
+#ifndef __extern_always_inline
+# define __MATH_INLINE __inline
+#else
+# define __MATH_INLINE __extern_always_inline
+#endif
+
+
+#if defined __USE_ISOC99 && defined __GNUC__
+
+/* Test for negative number.  Used in the signbit() macro.  */
+__MATH_INLINE int
+__NTH (__signbitf (float __x))
+{
+  return __builtin_signbitf (__x);
+}
+__MATH_INLINE int
+__NTH (__signbit (double __x))
+{
+  return __builtin_signbit (__x);
+}
+
+/* Leave it to the compiler to optimize these if __NO_MATH_ERRNO__.  */
+# if defined __riscv_flen && !defined __NO_MATH_ERRNO__
+
+#  if __riscv_xlen >= 64
+
+__MATH_INLINE long long int
+__NTH (llrintf (float __x))
+{
+  long long int __res;
+  __asm__ __volatile__ ("fcvt.l.s %0, %1" : "=r" (__res) : "f" (__x));
+  return __res;
+}
+
+__MATH_INLINE long long int
+__NTH (llroundf (float __x))
+{
+  long long int __res;
+  __asm__ __volatile__ ("fcvt.l.s %0, %1, rmm" : "=r" (__res) : "f" (__x));
+  return __res;
+}
+
+#  endif /* __riscv_xlen >= 64 */
+
+__MATH_INLINE long int
+__NTH (lrintf (float __x))
+{
+#ifdef __LP64__
+  return (long int) llrintf (__x);
+#else
+  long int __res;
+  __asm__ __volatile__ ("fcvt.w.s %0, %1" : "=r" (__res) : "f" (__x));
+  return __res;
+#endif
+}
+
+__MATH_INLINE long int
+__NTH (lroundf (float __x))
+{
+#ifdef __LP64__
+  return (long int) llroundf (__x);
+#else
+  long int __res;
+  __asm__ __volatile__ ("fcvt.w.s %0, %1, rmm" : "=r" (__res) : "f" (__x));
+  return __res;
+#endif
+}
+
+# endif /* __riscv_flen && !__NO_MATH_ERRNO__ */
+
+/* Leave it to the compiler to optimize these if __NO_MATH_ERRNO__.  */
+# if defined __riscv_flen && __riscv_flen >= 64 && !defined __NO_MATH_ERRNO__
+
+#  if __riscv_xlen >= 64
+
+__MATH_INLINE long long int
+__NTH (llrint (double __x))
+{
+  long long int __res;
+  __asm__ __volatile__ ("fcvt.l.d %0, %1" : "=r" (__res) : "f" (__x));
+  return __res;
+}
+
+__MATH_INLINE long long int
+__NTH (llround (double __x))
+{
+  long long int __res;
+  __asm__ __volatile__ ("fcvt.l.d %0, %1, rmm" : "=r" (__res) : "f" (__x));
+  return __res;
+}
+
+#  endif /* __riscv_xlen >= 64 */
+
+__MATH_INLINE long int
+__NTH (lrint (double __x))
+{
+#ifdef __LP64__
+  return (long int) llrint (__x);
+#else
+  long int __res;
+  __asm__ __volatile__ ("fcvt.w.d %0, %1" : "=r" (__res) : "f" (__x));
+  return __res;
+#endif
+}
+
+__MATH_INLINE long int
+__NTH (lround (double __x))
+{
+#ifdef __LP64__
+  return (long int) llround (__x);
+#else
+  long int __res;
+  __asm__ __volatile__ ("fcvt.w.d %0, %1, rmm" : "=r" (__res) : "f" (__x));
+  return __res;
+#endif
+}
+
+# endif /* __riscv_flen >= 64 && !__NO_MATH_ERRNO__ */
+
+#endif
diff --git a/sysdeps/riscv/fpu_control.h b/sysdeps/riscv/fpu_control.h
new file mode 100644
index 0000000000..ca0fa15a0c
--- /dev/null
+++ b/sysdeps/riscv/fpu_control.h
@@ -0,0 +1,77 @@ 
+/* FPU control word bits.  RISC-V version.
+   Copyright (C) 1996, 1997, 1998, 1999, 2000, 2006, 2008, 2017
+   Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Olaf Flebbe and Ralf Baechle.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Lesser General Public
+   License as published by the Free Software Foundation; either
+   version 2.1 of the License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Lesser General Public License for more details.
+
+   You should have received a copy of the GNU Lesser General Public
+   License along with the GNU C Library; if not, write to the Free
+   Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
+   02111-1307 USA.  */
+
+#ifndef _FPU_CONTROL_H
+#define _FPU_CONTROL_H
+
+#include <features.h>
+
+#ifndef __riscv_flen
+
+#define _FPU_RESERVED 0xffffffff
+#define _FPU_DEFAULT  0x00000000
+typedef unsigned int fpu_control_t;
+#define _FPU_GETCW(cw) (cw) = 0
+#define _FPU_SETCW(cw) do { } while (0)
+extern fpu_control_t __fpu_control;
+
+#else /* __riscv_flen */
+
+#define _FPU_RESERVED 0
+#define _FPU_DEFAULT  0
+#define _FPU_IEEE     _FPU_DEFAULT
+
+/* Type of the control word.  */
+typedef unsigned int fpu_control_t __attribute__ ((__mode__ (__SI__)));
+
+/* Macros for accessing the hardware control word.  */
+#define _FPU_GETCW(cw) __asm__ volatile ("frsr %0" : "=r" (cw))
+#define _FPU_SETCW(cw) __asm__ volatile ("fssr %z0" : : "rJ" (cw))
+
+/* Default control word set at startup.  */
+extern fpu_control_t __fpu_control;
+
+#define _FCLASS(x) ({ int res; \
+  if (sizeof (x) * 8 > __riscv_flen) __builtin_trap (); \
+  if (sizeof (x) == 4) asm ("fclass.s %0, %1" : "=r" (res) : "f" (x)); \
+  else if (sizeof (x) == 8) asm ("fclass.d %0, %1" : "=r" (res) : "f" (x)); \
+  else __builtin_trap (); \
+  res; })
+
+#define _FCLASS_MINF     (1<<0)
+#define _FCLASS_MNORM    (1<<1)
+#define _FCLASS_MSUBNORM (1<<2)
+#define _FCLASS_MZERO    (1<<3)
+#define _FCLASS_PZERO    (1<<4)
+#define _FCLASS_PSUBNORM (1<<5)
+#define _FCLASS_PNORM    (1<<6)
+#define _FCLASS_PINF     (1<<7)
+#define _FCLASS_SNAN     (1<<8)
+#define _FCLASS_QNAN     (1<<9)
+#define _FCLASS_ZERO     (_FCLASS_MZERO | _FCLASS_PZERO)
+#define _FCLASS_SUBNORM  (_FCLASS_MSUBNORM | _FCLASS_PSUBNORM)
+#define _FCLASS_NORM     (_FCLASS_MNORM | _FCLASS_PNORM)
+#define _FCLASS_INF      (_FCLASS_MINF | _FCLASS_PINF)
+#define _FCLASS_NAN      (_FCLASS_SNAN | _FCLASS_QNAN)
+
+#endif /* __riscv_flen */
+
+#endif	/* fpu_control.h */
diff --git a/sysdeps/riscv/ieee754.h b/sysdeps/riscv/ieee754.h
new file mode 100644
index 0000000000..f515aa4b36
--- /dev/null
+++ b/sysdeps/riscv/ieee754.h
@@ -0,0 +1,325 @@ 
+/* Copyright (C) 1992, 1995, 1996, 1999, 2002, 2003, 2017
+	Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Lesser General Public
+   License as published by the Free Software Foundation; either
+   version 2.1 of the License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Lesser General Public License for more details.
+
+   You should have received a copy of the GNU Lesser General Public
+   License along with the GNU C Library; if not, write to the Free
+   Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
+   02111-1307 USA.  */
+
+#ifndef _IEEE754_H
+
+#define _IEEE754_H 1
+#include <features.h>
+
+#include <endian.h>
+
+#include <float.h>
+
+__BEGIN_DECLS
+
+union ieee754_float
+  {
+    float f;
+
+    /* This is the IEEE 754 single-precision format.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:8;
+	unsigned int mantissa:23;
+#endif				/* Big endian.  */
+#if	__BYTE_ORDER == __LITTLE_ENDIAN
+	unsigned int mantissa:23;
+	unsigned int exponent:8;
+	unsigned int negative:1;
+#endif				/* Little endian.  */
+      } ieee;
+
+    /* This format makes it easier to see if a NaN is a signalling NaN.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:8;
+	unsigned int quiet_nan:1;
+	unsigned int mantissa:22;
+#endif				/* Big endian.  */
+#if	__BYTE_ORDER == __LITTLE_ENDIAN
+	unsigned int mantissa:22;
+	unsigned int quiet_nan:1;
+	unsigned int exponent:8;
+	unsigned int negative:1;
+#endif				/* Little endian.  */
+      } ieee_nan;
+  };
+
+#define IEEE754_FLOAT_BIAS	0x7f /* Added to exponent.  */
+
+
+union ieee754_double
+  {
+    double d;
+
+    /* This is the IEEE 754 double-precision format.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:11;
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa0:20;
+	unsigned int mantissa1:32;
+#endif				/* Big endian.  */
+#if	__BYTE_ORDER == __LITTLE_ENDIAN
+# if	__FLOAT_WORD_ORDER == __BIG_ENDIAN
+	unsigned int mantissa0:20;
+	unsigned int exponent:11;
+	unsigned int negative:1;
+	unsigned int mantissa1:32;
+# else
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa1:32;
+	unsigned int mantissa0:20;
+	unsigned int exponent:11;
+	unsigned int negative:1;
+# endif
+#endif				/* Little endian.  */
+      } ieee;
+
+    /* This format makes it easier to see if a NaN is a signalling NaN.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:11;
+	unsigned int quiet_nan:1;
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa0:19;
+	unsigned int mantissa1:32;
+#else
+# if	__FLOAT_WORD_ORDER == __BIG_ENDIAN
+	unsigned int mantissa0:19;
+	unsigned int quiet_nan:1;
+	unsigned int exponent:11;
+	unsigned int negative:1;
+	unsigned int mantissa1:32;
+# else
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa1:32;
+	unsigned int mantissa0:19;
+	unsigned int quiet_nan:1;
+	unsigned int exponent:11;
+	unsigned int negative:1;
+# endif
+#endif
+      } ieee_nan;
+  };
+
+#define IEEE754_DOUBLE_BIAS	0x3ff /* Added to exponent.  */
+
+#if LDBL_MANT_DIG == 113
+
+union ieee854_long_double
+  {
+    long double d;
+
+    /* This is the IEEE 854 quad-precision format.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:15;
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa0:16;
+	unsigned int mantissa1:32;
+	unsigned int mantissa2:32;
+	unsigned int mantissa3:32;
+#endif				/* Big endian.  */
+#if	__BYTE_ORDER == __LITTLE_ENDIAN
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa3:32;
+	unsigned int mantissa2:32;
+	unsigned int mantissa1:32;
+	unsigned int mantissa0:16;
+	unsigned int exponent:15;
+	unsigned int negative:1;
+#endif				/* Little endian.  */
+      } ieee;
+
+    /* This format makes it easier to see if a NaN is a signalling NaN.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:15;
+	unsigned int quiet_nan:1;
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa0:15;
+	unsigned int mantissa1:32;
+	unsigned int mantissa2:32;
+	unsigned int mantissa3:32;
+#endif				/* Big endian.  */
+#if	__BYTE_ORDER == __LITTLE_ENDIAN
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa3:32;
+	unsigned int mantissa2:32;
+	unsigned int mantissa1:32;
+	unsigned int mantissa0:15;
+	unsigned int quiet_nan:1;
+	unsigned int exponent:15;
+	unsigned int negative:1;
+#endif				/* Little endian.  */
+      } ieee_nan;
+  };
+
+#define IEEE854_LONG_DOUBLE_BIAS 0x3fff /* Added to exponent.  */
+
+#elif LDBL_MANT_DIG == 64
+
+union ieee854_long_double
+  {
+    long double d;
+
+    /* This is the IEEE 854 double-extended-precision format.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:15;
+	unsigned int empty:16;
+	unsigned int mantissa0:32;
+	unsigned int mantissa1:32;
+#endif
+#if	__BYTE_ORDER == __LITTLE_ENDIAN
+# if	__FLOAT_WORD_ORDER == __BIG_ENDIAN
+	unsigned int exponent:15;
+	unsigned int negative:1;
+	unsigned int empty:16;
+	unsigned int mantissa0:32;
+	unsigned int mantissa1:32;
+# else
+	unsigned int mantissa1:32;
+	unsigned int mantissa0:32;
+	unsigned int exponent:15;
+	unsigned int negative:1;
+	unsigned int empty:16;
+# endif
+#endif
+      } ieee;
+
+    /* This is for NaNs in the IEEE 854 double-extended-precision format.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:15;
+	unsigned int empty:16;
+	unsigned int one:1;
+	unsigned int quiet_nan:1;
+	unsigned int mantissa0:30;
+	unsigned int mantissa1:32;
+#endif
+#if	__BYTE_ORDER == __LITTLE_ENDIAN
+# if	__FLOAT_WORD_ORDER == __BIG_ENDIAN
+	unsigned int exponent:15;
+	unsigned int negative:1;
+	unsigned int empty:16;
+	unsigned int mantissa0:30;
+	unsigned int quiet_nan:1;
+	unsigned int one:1;
+	unsigned int mantissa1:32;
+# else
+	unsigned int mantissa1:32;
+	unsigned int mantissa0:30;
+	unsigned int quiet_nan:1;
+	unsigned int one:1;
+	unsigned int exponent:15;
+	unsigned int negative:1;
+	unsigned int empty:16;
+# endif
+#endif
+      } ieee_nan;
+  };
+
+#define IEEE854_LONG_DOUBLE_BIAS 0x3fff
+
+#elif LDBL_MANT_DIG == 53
+
+union ieee854_long_double
+  {
+    long double d;
+
+    /* This is the IEEE 754 double-precision format.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:11;
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa0:20;
+	unsigned int mantissa1:32;
+#endif				/* Big endian.  */
+#if	__BYTE_ORDER == __LITTLE_ENDIAN
+# if	__FLOAT_WORD_ORDER == __BIG_ENDIAN
+	unsigned int mantissa0:20;
+	unsigned int exponent:11;
+	unsigned int negative:1;
+	unsigned int mantissa1:32;
+# else
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa1:32;
+	unsigned int mantissa0:20;
+	unsigned int exponent:11;
+	unsigned int negative:1;
+# endif
+#endif				/* Little endian.  */
+      } ieee;
+
+    /* This format makes it easier to see if a NaN is a signalling NaN.  */
+    struct
+      {
+#if	__BYTE_ORDER == __BIG_ENDIAN
+	unsigned int negative:1;
+	unsigned int exponent:11;
+	unsigned int quiet_nan:1;
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa0:19;
+	unsigned int mantissa1:32;
+#else
+# if	__FLOAT_WORD_ORDER == __BIG_ENDIAN
+	unsigned int mantissa0:19;
+	unsigned int quiet_nan:1;
+	unsigned int exponent:11;
+	unsigned int negative:1;
+	unsigned int mantissa1:32;
+# else
+	/* Together these comprise the mantissa.  */
+	unsigned int mantissa1:32;
+	unsigned int mantissa0:19;
+	unsigned int quiet_nan:1;
+	unsigned int exponent:11;
+	unsigned int negative:1;
+# endif
+#endif
+      } ieee_nan;
+  };
+
+#define IEEE854_LONG_DOUBLE_BIAS	0x3ff /* Added to exponent.  */
+
+#endif /* LDBL_MANT_DIG == 53 */
+
+__END_DECLS
+
+#endif /* ieee754.h */
diff --git a/sysdeps/riscv/soft-fp/e_sqrtl.c b/sysdeps/riscv/soft-fp/e_sqrtl.c
new file mode 100644
index 0000000000..3241a877ac
--- /dev/null
+++ b/sysdeps/riscv/soft-fp/e_sqrtl.c
@@ -0,0 +1 @@ 
+#include "sysdeps/aarch64/soft-fp/e_sqrtl.c"
diff --git a/sysdeps/riscv/soft-fp/sfp-machine.h b/sysdeps/riscv/soft-fp/sfp-machine.h
new file mode 100644
index 0000000000..3068b5de78
--- /dev/null
+++ b/sysdeps/riscv/soft-fp/sfp-machine.h
@@ -0,0 +1,93 @@ 
+#include <fenv.h>
+#include <fpu_control.h>
+
+#if __riscv_xlen == 32
+
+#define _FP_W_TYPE_SIZE		32
+#define _FP_W_TYPE		unsigned long
+#define _FP_WS_TYPE		signed long
+#define _FP_I_TYPE		long
+
+#define _FP_MUL_MEAT_S(R,X,Y)				\
+  _FP_MUL_MEAT_1_wide(_FP_WFRACBITS_S,R,X,Y,umul_ppmm)
+#define _FP_MUL_MEAT_D(R,X,Y)				\
+  _FP_MUL_MEAT_2_wide(_FP_WFRACBITS_D,R,X,Y,umul_ppmm)
+#define _FP_MUL_MEAT_Q(R,X,Y)				\
+  _FP_MUL_MEAT_4_wide(_FP_WFRACBITS_Q,R,X,Y,umul_ppmm)
+
+#define _FP_DIV_MEAT_S(R,X,Y)	_FP_DIV_MEAT_1_udiv_norm(S,R,X,Y)
+#define _FP_DIV_MEAT_D(R,X,Y)	_FP_DIV_MEAT_2_udiv(D,R,X,Y)
+#define _FP_DIV_MEAT_Q(R,X,Y)	_FP_DIV_MEAT_4_udiv(Q,R,X,Y)
+
+#define _FP_NANFRAC_S		_FP_QNANBIT_S
+#define _FP_NANFRAC_D		_FP_QNANBIT_D, 0
+#define _FP_NANFRAC_Q		_FP_QNANBIT_Q, 0, 0, 0
+
+#else
+
+#define _FP_W_TYPE_SIZE		64
+#define _FP_W_TYPE		unsigned long long
+#define _FP_WS_TYPE		signed long long
+#define _FP_I_TYPE		long long
+
+#define _FP_MUL_MEAT_S(R,X,Y)					\
+  _FP_MUL_MEAT_1_imm(_FP_WFRACBITS_S,R,X,Y)
+#define _FP_MUL_MEAT_D(R,X,Y)					\
+  _FP_MUL_MEAT_1_wide(_FP_WFRACBITS_D,R,X,Y,umul_ppmm)
+#define _FP_MUL_MEAT_Q(R,X,Y)					\
+  _FP_MUL_MEAT_2_wide_3mul(_FP_WFRACBITS_Q,R,X,Y,umul_ppmm)
+
+#define _FP_DIV_MEAT_S(R,X,Y)	_FP_DIV_MEAT_1_imm(S,R,X,Y,_FP_DIV_HELP_imm)
+#define _FP_DIV_MEAT_D(R,X,Y)	_FP_DIV_MEAT_1_udiv_norm(D,R,X,Y)
+#define _FP_DIV_MEAT_Q(R,X,Y)	_FP_DIV_MEAT_2_udiv(Q,R,X,Y)
+
+#define _FP_NANFRAC_S		_FP_QNANBIT_S
+#define _FP_NANFRAC_D		_FP_QNANBIT_D
+#define _FP_NANFRAC_Q		_FP_QNANBIT_Q, 0
+
+#endif
+
+#define _FP_NANSIGN_S		0
+#define _FP_NANSIGN_D		0
+#define _FP_NANSIGN_Q		0
+
+#define _FP_KEEPNANFRACP 0
+#define _FP_QNANNEGATEDP 0
+
+#define _FP_CHOOSENAN(fs, wc, R, X, Y, OP)	\
+  do {						\
+    R##_s = _FP_NANSIGN_##fs;			\
+    _FP_FRAC_SET_##wc(R,_FP_NANFRAC_##fs);	\
+    R##_c = FP_CLS_NAN;				\
+  } while (0)
+
+#define _FP_DECL_EX		int _frm __attribute__ ((unused));
+#define FP_ROUNDMODE		_frm
+
+#define FP_RND_NEAREST		FE_TONEAREST
+#define FP_RND_ZERO		FE_TOWARDZERO
+#define FP_RND_PINF		FE_UPWARD
+#define FP_RND_MINF		FE_DOWNWARD
+
+#define FP_EX_INVALID		FE_INVALID
+#define FP_EX_OVERFLOW		FE_OVERFLOW
+#define FP_EX_UNDERFLOW		FE_UNDERFLOW
+#define FP_EX_DIVZERO		FE_DIVBYZERO
+#define FP_EX_INEXACT		FE_INEXACT
+
+#define _FP_TININESS_AFTER_ROUNDING 1
+
+#ifdef __riscv_flen
+#define FP_INIT_ROUNDMODE			\
+do {						\
+  __asm__ volatile ("frrm %0" : "=r" (_frm));	\
+} while (0)
+
+#define FP_HANDLE_EXCEPTIONS					\
+do {								\
+  if (__builtin_expect (_fex, 0))				\
+    __asm__ volatile ("csrs fflags, %0" : : "rK" (_fex));	\
+} while (0)
+#else
+#define FP_INIT_ROUNDMODE	_frm = FP_RND_NEAREST
+#endif