Fix PR fortran/93871 and re-implement degree-valued trigonometric intrinsics.

2020-04-01  Fritz Reese  <foreese@gcc.gnu.org>
	    Steven G. Kargl  <kargl@gcc.gnu.org>

gcc/fortran/ChangeLog

	PR fortran/93871
	* gfortran.h (GFC_ISYM_ACOSD, GFC_ISYM_ASIND, GFC_ISYM_ATAN2D,
	GFC_ISYM_ATAND, GFC_ISYM_COSD, GFC_ISYM_COTAND, GFC_ISYM_SIND,
	GFC_ISYM_TAND): New.
	* intrinsic.c (add_functions): Remove check for flag_dec_math.
	Give degree trig functions simplification and name resolution
	functions (e.g, gfc_simplify_atrigd () and gfc_resolve_atrigd ()).
	(do_simplify): Remove special casing of degree trig functions.
	* intrinsic.h (gfc_simplify_acosd, gfc_simplify_asind,
	gfc_simplify_atand, gfc_simplify_cosd, gfc_simplify_cotand,
	gfc_simplify_sind, gfc_simplify_tand, gfc_resolve_trigd2): Add new
	prototypes.
	(gfc_simplify_atrigd, gfc_simplify_trigd, gfc_resolve_cotan,
	resolve_atrigd): Remove prototypes of deleted functions.
	* iresolve.c (is_trig_resolved, copy_replace_function_shallow,
	gfc_resolve_cotan, get_radians, get_degrees, resolve_trig_call,
	gfc_resolve_atrigd, gfc_resolve_atan2d): Delete functions.
	(gfc_resolve_trigd, gfc_resolve_trigd2): Resolve to library functions.
	* simplify.c (rad2deg, deg2rad, gfc_simplify_acosd, gfc_simplify_asind,
	gfc_simplify_atand, gfc_simplify_atan2d, gfc_simplify_cosd,
	gfc_simplify_sind, gfc_simplify_tand, gfc_simplify_cotand): New
	functions.
	(gfc_simplify_atan2): Fix error message.
	(simplify_trig_call, gfc_simplify_trigd, gfc_simplify_atrigd,
	radians_f): Delete functions.
	* trans-intrinsic.c: Add LIB_FUNCTION decls for sind, cosd, tand.
	(rad2deg, gfc_conv_intrinsic_atrigd, gfc_conv_intrinsic_cotan,
	gfc_conv_intrinsic_cotand, gfc_conv_intrinsic_atan2d): New functions.
	(gfc_conv_intrinsic_function): Handle ACOSD, ASIND, ATAND, COTAN,
	COTAND, ATAN2D.
	* trigd_fe.inc: New file. Included by simplify.c to implement
	simplify_sind, simplify_cosd, simplify_tand with code common to the
	libgfortran implementation.

gcc/testsuite/ChangeLog

	PR fortran/93871
	* gfortran.dg/dec_math.f90: Extend coverage to real(10) and real(16).
	* gfortran.dg/dec_math_2.f90: New test.
	* gfortran.dg/dec_math_3.f90: Likewise.
	* gfortran.dg/dec_math_4.f90: Likewise.
	* gfortran.dg/dec_math_5.f90: Likewise.

libgfortran/ChangeLog

	PR fortran/93871
	* Makefile.am, Makefile.in: New make rule for intrinsics/trigd.c.
	* gfortran.map: New routines for {sind, cosd, tand}X{r4, r8, r10, r16}.
	* intrinsics/trigd.c, intrinsics/trigd_lib.inc, intrinsics/trigd.inc:
	New files. Defines native degree-valued trig functions.

1 files changed, 184 insertions, 116 deletions

diff --git a/gcc/fortran/simplify.c b/gcc/fortran/simplify.c
index 66ed925..f63f63c 100644
--- a/gcc/fortran/simplify.c
+++ b/gcc/fortran/simplify.c
@@ -1107,6 +1107,91 @@ gfc_simplify_asin (gfc_expr *x)
 }
 
 
+/* Convert radians to degrees, i.e., x * 180 / pi.  */
+
+static void
+rad2deg (mpfr_t x)
+{
+  mpfr_t tmp;
+
+  mpfr_init (tmp);
+  mpfr_const_pi (tmp, GFC_RND_MODE);
+  mpfr_mul_ui (x, x, 180, GFC_RND_MODE);
+  mpfr_div (x, x, tmp, GFC_RND_MODE);
+  mpfr_clear (tmp);
+}
+
+
+/* Simplify ACOSD(X) where the returned value has units of degree.  */
+
+gfc_expr *
+gfc_simplify_acosd (gfc_expr *x)
+{
+  gfc_expr *result;
+
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  if (mpfr_cmp_si (x->value.real, 1) > 0
+      || mpfr_cmp_si (x->value.real, -1) < 0)
+    {
+      gfc_error ("Argument of ACOSD at %L must be between -1 and 1",
+		 &x->where);
+      return &gfc_bad_expr;
+    }
+
+  result = gfc_get_constant_expr (x->ts.type, x->ts.kind, &x->where);
+  mpfr_acos (result->value.real, x->value.real, GFC_RND_MODE);
+  rad2deg (result->value.real);
+
+  return range_check (result, "ACOSD");
+}
+
+
+/* Simplify asind (x) where the returned value has units of degree. */
+
+gfc_expr *
+gfc_simplify_asind (gfc_expr *x)
+{
+  gfc_expr *result;
+
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  if (mpfr_cmp_si (x->value.real, 1) > 0
+      || mpfr_cmp_si (x->value.real, -1) < 0)
+    {
+      gfc_error ("Argument of ASIND at %L must be between -1 and 1",
+		 &x->where);
+      return &gfc_bad_expr;
+    }
+
+  result = gfc_get_constant_expr (x->ts.type, x->ts.kind, &x->where);
+  mpfr_asin (result->value.real, x->value.real, GFC_RND_MODE);
+  rad2deg (result->value.real);
+
+  return range_check (result, "ASIND");
+}
+
+
+/* Simplify atand (x) where the returned value has units of degree. */
+
+gfc_expr *
+gfc_simplify_atand (gfc_expr *x)
+{
+  gfc_expr *result;
+
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  result = gfc_get_constant_expr (x->ts.type, x->ts.kind, &x->where);
+  mpfr_atan (result->value.real, x->value.real, GFC_RND_MODE);
+  rad2deg (result->value.real);
+
+  return range_check (result, "ATAND");
+}
+
+
 gfc_expr *
 gfc_simplify_asinh (gfc_expr *x)
 {
@@ -1208,8 +1293,8 @@ gfc_simplify_atan2 (gfc_expr *y, gfc_expr *x)
 
   if (mpfr_zero_p (y->value.real) && mpfr_zero_p (x->value.real))
     {
-      gfc_error ("If first argument of ATAN2 %L is zero, then the "
-		 "second argument must not be zero", &x->where);
+      gfc_error ("If first argument of ATAN2 at %L is zero, then the "
+		 "second argument must not be zero", &y->where);
       return &gfc_bad_expr;
     }
 
@@ -1736,172 +1821,153 @@ gfc_simplify_conjg (gfc_expr *e)
   return range_check (result, "CONJG");
 }
 
-/* Return the simplification of the constant expression in icall, or NULL
-   if the expression is not constant.  */
 
-static gfc_expr *
-simplify_trig_call (gfc_expr *icall)
-{
-  gfc_isym_id func = icall->value.function.isym->id;
-  gfc_expr *x = icall->value.function.actual->expr;
-
-  /* The actual simplifiers will return NULL for non-constant x.  */
-  switch (func)
-    {
-    case GFC_ISYM_ACOS:
-      return gfc_simplify_acos (x);
-    case GFC_ISYM_ASIN:
-      return gfc_simplify_asin (x);
-    case GFC_ISYM_ATAN:
-      return gfc_simplify_atan (x);
-    case GFC_ISYM_COS:
-      return gfc_simplify_cos (x);
-    case GFC_ISYM_COTAN:
-      return gfc_simplify_cotan (x);
-    case GFC_ISYM_SIN:
-      return gfc_simplify_sin (x);
-    case GFC_ISYM_TAN:
-      return gfc_simplify_tan (x);
-    default:
-      gfc_internal_error ("in simplify_trig_call(): Bad intrinsic");
+/* Simplify atan2d (x) where the unit is degree.  */
+
+gfc_expr *
+gfc_simplify_atan2d (gfc_expr *y, gfc_expr *x)
+{
+  gfc_expr *result;
+
+  if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  if (mpfr_zero_p (y->value.real) && mpfr_zero_p (x->value.real))
+    {
+      gfc_error ("If first argument of ATAN2D at %L is zero, then the "
+		 "second argument must not be zero", &y->where);
+      return &gfc_bad_expr;
     }
+
+  result = gfc_get_constant_expr (x->ts.type, x->ts.kind, &x->where);
+  mpfr_atan2 (result->value.real, y->value.real, x->value.real, GFC_RND_MODE);
+  rad2deg (result->value.real);
+
+  return range_check (result, "ATAN2D");
 }
 
-/* Convert a floating-point number from radians to degrees.  */
 
-static void
-degrees_f (mpfr_t x, mpfr_rnd_t rnd_mode)
+gfc_expr *
+gfc_simplify_cos (gfc_expr *x)
 {
-  mpfr_t tmp;
-  mpfr_init (tmp);
+  gfc_expr *result;
 
-  /* Set x = x * 180.  */
-  mpfr_mul_ui (x, x, 180, rnd_mode);
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
 
-  /* Set x = x / pi.  */
-  mpfr_const_pi (tmp, rnd_mode);
-  mpfr_div (x, x, tmp, rnd_mode);
+  result = gfc_get_constant_expr (x->ts.type, x->ts.kind, &x->where);
 
-  mpfr_clear (tmp);
+  switch (x->ts.type)
+    {
+      case BT_REAL:
+	mpfr_cos (result->value.real, x->value.real, GFC_RND_MODE);
+	break;
+
+      case BT_COMPLEX:
+	gfc_set_model_kind (x->ts.kind);
+	mpc_cos (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+	break;
+
+      default:
+	gfc_internal_error ("in gfc_simplify_cos(): Bad type");
+    }
+
+  return range_check (result, "COS");
 }
 
-/* Convert a floating-point number from degrees to radians.  */
 
 static void
-radians_f (mpfr_t x, mpfr_rnd_t rnd_mode)
+deg2rad (mpfr_t x)
 {
-  mpfr_t tmp;
-  mpfr_init (tmp);
+  mpfr_t d2r;
 
-  /* Set x = x % 360 to avoid offsets with large angles.  */
-  mpfr_set_ui (tmp, 360, rnd_mode);
-  mpfr_fmod (tmp, x, tmp, rnd_mode);
+  mpfr_init (d2r);
+  mpfr_const_pi (d2r, GFC_RND_MODE);
+  mpfr_div_ui (d2r, d2r, 180, GFC_RND_MODE);
+  mpfr_mul (x, x, d2r, GFC_RND_MODE);
+  mpfr_clear (d2r);
+}
 
-  /* Set x = x * pi.  */
-  mpfr_const_pi (tmp, rnd_mode);
-  mpfr_mul (x, x, tmp, rnd_mode);
 
-  /* Set x = x / 180.  */
-  mpfr_div_ui (x, x, 180, rnd_mode);
-
-  mpfr_clear (tmp);
-}
+/* Simplification routines for SIND, COSD, TAND.  */
+#include "trigd_fe.inc"
 
 
-/* Convert argument to radians before calling a trig function.  */
+/* Simplify COSD(X) where X has the unit of degree.  */
 
 gfc_expr *
-gfc_simplify_trigd (gfc_expr *icall)
+gfc_simplify_cosd (gfc_expr *x)
 {
-  gfc_expr *arg;
-
-  arg = icall->value.function.actual->expr;
+  gfc_expr *result;
 
-  if (arg->ts.type != BT_REAL)
-    gfc_internal_error ("in gfc_simplify_trigd(): Bad type");
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
 
-  if (arg->expr_type == EXPR_CONSTANT)
-    /* Convert constant to radians before passing off to simplifier.  */
-    radians_f (arg->value.real, GFC_RND_MODE);
+  result = gfc_get_constant_expr (x->ts.type, x->ts.kind, &x->where);
+  mpfr_set (result->value.real, x->value.real, GFC_RND_MODE);
+  simplify_cosd (result->value.real);
 
-  /* Let the usual simplifier take over - we just simplified the arg.  */
-  return simplify_trig_call (icall);
+  return range_check (result, "COSD");
 }
 
-/* Convert result of an inverse trig function to degrees.  */
+
+/* Simplify SIND(X) where X has the unit of degree.  */
 
 gfc_expr *
-gfc_simplify_atrigd (gfc_expr *icall)
+gfc_simplify_sind (gfc_expr *x)
 {
   gfc_expr *result;
 
-  if (icall->value.function.actual->expr->ts.type != BT_REAL)
-    gfc_internal_error ("in gfc_simplify_atrigd(): Bad type");
-
-  /* See if another simplifier has work to do first.  */
-  result = simplify_trig_call (icall);
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
 
-  if (result && result->expr_type == EXPR_CONSTANT)
-    {
-      /* Convert constant to degrees after passing off to actual simplifier.  */
-      degrees_f (result->value.real, GFC_RND_MODE);
-      return result;
-    }
+  result = gfc_get_constant_expr (x->ts.type, x->ts.kind, &x->where);
+  mpfr_set (result->value.real, x->value.real, GFC_RND_MODE);
+  simplify_sind (result->value.real);
 
-  /* Let gfc_resolve_atrigd take care of the non-constant case.  */
-  return NULL;
+  return range_check (result, "SIND");
 }
 
-/* Convert the result of atan2 to degrees.  */
+
+/* Simplify TAND(X) where X has the unit of degree.  */
 
 gfc_expr *
-gfc_simplify_atan2d (gfc_expr *y, gfc_expr *x)
+gfc_simplify_tand (gfc_expr *x)
 {
   gfc_expr *result;
 
-  if (x->ts.type != BT_REAL || y->ts.type != BT_REAL)
-    gfc_internal_error ("in gfc_simplify_atan2d(): Bad type");
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
 
-  if (x->expr_type == EXPR_CONSTANT && y->expr_type == EXPR_CONSTANT)
-    {
-      result = gfc_simplify_atan2 (y, x);
-      if (result != NULL)
-	{
-	  degrees_f (result->value.real, GFC_RND_MODE);
-	  return result;
-	}
-    }
+  result = gfc_get_constant_expr (x->ts.type, x->ts.kind, &x->where);
+  mpfr_set (result->value.real, x->value.real, GFC_RND_MODE);
+  simplify_tand (result->value.real);
 
-  /* Let gfc_resolve_atan2d take care of the non-constant case.  */
-  return NULL;
+  return range_check (result, "TAND");
 }
 
+
+/* Simplify COTAND(X) where X has the unit of degree.  */
+
 gfc_expr *
-gfc_simplify_cos (gfc_expr *x)
+gfc_simplify_cotand (gfc_expr *x)
 {
   gfc_expr *result;
 
   if (x->expr_type != EXPR_CONSTANT)
     return NULL;
 
+  /* Implement COTAND = -TAND(x+90).
+     TAND offers correct exact values for multiples of 30 degrees.
+     This implementation is also compatible with the behavior of some legacy
+     compilers.  Keep this consistent with gfc_conv_intrinsic_cotand.  */
   result = gfc_get_constant_expr (x->ts.type, x->ts.kind, &x->where);
+  mpfr_set (result->value.real, x->value.real, GFC_RND_MODE);
+  mpfr_add_ui (result->value.real, result->value.real, 90, GFC_RND_MODE);
+  simplify_tand (result->value.real);
+  mpfr_neg (result->value.real, result->value.real, GFC_RND_MODE);
 
-  switch (x->ts.type)
-    {
-      case BT_REAL:
-	mpfr_cos (result->value.real, x->value.real, GFC_RND_MODE);
-	break;
-
-      case BT_COMPLEX:
-	gfc_set_model_kind (x->ts.kind);
-	mpc_cos (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
-	break;
-
-      default:
-	gfc_internal_error ("in gfc_simplify_cos(): Bad type");
-    }
-
-  return range_check (result, "COS");
+  return range_check (result, "COTAND");
 }
 
 
@@ -7778,6 +7844,8 @@ gfc_simplify_sum (gfc_expr *array, gfc_expr *dim, gfc_expr *mask)
 }
 
 
+/* Simplify COTAN(X) where X has the unit of radian.  */
+
 gfc_expr *
 gfc_simplify_cotan (gfc_expr *x)
 {
@@ -7799,8 +7867,8 @@ gfc_simplify_cotan (gfc_expr *x)
       /* There is no builtin mpc_cot, so compute cot = cos / sin.  */
       val = &result->value.complex;
       mpc_init2 (swp, mpfr_get_default_prec ());
-      mpc_cos (swp, x->value.complex, GFC_MPC_RND_MODE);
-      mpc_sin (*val, x->value.complex, GFC_MPC_RND_MODE);
+      mpc_sin_cos (*val, swp, x->value.complex, GFC_MPC_RND_MODE,
+		   GFC_MPC_RND_MODE);
       mpc_div (*val, swp, *val, GFC_MPC_RND_MODE);
       mpc_clear (swp);
       break;