re PR fortran/29383 (Fortran 2003/F95[TR15580:1999]: Floating point exception (IEEE) support)

	PR fortran/29383

gcc/fortran/
	* gfortran.h (gfc_simplify_ieee_selected_real_kind): New prototype.
	* libgfortran.h (GFC_FPE_*): Use simple integer values, valid in
	both C and Fortran.
	* expr.c (gfc_check_init_expr): Simplify IEEE_SELECTED_REAL_KIND.
	* simplify.c (gfc_simplify_ieee_selected_real_kind): New function.
	* module.c (mio_symbol): Keep track of symbols which came from
	intrinsic modules.
	(gfc_use_module): Keep track of the IEEE modules.
	* trans-decl.c (gfc_get_symbol_decl): Adjust code since
	we have new intrinsic modules.
	(gfc_build_builtin_function_decls): Build decls for
	ieee_procedure_entry and ieee_procedure_exit.
	(is_from_ieee_module, is_ieee_module_used, save_fp_state,
	restore_fp_state): New functions.
	(gfc_generate_function_code): Save and restore floating-point
	state on procedure entry/exit, when IEEE modules are used.
	* intrinsic.texi: Document the IEEE modules.

libgfortran/
	* configure.host: Add checks for IEEE support, rework priorities.
	* configure.ac: Define IEEE_SUPPORT, check for fpsetsticky and
	fpresetsticky.
	* configure: Regenerate.
	* Makefile.am: Build new ieee files, install IEEE_* modules.
	* Makefile.in: Regenerate.
	* gfortran.map (GFORTRAN_1.6): Add new symbols.
	* libgfortran.h (get_fpu_trap_exceptions, set_fpu_trap_exceptions,
	support_fpu_trap, set_fpu_except_flags, support_fpu_flag,
	support_fpu_rounding_mode, get_fpu_state, set_fpu_state): New
	prototypes.
	* config/fpu-*.h (get_fpu_trap_exceptions,
	set_fpu_trap_exceptions, support_fpu_trap, set_fpu_except_flags,
	support_fpu_flag, support_fpu_rounding_mode, get_fpu_state,
	set_fpu_state): New functions.
	* ieee/ieee_features.F90: New file.
	* ieee/ieee_exceptions.F90: New file.
	* ieee/ieee_arithmetic.F90: New file.
	* ieee/ieee_helper.c: New file.

gcc/testsuite/
	* lib/target-supports.exp (check_effective_target_fortran_ieee): 
	New function.
	* gfortran.dg/ieee/ieee.exp: New file.
	* gfortran.dg/ieee/ieee_1.F90: New file.
	* gfortran.dg/ieee/ieee_2.f90: New file.
	* gfortran.dg/ieee/ieee_3.f90: New file.
	* gfortran.dg/ieee/ieee_4.f90: New file.
	* gfortran.dg/ieee/ieee_5.f90: New file.
	* gfortran.dg/ieee/ieee_6.f90: New file.
	* gfortran.dg/ieee/ieee_7.f90: New file.
	* gfortran.dg/ieee/ieee_rounding_1.f90: New file.

From-SVN: r212102

5 files changed, 1066 insertions, 89 deletions

diff --git a/libgfortran/config/fpu-387.h b/libgfortran/config/fpu-387.h
index 7b56293..46720b2 100644
--- a/libgfortran/config/fpu-387.h
+++ b/libgfortran/config/fpu-387.h
@@ -23,6 +23,8 @@ a copy of the GCC Runtime Library Exception along with this program;
 see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
 <http://www.gnu.org/licenses/>.  */
 
+#include <assert.h>
+
 #ifndef __SSE_MATH__
 #include "cpuid.h"
 #endif
@@ -62,24 +64,122 @@ has_sse (void)
 
 #define _FPU_RC_MASK    0x3
 
+/* This structure corresponds to the layout of the block
+   written by FSTENV.  */
+typedef struct
+{
+  unsigned short int __control_word;
+  unsigned short int __unused1;
+  unsigned short int __status_word;
+  unsigned short int __unused2;
+  unsigned short int __tags;
+  unsigned short int __unused3;
+  unsigned int __eip;
+  unsigned short int __cs_selector;
+  unsigned int __opcode:11;
+  unsigned int __unused4:5;
+  unsigned int __data_offset;
+  unsigned short int __data_selector;
+  unsigned short int __unused5;
+  unsigned int __mxcsr;
+}
+my_fenv_t;
+
+
+/* Raise the supported floating-point exceptions from EXCEPTS.  Other
+   bits in EXCEPTS are ignored.  Code originally borrowed from
+   libatomic/config/x86/fenv.c.  */
+
+static void
+local_feraiseexcept (int excepts)
+{
+  if (excepts & _FPU_MASK_IM)
+    {
+      float f = 0.0f;
+#ifdef __SSE_MATH__
+      volatile float r __attribute__ ((unused));
+      __asm__ __volatile__ ("%vdivss\t{%0, %d0|%d0, %0}" : "+x" (f));
+      r = f; /* Needed to trigger exception.   */
+#else
+      __asm__ __volatile__ ("fdiv\t{%y0, %0|%0, %y0}" : "+t" (f));
+      /* No need for fwait, exception is triggered by emitted fstp.  */
+#endif
+    }
+  if (excepts & _FPU_MASK_DM)
+    {
+      my_fenv_t temp;
+      __asm__ __volatile__ ("fnstenv\t%0" : "=m" (temp));
+      temp.__status_word |= _FPU_MASK_DM;
+      __asm__ __volatile__ ("fldenv\t%0" : : "m" (temp));
+      __asm__ __volatile__ ("fwait");
+    }
+  if (excepts & _FPU_MASK_ZM)
+    {
+      float f = 1.0f, g = 0.0f;
+#ifdef __SSE_MATH__
+      volatile float r __attribute__ ((unused));
+      __asm__ __volatile__ ("%vdivss\t{%1, %d0|%d0, %1}" : "+x" (f) : "xm" (g));
+      r = f; /* Needed to trigger exception.   */
+#else
+      __asm__ __volatile__ ("fdivs\t%1" : "+t" (f) : "m" (g));
+      /* No need for fwait, exception is triggered by emitted fstp.  */
+#endif
+    }
+  if (excepts & _FPU_MASK_OM)
+    {
+      my_fenv_t temp;
+      __asm__ __volatile__ ("fnstenv\t%0" : "=m" (temp));
+      temp.__status_word |= _FPU_MASK_OM;
+      __asm__ __volatile__ ("fldenv\t%0" : : "m" (temp));
+      __asm__ __volatile__ ("fwait");
+    }
+  if (excepts & _FPU_MASK_UM)
+    {
+      my_fenv_t temp;
+      __asm__ __volatile__ ("fnstenv\t%0" : "=m" (temp));
+      temp.__status_word |= _FPU_MASK_UM;
+      __asm__ __volatile__ ("fldenv\t%0" : : "m" (temp));
+      __asm__ __volatile__ ("fwait");
+    }
+  if (excepts & _FPU_MASK_PM)
+    {
+      float f = 1.0f, g = 3.0f;
+#ifdef __SSE_MATH__
+      volatile float r __attribute__ ((unused));
+      __asm__ __volatile__ ("%vdivss\t{%1, %d0|%d0, %1}" : "+x" (f) : "xm" (g));
+      r = f; /* Needed to trigger exception.   */
+#else
+      __asm__ __volatile__ ("fdivs\t%1" : "+t" (f) : "m" (g));
+      /* No need for fwait, exception is triggered by emitted fstp.  */
+#endif
+    }
+}
+
 
 void
-set_fpu (void)
+set_fpu_trap_exceptions (int trap, int notrap)
 {
-  int excepts = 0;
+  int exc_set = 0, exc_clr = 0;
   unsigned short cw;
 
-  __asm__ __volatile__ ("fstcw\t%0" : "=m" (cw));
+  if (trap & GFC_FPE_INVALID) exc_set |= _FPU_MASK_IM;
+  if (trap & GFC_FPE_DENORMAL) exc_set |= _FPU_MASK_DM;
+  if (trap & GFC_FPE_ZERO) exc_set |= _FPU_MASK_ZM;
+  if (trap & GFC_FPE_OVERFLOW) exc_set |= _FPU_MASK_OM;
+  if (trap & GFC_FPE_UNDERFLOW) exc_set |= _FPU_MASK_UM;
+  if (trap & GFC_FPE_INEXACT) exc_set |= _FPU_MASK_PM;
+
+  if (notrap & GFC_FPE_INVALID) exc_clr |= _FPU_MASK_IM;
+  if (notrap & GFC_FPE_DENORMAL) exc_clr |= _FPU_MASK_DM;
+  if (notrap & GFC_FPE_ZERO) exc_clr |= _FPU_MASK_ZM;
+  if (notrap & GFC_FPE_OVERFLOW) exc_clr |= _FPU_MASK_OM;
+  if (notrap & GFC_FPE_UNDERFLOW) exc_clr |= _FPU_MASK_UM;
+  if (notrap & GFC_FPE_INEXACT) exc_clr |= _FPU_MASK_PM;
 
-  if (options.fpe & GFC_FPE_INVALID) excepts |= _FPU_MASK_IM;
-  if (options.fpe & GFC_FPE_DENORMAL) excepts |= _FPU_MASK_DM;
-  if (options.fpe & GFC_FPE_ZERO) excepts |= _FPU_MASK_ZM;
-  if (options.fpe & GFC_FPE_OVERFLOW) excepts |= _FPU_MASK_OM;
-  if (options.fpe & GFC_FPE_UNDERFLOW) excepts |= _FPU_MASK_UM;
-  if (options.fpe & GFC_FPE_INEXACT) excepts |= _FPU_MASK_PM;
+  __asm__ __volatile__ ("fstcw\t%0" : "=m" (cw));
 
-  cw |= _FPU_MASK_ALL;
-  cw &= ~excepts;
+  cw |= exc_clr;
+  cw &= ~exc_set;
 
   __asm__ __volatile__ ("fnclex\n\tfldcw\t%0" : : "m" (cw));
 
@@ -90,8 +190,8 @@ set_fpu (void)
       __asm__ __volatile__ ("%vstmxcsr\t%0" : "=m" (cw_sse));
 
       /* The SSE exception masks are shifted by 7 bits.  */
-      cw_sse |= _FPU_MASK_ALL << 7;
-      cw_sse &= ~(excepts << 7);
+      cw_sse |= (exc_clr << 7);
+      cw_sse &= ~(exc_set << 7);
 
       /* Clear stalled exception flags.  */
       cw_sse &= ~_FPU_EX_ALL;
@@ -100,6 +200,47 @@ set_fpu (void)
     }
 }
 
+void
+set_fpu (void)
+{
+  set_fpu_trap_exceptions (options.fpe, 0);
+}
+
+int
+get_fpu_trap_exceptions (void)
+{
+  int res = 0;
+  unsigned short cw;
+
+  __asm__ __volatile__ ("fstcw\t%0" : "=m" (cw));
+  cw &= _FPU_MASK_ALL;
+
+  if (has_sse())
+    {
+      unsigned int cw_sse;
+
+      __asm__ __volatile__ ("%vstmxcsr\t%0" : "=m" (cw_sse));
+
+      /* The SSE exception masks are shifted by 7 bits.  */
+      cw = cw | ((cw_sse >> 7) & _FPU_MASK_ALL);
+    }
+
+  if (~cw & _FPU_MASK_IM) res |= GFC_FPE_INVALID;
+  if (~cw & _FPU_MASK_DM) res |= GFC_FPE_DENORMAL;
+  if (~cw & _FPU_MASK_ZM) res |= GFC_FPE_ZERO;
+  if (~cw & _FPU_MASK_OM) res |= GFC_FPE_OVERFLOW;
+  if (~cw & _FPU_MASK_UM) res |= GFC_FPE_UNDERFLOW;
+  if (~cw & _FPU_MASK_PM) res |= GFC_FPE_INEXACT;
+
+  return res;
+}
+
+int
+support_fpu_trap (int flag __attribute__((unused)))
+{
+  return 1;
+}
+
 int
 get_fpu_except_flags (void)
 {
@@ -107,7 +248,7 @@ get_fpu_except_flags (void)
   int excepts;
   int result = 0;
 
-  __asm__ __volatile__ ("fnstsw\t%0" : "=a" (cw));
+  __asm__ __volatile__ ("fnstsw\t%0" : "=am" (cw));
   excepts = cw;
 
   if (has_sse())
@@ -131,6 +272,70 @@ get_fpu_except_flags (void)
 }
 
 void
+set_fpu_except_flags (int set, int clear)
+{
+  my_fenv_t temp;
+  int exc_set = 0, exc_clr = 0;
+
+  /* Translate from GFC_PE_* values to _FPU_MASK_* values.  */
+  if (set & GFC_FPE_INVALID)
+    exc_set |= _FPU_MASK_IM;
+  if (clear & GFC_FPE_INVALID)
+    exc_clr |= _FPU_MASK_IM;
+
+  if (set & GFC_FPE_DENORMAL)
+    exc_set |= _FPU_MASK_DM;
+  if (clear & GFC_FPE_DENORMAL)
+    exc_clr |= _FPU_MASK_DM;
+
+  if (set & GFC_FPE_ZERO)
+    exc_set |= _FPU_MASK_ZM;
+  if (clear & GFC_FPE_ZERO)
+    exc_clr |= _FPU_MASK_ZM;
+
+  if (set & GFC_FPE_OVERFLOW)
+    exc_set |= _FPU_MASK_OM;
+  if (clear & GFC_FPE_OVERFLOW)
+    exc_clr |= _FPU_MASK_OM;
+
+  if (set & GFC_FPE_UNDERFLOW)
+    exc_set |= _FPU_MASK_UM;
+  if (clear & GFC_FPE_UNDERFLOW)
+    exc_clr |= _FPU_MASK_UM;
+
+  if (set & GFC_FPE_INEXACT)
+    exc_set |= _FPU_MASK_PM;
+  if (clear & GFC_FPE_INEXACT)
+    exc_clr |= _FPU_MASK_PM;
+
+
+  /* Change the flags. This is tricky on 387 (unlike SSE), because we have
+     FNSTSW but no FLDSW instruction.  */
+  __asm__ __volatile__ ("fnstenv\t%0" : "=m" (temp));
+  temp.__status_word &= ~exc_clr;
+  __asm__ __volatile__ ("fldenv\t%0" : : "m" (temp));
+
+  /* Change the flags on SSE.  */
+
+  if (has_sse())
+  {
+    unsigned int cw_sse;
+
+    __asm__ __volatile__ ("%vstmxcsr\t%0" : "=m" (cw_sse));
+    cw_sse &= ~exc_clr;
+    __asm__ __volatile__ ("%vldmxcsr\t%0" : : "m" (cw_sse));
+  }
+
+  local_feraiseexcept (exc_set);
+}
+
+int
+support_fpu_flag (int flag __attribute__((unused)))
+{
+  return 1;
+}
+
+void
 set_fpu_rounding_mode (int round)
 {
   int round_mode;
@@ -213,3 +418,44 @@ get_fpu_rounding_mode (void)
       return GFC_FPE_INVALID; /* Should be unreachable.  */
     }
 }
+
+int
+support_fpu_rounding_mode (int mode __attribute__((unused)))
+{
+  return 1;
+}
+
+void
+get_fpu_state (void *state)
+{
+  my_fenv_t *envp = state;
+
+  /* Check we can actually store the FPU state in the allocated size.  */
+  assert (sizeof(my_fenv_t) <= (size_t) GFC_FPE_STATE_BUFFER_SIZE);
+
+  __asm__ __volatile__ ("fnstenv\t%0" : "=m" (*envp));
+
+  /* fnstenv has the side effect of masking all exceptions, so we need
+     to restore the control word after that.  */
+  __asm__ __volatile__ ("fldcw\t%0" : : "m" (envp->__control_word));
+
+  if (has_sse())
+    __asm__ __volatile__ ("%vstmxcsr\t%0" : "=m" (envp->__mxcsr));
+}
+
+void
+set_fpu_state (void *state)
+{
+  my_fenv_t *envp = state;
+
+  /* Check we can actually store the FPU state in the allocated size.  */
+  assert (sizeof(my_fenv_t) <= (size_t) GFC_FPE_STATE_BUFFER_SIZE);
+
+  /* glibc sources (sysdeps/x86_64/fpu/fesetenv.c) do something more
+     complex than this, but I think it suffices in our case.  */
+  __asm__ __volatile__ ("fldenv\t%0" : : "m" (*envp));
+
+  if (has_sse())
+    __asm__ __volatile__ ("%vldmxcsr\t%0" : : "m" (envp->__mxcsr));
+}
+
diff --git a/libgfortran/config/fpu-aix.h b/libgfortran/config/fpu-aix.h
index a05fab8..6b44ab7 100644
--- a/libgfortran/config/fpu-aix.h
+++ b/libgfortran/config/fpu-aix.h
@@ -33,15 +33,103 @@ see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
 #include <fpxcp.h>
 #endif
 
+#ifdef HAVE_FENV_H
+#include <fenv.h>
+#endif
+
+
 void
-set_fpu (void)
+set_fpu_trap_exceptions (int trap, int notrap)
 {
-  fptrap_t mode = 0;
+  fptrap_t mode_set = 0, mode_clr = 0;
 
-  if (options.fpe & GFC_FPE_INVALID)
 #ifdef TRP_INVALID
-    mode |= TRP_INVALID;
-#else
+  if (trap & GFC_FPE_INVALID)
+    mode_set |= TRP_INVALID;
+  if (notrap & GFC_FPE_INVALID)
+    mode_clr |= TRP_INVALID;
+#endif
+
+#ifdef TRP_DIV_BY_ZERO
+  if (trap & GFC_FPE_ZERO)
+    mode_set |= TRP_DIV_BY_ZERO;
+  if (notrap & GFC_FPE_ZERO)
+    mode_clr |= TRP_DIV_BY_ZERO;
+#endif
+
+#ifdef TRP_OVERFLOW
+  if (trap & GFC_FPE_OVERFLOW)
+    mode_set |= TRP_OVERFLOW;
+  if (notrap & GFC_FPE_OVERFLOW)
+    mode_clr |= TRP_OVERFLOW;
+#endif
+
+#ifdef TRP_UNDERFLOW
+  if (trap & GFC_FPE_UNDERFLOW)
+    mode_set |= TRP_UNDERFLOW;
+  if (notrap & GFC_FPE_UNDERFLOW)
+    mode_clr |= TRP_UNDERFLOW;
+#endif
+
+#ifdef TRP_INEXACT
+  if (trap & GFC_FPE_INEXACT)
+    mode_set |= TRP_INEXACT;
+  if (notrap & GFC_FPE_INEXACT)
+    mode_clr |= TRP_INEXACT;
+#endif
+
+  fp_trap (FP_TRAP_SYNC);
+  fp_enable (mode_set);
+  fp_disable (mode_clr);
+}
+
+
+int
+get_fpu_trap_exceptions (void)
+{
+  int res = 0;
+
+#ifdef TRP_INVALID
+  if (fp_is_enabled (TRP_INVALID))
+    res |= GFC_FPE_INVALID;
+#endif
+
+#ifdef TRP_DIV_BY_ZERO
+  if (fp_is_enabled (TRP_DIV_BY_ZERO))
+    res |= GFC_FPE_ZERO;
+#endif
+
+#ifdef TRP_OVERFLOW
+  if (fp_is_enabled (TRP_OVERFLOW))
+    res |= GFC_FPE_OVERFLOW;
+#endif
+
+#ifdef TRP_UNDERFLOW
+  if (fp_is_enabled (TRP_UNDERFLOW))
+    res |= GFC_FPE_UNDERFLOW;
+#endif
+
+#ifdef TRP_INEXACT
+  if (fp_is_enabled (TRP_INEXACT))
+    res |= GFC_FPE_INEXACT;
+#endif
+
+  return res;
+}
+
+
+int
+support_fpu_trap (int flag)
+{
+  return support_fpu_flag (flag);
+}
+
+
+void
+set_fpu (void)
+{
+#ifndef TRP_INVALID
+  if (options.fpe & GFC_FPE_INVALID)
     estr_write ("Fortran runtime warning: IEEE 'invalid operation' "
 	        "exception not supported.\n");
 #endif
@@ -50,43 +138,33 @@ set_fpu (void)
     estr_write ("Fortran runtime warning: Floating point 'denormal operand' "
 	        "exception not supported.\n");
 
+#ifndef TRP_DIV_BY_ZERO
   if (options.fpe & GFC_FPE_ZERO)
-#ifdef TRP_DIV_BY_ZERO
-    mode |= TRP_DIV_BY_ZERO;
-#else
     estr_write ("Fortran runtime warning: IEEE 'division by zero' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef TRP_OVERFLOW
   if (options.fpe & GFC_FPE_OVERFLOW)
-#ifdef TRP_OVERFLOW
-    mode |= TRP_OVERFLOW;
-#else
     estr_write ("Fortran runtime warning: IEEE 'overflow' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef TRP_UNDERFLOW
   if (options.fpe & GFC_FPE_UNDERFLOW)
-#ifdef TRP_UNDERFLOW
-    mode |= TRP_UNDERFLOW;
-#else
     estr_write ("Fortran runtime warning: IEEE 'underflow' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef TRP_INEXACT
   if (options.fpe & GFC_FPE_INEXACT)
-#ifdef TRP_INEXACT
-    mode |= TRP_INEXACT;
-#else
     estr_write ("Fortran runtime warning: IEEE 'inexact' "
 	        "exception not supported.\n");
 #endif
 
-  fp_trap(FP_TRAP_SYNC);
-  fp_enable(mode);
+  set_fpu_trap_exceptions (options.fpe, 0);
 }
 
-
 int
 get_fpu_except_flags (void)
 {
@@ -118,6 +196,98 @@ get_fpu_except_flags (void)
 }
 
 
+void
+set_fpu_except_flags (int set, int clear)
+{
+  int exc_set = 0, exc_clr = 0;
+
+#ifdef FP_INVALID
+  if (set & GFC_FPE_INVALID)
+    exc_set |= FP_INVALID;
+  else if (clear & GFC_FPE_INVALID)
+    exc_clr |= FP_INVALID;
+#endif
+
+#ifdef FP_DIV_BY_ZERO
+  if (set & GFC_FPE_ZERO)
+    exc_set |= FP_DIV_BY_ZERO;
+  else if (clear & GFC_FPE_ZERO)
+    exc_clr |= FP_DIV_BY_ZERO;
+#endif
+
+#ifdef FP_OVERFLOW
+  if (set & GFC_FPE_OVERFLOW)
+    exc_set |= FP_OVERFLOW;
+  else if (clear & GFC_FPE_OVERFLOW)
+    exc_clr |= FP_OVERFLOW;
+#endif
+
+#ifdef FP_UNDERFLOW
+  if (set & GFC_FPE_UNDERFLOW)
+    exc_set |= FP_UNDERFLOW;
+  else if (clear & GFC_FPE_UNDERFLOW)
+    exc_clr |= FP_UNDERFLOW;
+#endif
+
+/* AIX does not have FP_DENORMAL.  */
+
+#ifdef FP_INEXACT
+  if (set & GFC_FPE_INEXACT)
+    exc_set |= FP_INEXACT;
+  else if (clear & GFC_FPE_INEXACT)
+    exc_clr |= FP_INEXACT;
+#endif
+
+  fp_clr_flag (exc_clr);
+  fp_set_flag (exc_set);
+}
+
+
+int
+support_fpu_flag (int flag)
+{
+  if (flag & GFC_FPE_INVALID)
+  {
+#ifndef FP_INVALID
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_ZERO)
+  {
+#ifndef FP_DIV_BY_ZERO
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_OVERFLOW)
+  {
+#ifndef FP_OVERFLOW
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_UNDERFLOW)
+  {
+#ifndef FP_UNDERFLOW
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_DENORMAL)
+  {
+    /* AIX does not support denormal flag.  */
+    return 0;
+  }
+  else if (flag & GFC_FPE_INEXACT)
+  {
+#ifndef FP_INEXACT
+    return 0;
+#endif
+  }
+
+  return 1;
+}
+
+
+
+
 int
 get_fpu_rounding_mode (void)
 {
@@ -188,3 +358,60 @@ set_fpu_rounding_mode (int mode)
 
   fesetround (rnd_mode);
 }
+
+
+int
+support_fpu_rounding_mode (int mode)
+{
+  switch (mode)
+    {
+      case GFC_FPE_TONEAREST:
+#ifdef FE_TONEAREST
+	return 1;
+#else
+	return 0;
+#endif
+
+#ifdef FE_UPWARD
+	return 1;
+#else
+	return 0;
+#endif
+
+#ifdef FE_DOWNWARD
+	return 1;
+#else
+	return 0;
+#endif
+
+#ifdef FE_TOWARDZERO
+	return 1;
+#else
+	return 0;
+#endif
+
+      default:
+	return 0;
+    }
+}
+
+
+
+void
+get_fpu_state (void *state)
+{
+  /* Check we can actually store the FPU state in the allocated size.  */
+  assert (sizeof(fenv_t) <= GFC_FPE_STATE_BUFFER_SIZE);
+
+  fegetenv (state);
+}
+
+void
+set_fpu_state (void *state)
+{
+  /* Check we can actually store the FPU state in the allocated size.  */
+  assert (sizeof(fenv_t) <= GFC_FPE_STATE_BUFFER_SIZE);
+
+  fesetenv (state);
+}
+
diff --git a/libgfortran/config/fpu-generic.h b/libgfortran/config/fpu-generic.h
index d9be4d9..bbad875 100644
--- a/libgfortran/config/fpu-generic.h
+++ b/libgfortran/config/fpu-generic.h
@@ -51,6 +51,12 @@ set_fpu (void)
 	        "exception not supported.\n");
 }
 
+void
+set_fpu_trap_exceptions (int trap __attribute__((unused)),
+			 int notrap __attribute__((unused)))
+{
+}
+
 int
 get_fpu_except_flags (void)
 {
diff --git a/libgfortran/config/fpu-glibc.h b/libgfortran/config/fpu-glibc.h
index cf21684..695b9d3 100644
--- a/libgfortran/config/fpu-glibc.h
+++ b/libgfortran/config/fpu-glibc.h
@@ -27,63 +27,141 @@ see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
    feenableexcept function in fenv.h to set individual exceptions
    (there's nothing to do that in C99).  */
 
+#include <assert.h>
+
 #ifdef HAVE_FENV_H
 #include <fenv.h>
 #endif
 
-void set_fpu (void)
-{
-  if (FE_ALL_EXCEPT != 0)
-    fedisableexcept (FE_ALL_EXCEPT);
 
-  if (options.fpe & GFC_FPE_INVALID)
+void set_fpu_trap_exceptions (int trap, int notrap)
+{
 #ifdef FE_INVALID
+  if (trap & GFC_FPE_INVALID)
     feenableexcept (FE_INVALID);
-#else
+  if (notrap & GFC_FPE_INVALID)
+    fedisableexcept (FE_INVALID);
+#endif
+
+/* glibc does never have a FE_DENORMAL.  */
+#ifdef FE_DENORMAL
+  if (trap & GFC_FPE_DENORMAL)
+    feenableexcept (FE_DENORMAL);
+  if (notrap & GFC_FPE_DENORMAL)
+    fedisableexcept (FE_DENORMAL);
+#endif
+
+#ifdef FE_DIVBYZERO
+  if (trap & GFC_FPE_ZERO)
+    feenableexcept (FE_DIVBYZERO);
+  if (notrap & GFC_FPE_ZERO)
+    fedisableexcept (FE_DIVBYZERO);
+#endif
+
+#ifdef FE_OVERFLOW
+  if (trap & GFC_FPE_OVERFLOW)
+    feenableexcept (FE_OVERFLOW);
+  if (notrap & GFC_FPE_OVERFLOW)
+    fedisableexcept (FE_OVERFLOW);
+#endif
+
+#ifdef FE_UNDERFLOW
+  if (trap & GFC_FPE_UNDERFLOW)
+    feenableexcept (FE_UNDERFLOW);
+  if (notrap & GFC_FPE_UNDERFLOW)
+    fedisableexcept (FE_UNDERFLOW);
+#endif
+
+#ifdef FE_INEXACT
+  if (trap & GFC_FPE_INEXACT)
+    feenableexcept (FE_INEXACT);
+  if (notrap & GFC_FPE_INEXACT)
+    fedisableexcept (FE_INEXACT);
+#endif
+}
+
+
+int
+get_fpu_trap_exceptions (void)
+{
+  int exceptions = fegetexcept ();
+  int res = 0;
+
+#ifdef FE_INVALID
+  if (exceptions & FE_INVALID) res |= GFC_FPE_INVALID;
+#endif
+
+#ifdef FE_DENORMAL
+  if (exceptions & FE_DENORMAL) res |= GFC_FPE_DENORMAL;
+#endif
+
+#ifdef FE_DIVBYZERO
+  if (exceptions & FE_DIVBYZERO) res |= GFC_FPE_ZERO;
+#endif
+
+#ifdef FE_OVERFLOW
+  if (exceptions & FE_OVERFLOW) res |= GFC_FPE_OVERFLOW;
+#endif
+
+#ifdef FE_UNDERFLOW
+  if (exceptions & FE_UNDERFLOW) res |= GFC_FPE_UNDERFLOW;
+#endif
+
+#ifdef FE_INEXACT
+  if (exceptions & FE_INEXACT) res |= GFC_FPE_INEXACT;
+#endif
+
+  return res;
+}
+
+
+int
+support_fpu_trap (int flag)
+{
+  return support_fpu_flag (flag);
+}
+
+
+void set_fpu (void)
+{
+#ifndef FE_INVALID
+  if (options.fpe & GFC_FPE_INVALID)
     estr_write ("Fortran runtime warning: IEEE 'invalid operation' "
 	        "exception not supported.\n");
 #endif
 
 /* glibc does never have a FE_DENORMAL.  */
+#ifndef FE_DENORMAL
   if (options.fpe & GFC_FPE_DENORMAL)
-#ifdef FE_DENORMAL
-    feenableexcept (FE_DENORMAL);
-#else
     estr_write ("Fortran runtime warning: Floating point 'denormal operand' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef FE_DIVBYZERO
   if (options.fpe & GFC_FPE_ZERO)
-#ifdef FE_DIVBYZERO
-    feenableexcept (FE_DIVBYZERO);
-#else
     estr_write ("Fortran runtime warning: IEEE 'division by zero' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef FE_OVERFLOW
   if (options.fpe & GFC_FPE_OVERFLOW)
-#ifdef FE_OVERFLOW
-    feenableexcept (FE_OVERFLOW);
-#else
     estr_write ("Fortran runtime warning: IEEE 'overflow' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef FE_UNDERFLOW
   if (options.fpe & GFC_FPE_UNDERFLOW)
-#ifdef FE_UNDERFLOW
-    feenableexcept (FE_UNDERFLOW);
-#else
     estr_write ("Fortran runtime warning: IEEE 'underflow' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef FE_INEXACT
   if (options.fpe & GFC_FPE_INEXACT)
-#ifdef FE_INEXACT
-    feenableexcept (FE_INEXACT);
-#else
     estr_write ("Fortran runtime warning: IEEE 'inexact' "
 	        "exception not supported.\n");
 #endif
+
+  set_fpu_trap_exceptions (options.fpe, 0);
 }
 
 
@@ -129,6 +207,102 @@ get_fpu_except_flags (void)
 }
 
 
+void
+set_fpu_except_flags (int set, int clear)
+{
+  int exc_set = 0, exc_clr = 0;
+
+#ifdef FE_INVALID
+  if (set & GFC_FPE_INVALID)
+    exc_set |= FE_INVALID;
+  else if (clear & GFC_FPE_INVALID)
+    exc_clr |= FE_INVALID;
+#endif
+
+#ifdef FE_DIVBYZERO
+  if (set & GFC_FPE_ZERO)
+    exc_set |= FE_DIVBYZERO;
+  else if (clear & GFC_FPE_ZERO)
+    exc_clr |= FE_DIVBYZERO;
+#endif
+
+#ifdef FE_OVERFLOW
+  if (set & GFC_FPE_OVERFLOW)
+    exc_set |= FE_OVERFLOW;
+  else if (clear & GFC_FPE_OVERFLOW)
+    exc_clr |= FE_OVERFLOW;
+#endif
+
+#ifdef FE_UNDERFLOW
+  if (set & GFC_FPE_UNDERFLOW)
+    exc_set |= FE_UNDERFLOW;
+  else if (clear & GFC_FPE_UNDERFLOW)
+    exc_clr |= FE_UNDERFLOW;
+#endif
+
+#ifdef FE_DENORMAL
+  if (set & GFC_FPE_DENORMAL)
+    exc_set |= FE_DENORMAL;
+  else if (clear & GFC_FPE_DENORMAL)
+    exc_clr |= FE_DENORMAL;
+#endif
+
+#ifdef FE_INEXACT
+  if (set & GFC_FPE_INEXACT)
+    exc_set |= FE_INEXACT;
+  else if (clear & GFC_FPE_INEXACT)
+    exc_clr |= FE_INEXACT;
+#endif
+
+  feclearexcept (exc_clr);
+  feraiseexcept (exc_set);
+}
+
+
+int
+support_fpu_flag (int flag)
+{
+  if (flag & GFC_FPE_INVALID)
+  {
+#ifndef FE_INVALID
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_ZERO)
+  {
+#ifndef FE_DIVBYZERO
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_OVERFLOW)
+  {
+#ifndef FE_OVERFLOW
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_UNDERFLOW)
+  {
+#ifndef FE_UNDERFLOW
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_DENORMAL)
+  {
+#ifndef FE_DENORMAL
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_INEXACT)
+  {
+#ifndef FE_INEXACT
+    return 0;
+#endif
+  }
+
+  return 1;
+}
+
+
 int
 get_fpu_rounding_mode (void)
 {
@@ -199,3 +373,60 @@ set_fpu_rounding_mode (int mode)
 
   fesetround (rnd_mode);
 }
+
+
+int
+support_fpu_rounding_mode (int mode)
+{
+  switch (mode)
+    {
+      case GFC_FPE_TONEAREST:
+#ifdef FE_TONEAREST
+	return 1;
+#else
+	return 0;
+#endif
+
+#ifdef FE_UPWARD
+	return 1;
+#else
+	return 0;
+#endif
+
+#ifdef FE_DOWNWARD
+	return 1;
+#else
+	return 0;
+#endif
+
+#ifdef FE_TOWARDZERO
+	return 1;
+#else
+	return 0;
+#endif
+
+      default:
+	return 0;
+    }
+}
+
+
+void
+get_fpu_state (void *state)
+{
+  /* Check we can actually store the FPU state in the allocated size.  */
+  assert (sizeof(fenv_t) <= GFC_FPE_STATE_BUFFER_SIZE);
+
+  fegetenv (state);
+}
+
+
+void
+set_fpu_state (void *state)
+{
+  /* Check we can actually store the FPU state in the allocated size.  */
+  assert (sizeof(fenv_t) <= GFC_FPE_STATE_BUFFER_SIZE);
+
+  fesetenv (state);
+}
+
diff --git a/libgfortran/config/fpu-sysv.h b/libgfortran/config/fpu-sysv.h
index e7ba88f..0105cf7 100644
--- a/libgfortran/config/fpu-sysv.h
+++ b/libgfortran/config/fpu-sysv.h
@@ -25,73 +25,174 @@ see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
 
 /* FPU-related code for SysV platforms with fpsetmask().  */
 
+/* BSD and Solaris systems have slightly different types and functions
+   naming.  We deal with these here, to simplify the code below.  */
+
+#if HAVE_FP_EXCEPT
+# define FP_EXCEPT_TYPE fp_except
+#elif HAVE_FP_EXCEPT_T
+# define FP_EXCEPT_TYPE fp_except_t
+#else
+  choke me
+#endif
+
+#if HAVE_FP_RND
+# define FP_RND_TYPE fp_rnd
+#elif HAVE_FP_RND_T
+# define FP_RND_TYPE fp_rnd_t
+#else
+  choke me
+#endif
+
+#if HAVE_FPSETSTICKY
+# define FPSETSTICKY fpsetsticky
+#elif HAVE_FPRESETSTICKY
+# define FPSETSTICKY fpresetsticky
+#else
+  choke me
+#endif
+
+
 void
-set_fpu (void)
+set_fpu_trap_exceptions (int trap, int notrap)
 {
-  int cw = 0;
+  FP_EXCEPT_TYPE cw = fpgetmask();
 
-  if (options.fpe & GFC_FPE_INVALID)
 #ifdef FP_X_INV
+  if (trap & GFC_FPE_INVALID)
     cw |= FP_X_INV;
-#else
+  if (notrap & GFC_FPE_INVALID)
+    cw &= ~FP_X_INV;
+#endif
+
+#ifdef FP_X_DNML
+  if (trap & GFC_FPE_DENORMAL)
+    cw |= FP_X_DNML;
+  if (notrap & GFC_FPE_DENORMAL)
+    cw &= ~FP_X_DNML;
+#endif
+
+#ifdef FP_X_DZ
+  if (trap & GFC_FPE_ZERO)
+    cw |= FP_X_DZ;
+  if (notrap & GFC_FPE_ZERO)
+    cw &= ~FP_X_DZ;
+#endif
+
+#ifdef FP_X_OFL
+  if (trap & GFC_FPE_OVERFLOW)
+    cw |= FP_X_OFL;
+  if (notrap & GFC_FPE_OVERFLOW)
+    cw &= ~FP_X_OFL;
+#endif
+
+#ifdef FP_X_UFL
+  if (trap & GFC_FPE_UNDERFLOW)
+    cw |= FP_X_UFL;
+  if (notrap & GFC_FPE_UNDERFLOW)
+    cw &= ~FP_X_UFL;
+#endif
+
+#ifdef FP_X_IMP
+  if (trap & GFC_FPE_INEXACT)
+    cw |= FP_X_IMP;
+  if (notrap & GFC_FPE_INEXACT)
+    cw &= ~FP_X_IMP;
+#endif
+
+  fpsetmask(cw);
+}
+
+
+int
+get_fpu_trap_exceptions (void)
+{
+  int res = 0;
+  FP_EXCEPT_TYPE cw = fpgetmask();
+
+#ifdef FP_X_INV
+  if (cw & FP_X_INV) res |= GFC_FPE_INVALID;
+#endif
+
+#ifdef FP_X_DNML
+  if (cw & FP_X_DNML) res |= GFC_FPE_DENORMAL;
+#endif
+
+#ifdef FP_X_DZ
+  if (cw & FP_X_DZ) res |= GFC_FPE_ZERO;
+#endif
+
+#ifdef FP_X_OFL
+  if (cw & FP_X_OFL) res |= GFC_FPE_OVERFLOW;
+#endif
+
+#ifdef FP_X_UFL
+  if (cw & FP_X_UFL) res |= GFC_FPE_UNDERFLOW;
+#endif
+
+#ifdef FP_X_IMP
+  if (cw & FP_X_IMP) res |= GFC_FPE_INEXACT;
+#endif
+
+  return res;
+}
+
+
+int
+support_fpu_trap (int flag)
+{
+  return support_fpu_flag (flag);
+}
+
+
+void
+set_fpu (void)
+{
+#ifndef FP_X_INV
+  if (options.fpe & GFC_FPE_INVALID)
     estr_write ("Fortran runtime warning: IEEE 'invalid operation' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef FP_X_DNML
   if (options.fpe & GFC_FPE_DENORMAL)
-#ifdef FP_X_DNML
-    cw |= FP_X_DNML;
-#else
     estr_write ("Fortran runtime warning: Floating point 'denormal operand' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef FP_X_DZ
   if (options.fpe & GFC_FPE_ZERO)
-#ifdef FP_X_DZ
-    cw |= FP_X_DZ;
-#else
     estr_write ("Fortran runtime warning: IEEE 'division by zero' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef FP_X_OFL
   if (options.fpe & GFC_FPE_OVERFLOW)
-#ifdef FP_X_OFL
-    cw |= FP_X_OFL;
-#else
     estr_write ("Fortran runtime warning: IEEE 'overflow' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef FP_X_UFL
   if (options.fpe & GFC_FPE_UNDERFLOW)
-#ifdef FP_X_UFL
-    cw |= FP_X_UFL;
-#else
     estr_write ("Fortran runtime warning: IEEE 'underflow' "
 	        "exception not supported.\n");
 #endif
 
+#ifndef FP_X_IMP
   if (options.fpe & GFC_FPE_INEXACT)
-#ifdef FP_X_IMP
-    cw |= FP_X_IMP;
-#else
     estr_write ("Fortran runtime warning: IEEE 'inexact' "
 	        "exception not supported.\n");
 #endif
 
-  fpsetmask(cw);
+  set_fpu_trap_exceptions (options.fpe, 0);
 }
 
+
 int
 get_fpu_except_flags (void)
 {
   int result;
-#if HAVE_FP_EXCEPT
-  fp_except set_excepts;
-#elif HAVE_FP_EXCEPT_T
-  fp_except_t set_excepts;
-#else
-  choke me
-#endif
+  FP_EXCEPT_TYPE set_excepts;
 
   result = 0;
   set_excepts = fpgetsticky ();
@@ -130,6 +231,103 @@ get_fpu_except_flags (void)
 }
 
 
+void
+set_fpu_except_flags (int set, int clear)
+{
+  FP_EXCEPT_TYPE flags;
+
+  flags = fpgetsticky ();
+
+#ifdef FP_X_INV
+  if (set & GFC_FPE_INVALID)
+    flags |= FP_X_INV;
+  if (clear & GFC_FPE_INVALID)
+    flags &= ~FP_X_INV;
+#endif
+
+#ifdef FP_X_DZ
+  if (set & GFC_FPE_ZERO)
+    flags |= FP_X_DZ;
+  if (clear & GFC_FPE_ZERO)
+    flags &= ~FP_X_DZ;
+#endif
+
+#ifdef FP_X_OFL
+  if (set & GFC_FPE_OVERFLOW)
+    flags |= FP_X_OFL;
+  if (clear & GFC_FPE_OVERFLOW)
+    flags &= ~FP_X_OFL;
+#endif
+
+#ifdef FP_X_UFL
+  if (set & GFC_FPE_UNDERFLOW)
+    flags |= FP_X_UFL;
+  if (clear & GFC_FPE_UNDERFLOW)
+    flags &= ~FP_X_UFL;
+#endif
+
+#ifdef FP_X_DNML
+  if (set & GFC_FPE_DENORMAL)
+    flags |= FP_X_DNML;
+  if (clear & GFC_FPE_DENORMAL)
+    flags &= ~FP_X_DNML;
+#endif
+
+#ifdef FP_X_IMP
+  if (set & GFC_FPE_INEXACT)
+    flags |= FP_X_IMP;
+  if (clear & GFC_FPE_INEXACT)
+    flags &= ~FP_X_IMP;
+#endif
+
+  FPSETSTICKY (flags);
+}
+
+
+int
+support_fpu_flag (int flag)
+{
+  if (flag & GFC_FPE_INVALID)
+  {
+#ifndef FP_X_INV
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_ZERO)
+  {
+#ifndef FP_X_DZ
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_OVERFLOW)
+  {
+#ifndef FP_X_OFL
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_UNDERFLOW)
+  {
+#ifndef FP_X_UFL
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_DENORMAL)
+  {
+#ifndef FP_X_DNML
+    return 0;
+#endif
+  }
+  else if (flag & GFC_FPE_INEXACT)
+  {
+#ifndef FP_X_IMP
+    return 0;
+#endif
+  }
+
+  return 1;
+}
+
+
 int
 get_fpu_rounding_mode (void)
 {
@@ -163,13 +361,7 @@ get_fpu_rounding_mode (void)
 void
 set_fpu_rounding_mode (int mode)
 {
-#if HAVE_FP_RND
-  fp_rnd rnd_mode;
-#elif HAVE_FP_RND_T
-  fp_rnd_t rnd_mode;
-#else
-  choke me
-#endif
+  FP_RND_TYPE rnd_mode;
 
   switch (mode)
     {
@@ -201,3 +393,78 @@ set_fpu_rounding_mode (int mode)
     }
   fpsetround (rnd_mode);
 }
+
+
+int
+support_fpu_rounding_mode (int mode)
+{
+  switch (mode)
+    {
+      case GFC_FPE_TONEAREST:
+#ifdef FP_RN
+	return 1;
+#else
+	return 0;
+#endif
+
+      case GFC_FPE_UPWARD:
+#ifdef FP_RP
+	return 1;
+#else
+	return 0;
+#endif
+
+      case GFC_FPE_DOWNWARD:
+#ifdef FP_RM
+	return 1;
+#else
+	return 0;
+#endif
+
+      case GFC_FPE_TOWARDZERO:
+#ifdef FP_RZ
+	return 1;
+#else
+	return 0;
+#endif
+
+      default:
+	return 0;
+    }
+}
+
+
+typedef struct
+{
+  FP_EXCEPT_TYPE mask;
+  FP_EXCEPT_TYPE sticky;
+  FP_RND_TYPE round;
+} fpu_state_t;
+
+
+void
+get_fpu_state (void *s)
+{
+  fpu_state_t *state = s;
+
+  /* Check we can actually store the FPU state in the allocated size.  */
+  assert (sizeof(fpu_state_t) <= GFC_FPE_STATE_BUFFER_SIZE);
+
+  state->mask = fpgetmask ();
+  state->sticky = fpgetsticky ();
+  state->round = fpgetround ();
+}
+
+void
+set_fpu_state (void *s)
+{
+  fpu_state_t *state = s;
+
+  /* Check we can actually store the FPU state in the allocated size.  */
+  assert (sizeof(fpu_state_t) <= GFC_FPE_STATE_BUFFER_SIZE);
+
+  fpsetmask (state->mask);
+  FPSETSTICKY (state->sticky);
+  fpsetround (state->round);
+}
+