re PR middle-end/30251 (Evaluate bessel functions at compile-time)

	PR middle-end/30251
	* builtins.c (do_mpfr_bessel_n): New.
	(fold_builtin_1): Handle BUILT_IN_J0 and BUILT_IN_J1.
	(fold_builtin_2): Handle BUILT_IN_JN.

testsuite:
	* gcc.dg/torture/builtin-math-4.c: New test.

From-SVN: r124818

4 files changed, 201 insertions, 0 deletions

diff --git a/gcc/ChangeLog b/gcc/ChangeLog
index 371c5e3..6e3c321 100644
--- a/gcc/ChangeLog
+++ b/gcc/ChangeLog
@@ -1,3 +1,10 @@
+2007-05-17  Kaveh R. Ghazi  <ghazi@caip.rutgers.edu>
+
+	PR middle-end/30251
+	* builtins.c (do_mpfr_bessel_n): New.
+	(fold_builtin_1): Handle BUILT_IN_J0 and BUILT_IN_J1.
+	(fold_builtin_2): Handle BUILT_IN_JN.
+
 2007-05-17  Danny Smith  <dannysmith@users.sourceforge.net>
 
 	PR target/31965
diff --git a/gcc/builtins.c b/gcc/builtins.c
index 2de66dd..221b5a3 100644
--- a/gcc/builtins.c
+++ b/gcc/builtins.c
@@ -231,6 +231,11 @@ static tree do_mpfr_arg2 (tree, tree, tree,
 static tree do_mpfr_arg3 (tree, tree, tree, tree,
 			  int (*)(mpfr_ptr, mpfr_srcptr, mpfr_srcptr, mpfr_srcptr, mp_rnd_t));
 static tree do_mpfr_sincos (tree, tree, tree);
+#if MPFR_VERSION >= MPFR_VERSION_NUM(2,3,0)
+static tree do_mpfr_bessel_n (tree, tree, tree,
+			      int (*)(mpfr_ptr, long, mpfr_srcptr, mp_rnd_t),
+			      const REAL_VALUE_TYPE *, bool);
+#endif
 
 /* Return true if NODE should be considered for inline expansion regardless
    of the optimization level.  This means whenever a function is invoked with
@@ -9766,6 +9771,20 @@ fold_builtin_1 (tree fndecl, tree arg0, bool ignore)
 			     &dconstm1, NULL, false);
     break;
 
+#if MPFR_VERSION >= MPFR_VERSION_NUM(2,3,0)
+    CASE_FLT_FN (BUILT_IN_J0):
+      if (validate_arg (arg0, REAL_TYPE))
+	return do_mpfr_arg1 (arg0, type, mpfr_j0,
+			     NULL, NULL, 0);
+    break;
+
+    CASE_FLT_FN (BUILT_IN_J1):
+      if (validate_arg (arg0, REAL_TYPE))
+	return do_mpfr_arg1 (arg0, type, mpfr_j1,
+			     NULL, NULL, 0);
+    break;
+#endif
+
     CASE_FLT_FN (BUILT_IN_NAN):
     case BUILT_IN_NAND32:
     case BUILT_IN_NAND64:
@@ -9876,6 +9895,13 @@ fold_builtin_2 (tree fndecl, tree arg0, tree arg1, bool ignore)
 
   switch (fcode)
     {
+#if MPFR_VERSION >= MPFR_VERSION_NUM(2,3,0)
+    CASE_FLT_FN (BUILT_IN_JN):
+      if (validate_arg (arg0, INTEGER_TYPE)
+	  && validate_arg (arg1, REAL_TYPE))
+	return do_mpfr_bessel_n (arg0, arg1, type, mpfr_jn, NULL, 0);
+    break;
+#endif
 
     CASE_FLT_FN (BUILT_IN_ATAN2):
       if (validate_arg (arg0, REAL_TYPE)
@@ -12505,3 +12531,50 @@ do_mpfr_sincos (tree arg, tree arg_sinp, tree arg_cosp)
     }
   return result;
 }
+
+#if MPFR_VERSION >= MPFR_VERSION_NUM(2,3,0)
+/* If argument ARG1 is an INTEGER_CST and ARG2 is a REAL_CST, call the
+   two-argument mpfr order N Bessel function FUNC on them and return
+   the resulting value as a tree with type TYPE.  The mpfr precision
+   is set to the precision of TYPE.  We assume that function FUNC
+   returns zero if the result could be calculated exactly within the
+   requested precision.  */
+static tree
+do_mpfr_bessel_n (tree arg1, tree arg2, tree type,
+		  int (*func)(mpfr_ptr, long, mpfr_srcptr, mp_rnd_t),
+		  const REAL_VALUE_TYPE *min, bool inclusive)
+{
+  tree result = NULL_TREE;
+
+  STRIP_NOPS (arg1);
+  STRIP_NOPS (arg2);
+
+  /* To proceed, MPFR must exactly represent the target floating point
+     format, which only happens when the target base equals two.  */
+  if (REAL_MODE_FORMAT (TYPE_MODE (type))->b == 2
+      && host_integerp (arg1, 0)
+      && TREE_CODE (arg2) == REAL_CST && !TREE_OVERFLOW (arg2))
+    {
+      const HOST_WIDE_INT n = tree_low_cst(arg1, 0);
+      const REAL_VALUE_TYPE *const ra = &TREE_REAL_CST (arg2);
+
+      if (n == (long)n
+	  && !real_isnan (ra) && !real_isinf (ra)
+	  && (!min || real_compare (inclusive ? GE_EXPR: GT_EXPR , ra, min)))
+        {
+	  const int prec = REAL_MODE_FORMAT (TYPE_MODE (type))->p;
+	  int inexact;
+	  mpfr_t m;
+
+	  mpfr_init2 (m, prec);
+	  mpfr_from_real (m, ra, GMP_RNDN);
+	  mpfr_clear_flags ();
+	  inexact = func (m, n, m, GMP_RNDN);
+	  result = do_mpfr_ckconv (m, type, inexact);
+	  mpfr_clear (m);
+	}
+    }
+  
+  return result;
+}
+#endif
diff --git a/gcc/testsuite/ChangeLog b/gcc/testsuite/ChangeLog
index 77efbca..c2be740 100644
--- a/gcc/testsuite/ChangeLog
+++ b/gcc/testsuite/ChangeLog
@@ -1,3 +1,7 @@
+2007-05-17  Kaveh R. Ghazi  <ghazi@caip.rutgers.edu>
+
+	* gcc.dg/torture/builtin-math-4.c: New test.
+
 2007-05-17  Janis Johnson  <janis187@us.ibm.com>
 	    Manuel Lopez-Ibanez  <manu@gcc.gnu.org>
 
diff --git a/gcc/testsuite/gcc.dg/torture/builtin-math-4.c b/gcc/testsuite/gcc.dg/torture/builtin-math-4.c
new file mode 100644
index 0000000..294dd93
--- /dev/null
+++ b/gcc/testsuite/gcc.dg/torture/builtin-math-4.c
@@ -0,0 +1,117 @@
+/* Copyright (C) 2007  Free Software Foundation.
+
+   Verify that built-in math function constant folding of constant
+   arguments is correctly performed by the compiler.  This testcase is
+   for functionality that was available as of mpfr-2.3.0.
+
+   Origin: Kaveh R. Ghazi,  April 23, 2007.  */
+
+/* { dg-do link } */
+/* Expect failures at least until mpfr-2.3.0 is released. */
+/* { dg-xfail-if "mpfr-2.3.0" { *-*-* } } */
+
+/* All references to link_error should go away at compile-time.  */
+extern void link_error(int);
+
+/* Return TRUE if the sign of X != sign of Y.  This is important when
+   comparing signed zeros.  */
+#define CKSGN_F(X,Y) \
+  (__builtin_copysignf(1.0F,(X)) != __builtin_copysignf(1.0F,(Y)))
+#define CKSGN(X,Y) \
+  (__builtin_copysign(1.0,(X)) != __builtin_copysign(1.0,(Y)))
+#define CKSGN_L(X,Y) \
+  (__builtin_copysignl(1.0L,(X)) != __builtin_copysignl(1.0L,(Y)))
+
+/* Test that FUNC(ARG) == (RES).  */
+#define TESTIT(FUNC,ARG,RES) do { \
+  if (__builtin_##FUNC##f(ARG##F) != RES##F \
+      || CKSGN_F(__builtin_##FUNC##f(ARG##F),RES##F)) \
+    link_error(__LINE__); \
+  if (__builtin_##FUNC(ARG) != RES \
+      || CKSGN(__builtin_##FUNC(ARG),RES)) \
+    link_error(__LINE__); \
+  if (__builtin_##FUNC##l(ARG##L) != RES##L \
+      || CKSGN_L(__builtin_##FUNC##l(ARG##L),RES##L)) \
+    link_error(__LINE__); \
+  } while (0)
+
+/* Range test, check that (LOW) < FUNC(ARG) < (HI).  */
+#define TESTIT_R(FUNC,ARG,LOW,HI) do { \
+  if (__builtin_##FUNC##f(ARG) <= (LOW) || __builtin_##FUNC##f(ARG) >= (HI)) \
+    link_error(__LINE__); \
+  if (__builtin_##FUNC(ARG) <= (LOW) || __builtin_##FUNC(ARG) >= (HI)) \
+    link_error(__LINE__); \
+  if (__builtin_##FUNC##l(ARG) <= (LOW) || __builtin_##FUNC##l(ARG) >= (HI)) \
+    link_error(__LINE__); \
+  } while (0)
+
+/* Test that FUNC(ARG1, ARG2) == (RES).  */
+#define TESTIT2(FUNC,ARG1,ARG2,RES) do { \
+  if (__builtin_##FUNC##f(ARG1, ARG2##F) != RES##F \
+      || CKSGN_F(__builtin_##FUNC##f(ARG1,ARG2##F),RES##F)) \
+    link_error(__LINE__); \
+  if (__builtin_##FUNC(ARG1, ARG2) != RES \
+      || CKSGN(__builtin_##FUNC(ARG1,ARG2),RES)) \
+    link_error(__LINE__); \
+  if (__builtin_##FUNC##l(ARG1, ARG2##L) != RES##L \
+      || CKSGN_L(__builtin_##FUNC##l(ARG1,ARG2##L),RES##L)) \
+    link_error(__LINE__); \
+  } while (0)
+
+/* Range test, check that (LOW) < FUNC(ARG1,ARG2) < (HI).  */
+#define TESTIT2_R(FUNC,ARG1,ARG2,LOW,HI) do { \
+  if (__builtin_##FUNC##f(ARG1, ARG2##F) <= (LOW) \
+      || __builtin_##FUNC##f(ARG1, ARG2##F) >= (HI)) \
+    link_error(__LINE__); \
+  if (__builtin_##FUNC(ARG1, ARG2) <= (LOW) \
+      || __builtin_##FUNC(ARG1, ARG2) >= (HI)) \
+    link_error(__LINE__); \
+  if (__builtin_##FUNC##l(ARG1, ARG2##L) <= (LOW) \
+      || __builtin_##FUNC##l(ARG1, ARG2##L) >= (HI)) \
+    link_error(__LINE__); \
+  } while (0)
+
+int main (void)
+{
+  TESTIT (j0, 0.0, 1.0); /* j0(0) == 1 */
+  TESTIT (j0, -0.0, 1.0); /* j0(-0) == 1 */
+  TESTIT_R (j0, 1.0, 0.765, 0.766); /* j0(1) == 0.7651... */
+  TESTIT_R (j0, -1.0, 0.765, 0.766); /* j0(-1) == 0.7651... */
+
+  TESTIT (j1, 0.0, 0.0); /* j1(0) == 0 */
+  TESTIT (j1, -0.0, -0.0); /* j1(-0) == -0 */
+  TESTIT_R (j1, 1.0, 0.44, 0.45); /* j1(1) == 0.440... */
+  TESTIT_R (j1, -1.0, -0.45, -0.44); /* j1(-1) == -0.440... */
+
+  TESTIT2 (jn, 5, 0.0, 0.0); /* jn(5,0) == 0 */
+  TESTIT2 (jn, 5, -0.0, -0.0); /* jn(5,-0) == -0 */
+  TESTIT2 (jn, 6, 0.0, 0.0); /* jn(6,0) == 0 */
+  TESTIT2 (jn, 6, -0.0, 0.0); /* jn(6,-0) == 0 */
+
+  TESTIT2 (jn, -5, 0.0, -0.0); /* jn(-5,0) == -0 */
+  TESTIT2 (jn, -5, -0.0, 0.0); /* jn(-5,-0) == 0 */
+  TESTIT2 (jn, -6, 0.0, 0.0); /* jn(-6,0) == 0 */
+  TESTIT2 (jn, -6, -0.0, 0.0); /* jn(-6,-0) == 0 */
+
+  TESTIT2_R (jn, 2, 1.0, 0.11, 0.12); /* jn(2,1) == 0.114... */
+  TESTIT2_R (jn, 2, -1.0, 0.11, 0.12); /* jn(2,-1) == 0.114... */
+  TESTIT2_R (jn, 3, 5.0, 0.36, 0.37); /* jn(3,5) == 0.364... */
+  TESTIT2_R (jn, 3, -5.0, -0.37, -0.36); /* jn(3,-5) == -0.364... */
+
+  TESTIT2_R (jn, -2, 1.0, 0.11, 0.12); /* jn(-2,1) == 0.114... */
+  TESTIT2_R (jn, -2, -1.0, 0.11, 0.12); /* jn(-2,-1) == 0.114... */
+  TESTIT2_R (jn, -3, 5.0, -0.37, -0.36); /* jn(-3,5) == -0.364... */
+  TESTIT2_R (jn, -3, -5.0, 0.36, 0.37); /* jn(-3,-5) == 0.364... */
+
+  TESTIT2_R (jn, 4, 3.5, 0.20, 0.21); /* jn(4,3.5) == 0.204... */
+  TESTIT2_R (jn, 4, -3.5, 0.20, 0.21); /* jn(4,-3.5) == 0.204... */
+  TESTIT2_R (jn, 5, 4.6, 0.20, 0.21); /* jn(5,4.6) == 0.207... */
+  TESTIT2_R (jn, 5, -4.6, -0.21, -0.20); /* jn(5,-4.6) == -0.207... */
+
+  TESTIT2_R (jn, -4, 3.5, 0.20, 0.21); /* jn(-4,3.5) == 0.204... */
+  TESTIT2_R (jn, -4, -3.5, 0.20, 0.21); /* jn(-4,-3.5) == 0.204... */
+  TESTIT2_R (jn, -5, 4.6, -0.21, -0.20); /* jn(-5,4.6) == -0.207... */
+  TESTIT2_R (jn, -5, -4.6, 0.20, 0.21); /* jn(-5,-4.6) == 0.207... */
+
+  return 0;
+}