re PR fortran/25104 ([F2003] Non-initialization expr. as case-selector)

2009-06-07  Daniel Franke  <franke.daniel@gmail.com>

        PR fortran/25104
        PR fortran/29962
        * intrinsic.h (gfc_simplify_dot_product): New prototype.
        (gfc_simplify_matmul): Likewise.
        (gfc_simplify_transpose): Likewise.
        * intrinsic.c (add_functions): Added new simplifier callbacks.
        * simplify.c (init_result_expr): New.
        (compute_dot_product): New.
        (gfc_simplify_dot_product): New.
        (gfc_simplify_matmul): New.
        (gfc_simplify_transpose): New.
        * expr.c (check_transformational): Allow transformational
        * intrinsics
        with simplifier in initialization expression.


2009-06-07  Daniel Franke  <franke.daniel@gmail.com>

        PR fortran/25104
        PR fortran/29962
        * gfortran.dg/dot_product_1.f03: New.
        * gfortran.dg/matmul_8.f03: New.
        * gfortran.dg/transpose_3.f03: New.

From-SVN: r148243

gcc/fortran/ChangeLog

+2009-06-07  Daniel Franke  <franke.daniel@gmail.com>
+
+	PR fortran/25104
+	PR fortran/29962
+	* intrinsic.h (gfc_simplify_dot_product): New prototype.
+	(gfc_simplify_matmul): Likewise.
+	(gfc_simplify_transpose): Likewise.
+	* intrinsic.c (add_functions): Added new simplifier callbacks.
+	* simplify.c (init_result_expr): New.
+	(compute_dot_product): New.
+	(gfc_simplify_dot_product): New.
+	(gfc_simplify_matmul): New.
+	(gfc_simplify_transpose): New.
+	* expr.c (check_transformational): Allow transformational intrinsics
+	with simplifier in initialization expression.
+
 2009-06-06  Daniel Franke  <franke.daniel@gmail.com>
 
 	PR fortran/37203

gcc/fortran/expr.c

@@ -2127,8 +2127,15 @@ check_transformational (gfc_expr *e)
     "selected_real_kind", "transfer", "trim", NULL
   };
 
+  static const char * const trans_func_f2003[] =  {
+    "dot_product", "matmul", "null", "pack", "repeat",
+    "reshape", "selected_char_kind", "selected_int_kind",
+    "selected_real_kind", "transfer", "transpose", "trim", NULL
+  };
+
   int i;
   const char *name;
+  const char *const *functions;
 
   if (!e->value.function.isym
       || !e->value.function.isym->transformational)
@@ -2136,31 +2143,23 @@ check_transformational (gfc_expr *e)
 
   name = e->symtree->n.sym->name;
 
+  functions = (gfc_option.allow_std & GFC_STD_F2003) 
+		? trans_func_f2003 : trans_func_f95;
+
   /* NULL() is dealt with below.  */
   if (strcmp ("null", name) == 0)
     return MATCH_NO;
 
-  for (i = 0; trans_func_f95[i]; i++)
-    if (strcmp (trans_func_f95[i], name) == 0)
-      break;
+  for (i = 0; functions[i]; i++)
+    if (strcmp (functions[i], name) == 0)
+       break;
 
-  /* FIXME, F2003: implement translation of initialization
-     expressions before enabling this check. For F95, error
-     out if the transformational function is not in the list.  */
-#if 0
-  if (trans_func_f95[i] == NULL
-      && gfc_notify_std (GFC_STD_F2003, 
-			 "transformational intrinsic '%s' at %L is not permitted "
-			 "in an initialization expression", name, &e->where) == FAILURE)
-    return MATCH_ERROR;
-#else
-  if (trans_func_f95[i] == NULL)
+  if (functions[i] == NULL)
     {
       gfc_error("transformational intrinsic '%s' at %L is not permitted "
 		"in an initialization expression", name, &e->where);
       return MATCH_ERROR;
     }
-#endif
 
   return check_init_expr_arguments (e);
 }

gcc/fortran/intrinsic.c

@@ -1499,7 +1499,7 @@ add_functions (void)
   make_generic ("dim", GFC_ISYM_DIM, GFC_STD_F77);
 
   add_sym_2 ("dot_product", GFC_ISYM_DOT_PRODUCT, CLASS_TRANSFORMATIONAL, ACTUAL_NO, BT_REAL, dr,
-	     GFC_STD_F95, gfc_check_dot_product, NULL, gfc_resolve_dot_product,
+	     GFC_STD_F95, gfc_check_dot_product, gfc_simplify_dot_product, gfc_resolve_dot_product,
 	     va, BT_REAL, dr, REQUIRED, vb, BT_REAL, dr, REQUIRED);
 
   make_generic ("dot_product", GFC_ISYM_DOT_PRODUCT, GFC_STD_F95);
@@ -2034,7 +2034,7 @@ add_functions (void)
   make_generic ("malloc", GFC_ISYM_MALLOC, GFC_STD_GNU);
 
   add_sym_2 ("matmul", GFC_ISYM_MATMUL, CLASS_TRANSFORMATIONAL, ACTUAL_NO, BT_REAL, dr, GFC_STD_F95,
-	     gfc_check_matmul, NULL, gfc_resolve_matmul,
+	     gfc_check_matmul, gfc_simplify_matmul, gfc_resolve_matmul,
 	     ma, BT_REAL, dr, REQUIRED, mb, BT_REAL, dr, REQUIRED);
 
   make_generic ("matmul", GFC_ISYM_MATMUL, GFC_STD_F95);
@@ -2535,7 +2535,7 @@ add_functions (void)
   make_generic ("transfer", GFC_ISYM_TRANSFER, GFC_STD_F95);
 
   add_sym_1 ("transpose", GFC_ISYM_TRANSPOSE, CLASS_TRANSFORMATIONAL, ACTUAL_NO, BT_REAL, dr, GFC_STD_F95,
-	     gfc_check_transpose, NULL, gfc_resolve_transpose,
+	     gfc_check_transpose, gfc_simplify_transpose, gfc_resolve_transpose,
 	     m, BT_REAL, dr, REQUIRED);
 
   make_generic ("transpose", GFC_ISYM_TRANSPOSE, GFC_STD_F95);

gcc/fortran/intrinsic.h

@@ -229,6 +229,7 @@ gfc_expr *gfc_simplify_dble (gfc_expr *);
 gfc_expr *gfc_simplify_digits (gfc_expr *);
 gfc_expr *gfc_simplify_dim (gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_dprod (gfc_expr *, gfc_expr *);
+gfc_expr *gfc_simplify_dot_product (gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_epsilon (gfc_expr *);
 gfc_expr *gfc_simplify_erf (gfc_expr *);
 gfc_expr *gfc_simplify_erfc (gfc_expr *);
@@ -271,6 +272,7 @@ gfc_expr *gfc_simplify_llt (gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_log (gfc_expr *);
 gfc_expr *gfc_simplify_log10 (gfc_expr *);
 gfc_expr *gfc_simplify_logical (gfc_expr *, gfc_expr *);
+gfc_expr *gfc_simplify_matmul (gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_merge (gfc_expr *, gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_min (gfc_expr *);
 gfc_expr *gfc_simplify_minval (gfc_expr *, gfc_expr*, gfc_expr*);
@@ -318,6 +320,7 @@ gfc_expr *gfc_simplify_tanh (gfc_expr *);
 gfc_expr *gfc_simplify_tiny (gfc_expr *);
 gfc_expr *gfc_simplify_trailz (gfc_expr *);
 gfc_expr *gfc_simplify_transfer (gfc_expr *, gfc_expr *, gfc_expr *);
+gfc_expr *gfc_simplify_transpose (gfc_expr *);
 gfc_expr *gfc_simplify_trim (gfc_expr *);
 gfc_expr *gfc_simplify_ubound (gfc_expr *, gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_verify (gfc_expr *, gfc_expr *, gfc_expr *, gfc_expr *);

gcc/fortran/simplify.c

@@ -255,6 +255,138 @@ is_constant_array_expr (gfc_expr *e)
 }
 
 
+/* Initialize a transformational result expression with a given value.  */
+
+static void
+init_result_expr (gfc_expr *e, int init, gfc_expr *array)
+{
+  if (e && e->expr_type == EXPR_ARRAY)
+    {
+      gfc_constructor *ctor = e->value.constructor;
+      while (ctor)
+	{
+	  init_result_expr (ctor->expr, init, array);
+	  ctor = ctor->next;
+	}
+    }
+  else if (e && e->expr_type == EXPR_CONSTANT)
+    {
+      int i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
+      int length;
+      gfc_char_t *string;
+
+      switch (e->ts.type)
+	{
+	  case BT_LOGICAL:
+	    e->value.logical = (init ? 1 : 0);
+	    break;
+
+	  case BT_INTEGER:
+	    if (init == INT_MIN)
+	      mpz_set (e->value.integer, gfc_integer_kinds[i].min_int);
+	    else if (init == INT_MAX)
+	      mpz_set (e->value.integer, gfc_integer_kinds[i].huge);
+	    else
+	      mpz_set_si (e->value.integer, init);
+	    break;
+
+	  case BT_REAL:
+	    if (init == INT_MIN)
+	      {
+		mpfr_set (e->value.real, gfc_real_kinds[i].huge, GFC_RND_MODE);
+		mpfr_neg (e->value.real, e->value.real, GFC_RND_MODE);
+	      }
+	    else if (init == INT_MAX)
+	      mpfr_set (e->value.real, gfc_real_kinds[i].huge, GFC_RND_MODE);
+	    else
+	      mpfr_set_si (e->value.real, init, GFC_RND_MODE);
+	    break;
+
+	  case BT_COMPLEX:
+	    mpfr_set_si (e->value.complex.r, init, GFC_RND_MODE);
+	    mpfr_set_si (e->value.complex.i, 0, GFC_RND_MODE);
+	    break;
+
+	  case BT_CHARACTER:
+	    if (init == INT_MIN)
+	      {
+		gfc_expr *len = gfc_simplify_len (array, NULL);
+		gfc_extract_int (len, &length);
+		string = gfc_get_wide_string (length + 1);
+		gfc_wide_memset (string, 0, length);
+	      }
+	    else if (init == INT_MAX)
+	      {
+		gfc_expr *len = gfc_simplify_len (array, NULL);
+		gfc_extract_int (len, &length);
+		string = gfc_get_wide_string (length + 1);
+		gfc_wide_memset (string, 255, length);
+	      }
+	    else
+	      {
+		length = 0;
+		string = gfc_get_wide_string (1);
+	      }
+
+	    string[length] = '\0';
+	    e->value.character.length = length;
+	    e->value.character.string = string;
+	    break;
+
+	  default:
+	    gcc_unreachable();
+	}
+    }
+  else
+    gcc_unreachable();
+}
+
+
+/* Helper function for gfc_simplify_dot_product() and gfc_simplify_matmul.  */
+
+static gfc_expr *
+compute_dot_product (gfc_constructor *ctor_a, int stride_a,
+		     gfc_constructor *ctor_b, int stride_b)
+{
+  gfc_expr *result;
+  gfc_expr *a = ctor_a->expr, *b = ctor_b->expr;
+
+  gcc_assert (gfc_compare_types (&a->ts, &b->ts));
+
+  result = gfc_constant_result (a->ts.type, a->ts.kind, &a->where);
+  init_result_expr (result, 0, NULL);
+
+  while (ctor_a && ctor_b)
+    {
+      /* Copying of expressions is required as operands are free'd
+	 by the gfc_arith routines.  */
+      switch (result->ts.type)
+	{
+	  case BT_LOGICAL:
+	    result = gfc_or (result,
+			     gfc_and (gfc_copy_expr (ctor_a->expr),
+				      gfc_copy_expr (ctor_b->expr)));
+	    break;
+
+	  case BT_INTEGER:
+	  case BT_REAL:
+	  case BT_COMPLEX:
+	    result = gfc_add (result,
+			      gfc_multiply (gfc_copy_expr (ctor_a->expr),
+					    gfc_copy_expr (ctor_b->expr)));
+	    break;
+
+	  default:
+	    gcc_unreachable();
+	}
+
+      ADVANCE (ctor_a, stride_a);
+      ADVANCE (ctor_b, stride_b);
+    }
+
+  return result;
+}
+
 /********************** Simplification functions *****************************/
 
 gfc_expr *
@@ -1210,6 +1342,32 @@ gfc_simplify_dim (gfc_expr *x, gfc_expr *y)
 }
 
 
+gfc_expr*
+gfc_simplify_dot_product (gfc_expr *vector_a, gfc_expr *vector_b)
+{
+  gfc_expr *result;
+
+  if (!is_constant_array_expr (vector_a)
+      || !is_constant_array_expr (vector_b))
+    return NULL;
+
+  gcc_assert (vector_a->rank == 1);
+  gcc_assert (vector_b->rank == 1);
+  gcc_assert (gfc_compare_types (&vector_a->ts, &vector_b->ts));
+
+  if (vector_a->value.constructor && vector_b->value.constructor)
+    return compute_dot_product (vector_a->value.constructor, 1,
+			        vector_b->value.constructor, 1);
+
+  /* Zero sized array ...  */
+  result = gfc_constant_result (vector_a->ts.type,
+				vector_a->ts.kind,
+				&vector_a->where);
+  init_result_expr (result, 0, NULL);
+  return result;
+}
+
+
 gfc_expr *
 gfc_simplify_dprod (gfc_expr *x, gfc_expr *y)
 {
@@ -2856,6 +3014,84 @@ gfc_simplify_logical (gfc_expr *e, gfc_expr *k)
 }
 
 
+gfc_expr*
+gfc_simplify_matmul (gfc_expr *matrix_a, gfc_expr *matrix_b)
+{
+  gfc_expr *result;
+  gfc_constructor *ma_ctor, *mb_ctor;
+  int row, result_rows, col, result_columns, stride_a, stride_b;
+
+  if (!is_constant_array_expr (matrix_a)
+      || !is_constant_array_expr (matrix_b))
+    return NULL;
+
+  gcc_assert (gfc_compare_types (&matrix_a->ts, &matrix_b->ts));
+  result = gfc_start_constructor (matrix_a->ts.type,
+				  matrix_a->ts.kind,
+				  &matrix_a->where);
+
+  if (matrix_a->rank == 1 && matrix_b->rank == 2)
+    {
+      result_rows = 1;
+      result_columns = mpz_get_si (matrix_b->shape[0]);
+      stride_a = 1;
+      stride_b = mpz_get_si (matrix_b->shape[0]);
+
+      result->rank = 1;
+      result->shape = gfc_get_shape (result->rank);
+      mpz_init_set_si (result->shape[0], result_columns);
+    }
+  else if (matrix_a->rank == 2 && matrix_b->rank == 1)
+    {
+      result_rows = mpz_get_si (matrix_b->shape[0]);
+      result_columns = 1;
+      stride_a = mpz_get_si (matrix_a->shape[0]);
+      stride_b = 1;
+
+      result->rank = 1;
+      result->shape = gfc_get_shape (result->rank);
+      mpz_init_set_si (result->shape[0], result_rows);
+    }
+  else if (matrix_a->rank == 2 && matrix_b->rank == 2)
+    {
+      result_rows = mpz_get_si (matrix_a->shape[0]);
+      result_columns = mpz_get_si (matrix_b->shape[1]);
+      stride_a = mpz_get_si (matrix_a->shape[1]);
+      stride_b = mpz_get_si (matrix_b->shape[0]);
+
+      result->rank = 2;
+      result->shape = gfc_get_shape (result->rank);
+      mpz_init_set_si (result->shape[0], result_rows);
+      mpz_init_set_si (result->shape[1], result_columns);
+    }
+  else
+    gcc_unreachable();
+
+  ma_ctor = matrix_a->value.constructor;
+  mb_ctor = matrix_b->value.constructor;
+
+  for (col = 0; col < result_columns; ++col)
+    {
+      ma_ctor = matrix_a->value.constructor;
+
+      for (row = 0; row < result_rows; ++row)
+	{
+	  gfc_expr *e;
+	  e = compute_dot_product (ma_ctor, stride_a,
+				   mb_ctor, 1);
+
+	  gfc_append_constructor (result, e);
+
+	  ADVANCE (ma_ctor, 1);
+	}
+
+      ADVANCE (mb_ctor, stride_b);
+    }
+
+  return result;
+}
+
+
 gfc_expr *
 gfc_simplify_merge (gfc_expr *tsource, gfc_expr *fsource, gfc_expr *mask)
 {
@@ -4761,6 +4997,47 @@ gfc_simplify_transfer (gfc_expr *source, gfc_expr *mold, gfc_expr *size)
 
 
 gfc_expr *
+gfc_simplify_transpose (gfc_expr *matrix)
+{
+  int i, matrix_rows;
+  gfc_expr *result;
+  gfc_constructor *matrix_ctor;
+
+  if (!is_constant_array_expr (matrix))
+    return NULL;
+
+  gcc_assert (matrix->rank == 2);
+
+  result = gfc_start_constructor (matrix->ts.type, matrix->ts.kind, &matrix->where);
+  result->rank = 2;
+  result->shape = gfc_get_shape (result->rank);
+  mpz_set (result->shape[0], matrix->shape[1]);
+  mpz_set (result->shape[1], matrix->shape[0]);
+
+  if (matrix->ts.type == BT_CHARACTER)
+    result->ts.cl = matrix->ts.cl;
+
+  matrix_rows = mpz_get_si (matrix->shape[0]);
+  matrix_ctor = matrix->value.constructor;
+  for (i = 0; i < matrix_rows; ++i)
+    {
+      gfc_constructor *column_ctor = matrix_ctor;
+      while (column_ctor)
+	{
+	  gfc_append_constructor (result, 
+				  gfc_copy_expr (column_ctor->expr));
+
+	  ADVANCE (column_ctor, matrix_rows);
+	}
+
+      ADVANCE (matrix_ctor, 1);
+    }
+
+  return result;
+}
+
+
+gfc_expr *
 gfc_simplify_trim (gfc_expr *e)
 {
   gfc_expr *result;

gcc/testsuite/ChangeLog

+2009-06-07  Daniel Franke  <franke.daniel@gmail.com>
+
+	PR fortran/25104
+	PR fortran/29962
+	* gfortran.dg/dot_product_1.f03: New.
+	* gfortran.dg/matmul_8.f03: New.
+	* gfortran.dg/transpose_3.f03: New.
+
 2009-06-06  Ian Lance Taylor  <iant@google.com>
 
 	* gcc.dg/Wunused-label-1.c: New test case.

gcc/testsuite/gfortran.dg/dot_product_1.f03

+! { dg-do "run" }
+! Transformational intrinsic DOT_PRODUCT as initialization expression.
+
+  INTEGER, PARAMETER :: n = 10
+  INTEGER, PARAMETER :: a(n) = 1
+  INTEGER, PARAMETER :: p = DOT_PRODUCT(a, a)
+  INTEGER, PARAMETER :: e = DOT_PRODUCT(SHAPE(1), SHAPE(1))
+
+  IF (p /= n) CALL abort()
+  IF (e /= 0) CALL abort()
+END

gcc/testsuite/gfortran.dg/matmul_8.f03

+! { dg-do "run" }
+! Transformational intrinsic MATMUL as initialization expression.
+
+  REAL, PARAMETER :: PI = 3.141592654, theta = PI/6.0
+
+  REAL, PARAMETER :: unity(2,2) = RESHAPE([1.0, 0.0, 0.0, 1.0], [2, 2])
+  REAL, PARAMETER :: m1(2,2)    = RESHAPE([COS(theta), SIN(theta), -SIN(theta), COS(theta)], [2, 2])
+  REAL, PARAMETER :: m2(2,2)    = RESHAPE([COS(theta), -SIN(theta), SIN(theta), COS(theta)], [2, 2])
+  REAL, PARAMETER :: m(2,2)     = MATMUL(m1, m2)
+
+  IF (ANY(ABS(m - unity) > EPSILON(0.0))) CALL abort()
+END

gcc/testsuite/gfortran.dg/transpose_3.f03

+! { dg-do "run" }
+! Transformational intrinsic TRANSPOSE as initialization expression.
+
+  INTEGER, PARAMETER :: n = 10
+  INTEGER, PARAMETER :: a(n,1) = RESHAPE([ (i, i = 1, n) ], [n, 1])
+  INTEGER, PARAMETER :: b(1,n) = TRANSPOSE(a)
+  INTEGER, PARAMETER :: c(n,1) = TRANSPOSE(b)
+
+  IF (ANY(c /= a)) CALL abort()
+END