OSDN Git Service

2010-02-20 Paul Thomas <pault@gcc.gnu.org>
[pf3gnuchains/gcc-fork.git] / gcc / fortran / simplify.c
index db28d36..8768cb6 100644 (file)
@@ -214,26 +214,6 @@ convert_mpz_to_signed (mpz_t x, int bitsize)
     }
 }
 
-/* Helper function to convert to/from mpfr_t & mpc_t and call the
-   supplied mpc function on the respective values.  */
-
-#ifdef HAVE_mpc
-static void
-call_mpc_func (mpfr_ptr result_re, mpfr_ptr result_im,
-              mpfr_srcptr input_re, mpfr_srcptr input_im,
-              int (*func)(mpc_ptr, mpc_srcptr, mpc_rnd_t))
-{
-  mpc_t c;
-  mpc_init2 (c, mpfr_get_default_prec());
-  mpc_set_fr_fr (c, input_re, input_im, GFC_MPC_RND_MODE);
-  func (c, c, GFC_MPC_RND_MODE);
-  mpfr_set (result_re, mpc_realref (c), GFC_RND_MODE);
-  mpfr_set (result_im, mpc_imagref (c), GFC_RND_MODE);
-  mpc_clear (c);
-}
-#endif
-
-
 /* Test that the expression is an constant array.  */
 
 static bool
@@ -303,8 +283,7 @@ init_result_expr (gfc_expr *e, int init, gfc_expr *array)
            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);
+           mpc_set_si (e->value.complex, init, GFC_MPC_RND_MODE);
            break;
 
          case BT_CHARACTER:
@@ -387,6 +366,246 @@ compute_dot_product (gfc_constructor *ctor_a, int stride_a,
   return result;
 }
 
+
+/* Build a result expression for transformational intrinsics, 
+   depending on DIM. */
+
+static gfc_expr *
+transformational_result (gfc_expr *array, gfc_expr *dim, bt type,
+                        int kind, locus* where)
+{
+  gfc_expr *result;
+  int i, nelem;
+
+  if (!dim || array->rank == 1)
+    return gfc_constant_result (type, kind, where);
+
+  result = gfc_start_constructor (type, kind, where);
+  result->shape = gfc_copy_shape_excluding (array->shape, array->rank, dim);
+  result->rank = array->rank - 1;
+
+  /* gfc_array_size() would count the number of elements in the constructor,
+     we have not built those yet.  */
+  nelem = 1;
+  for  (i = 0; i < result->rank; ++i)
+    nelem *= mpz_get_ui (result->shape[i]);
+
+  for (i = 0; i < nelem; ++i)
+    {
+      gfc_expr *e = gfc_constant_result (type, kind, where);
+      gfc_append_constructor (result, e);
+    }
+
+  return result;
+}
+
+
+typedef gfc_expr* (*transformational_op)(gfc_expr*, gfc_expr*);
+
+/* Wrapper function, implements 'op1 += 1'. Only called if MASK
+   of COUNT intrinsic is .TRUE..
+
+   Interface and implimentation mimics arith functions as
+   gfc_add, gfc_multiply, etc.  */
+
+static gfc_expr* gfc_count (gfc_expr *op1, gfc_expr *op2)
+{
+  gfc_expr *result;
+
+  gcc_assert (op1->ts.type == BT_INTEGER);
+  gcc_assert (op2->ts.type == BT_LOGICAL);
+  gcc_assert (op2->value.logical);
+
+  result = gfc_copy_expr (op1);
+  mpz_add_ui (result->value.integer, result->value.integer, 1);
+
+  gfc_free_expr (op1);
+  gfc_free_expr (op2);
+  return result;
+}
+
+
+/* Transforms an ARRAY with operation OP, according to MASK, to a
+   scalar RESULT. E.g. called if
+
+     REAL, PARAMETER :: array(n, m) = ...
+     REAL, PARAMETER :: s = SUM(array)
+
+  where OP == gfc_add().  */
+
+static gfc_expr *
+simplify_transformation_to_scalar (gfc_expr *result, gfc_expr *array, gfc_expr *mask,
+                                  transformational_op op)
+{
+  gfc_expr *a, *m;
+  gfc_constructor *array_ctor, *mask_ctor;
+
+  /* Shortcut for constant .FALSE. MASK.  */
+  if (mask
+      && mask->expr_type == EXPR_CONSTANT
+      && !mask->value.logical)
+    return result;
+
+  array_ctor = array->value.constructor;
+  mask_ctor = NULL;
+  if (mask && mask->expr_type == EXPR_ARRAY)
+    mask_ctor = mask->value.constructor;
+
+  while (array_ctor)
+    {
+      a = array_ctor->expr;
+      array_ctor = array_ctor->next;
+
+      /* A constant MASK equals .TRUE. here and can be ignored.  */
+      if (mask_ctor)
+       {
+         m = mask_ctor->expr;
+         mask_ctor = mask_ctor->next;
+         if (!m->value.logical)
+           continue;
+       }
+
+      result = op (result, gfc_copy_expr (a));
+    }
+
+  return result;
+}
+
+/* Transforms an ARRAY with operation OP, according to MASK, to an
+   array RESULT. E.g. called if
+
+     REAL, PARAMETER :: array(n, m) = ...
+     REAL, PARAMETER :: s(n) = PROD(array, DIM=1)
+
+  where OP == gfc_multiply().  */
+
+static gfc_expr *
+simplify_transformation_to_array (gfc_expr *result, gfc_expr *array, gfc_expr *dim,
+                                 gfc_expr *mask, transformational_op op)
+{
+  mpz_t size;
+  int done, i, n, arraysize, resultsize, dim_index, dim_extent, dim_stride;
+  gfc_expr **arrayvec, **resultvec, **base, **src, **dest;
+  gfc_constructor *array_ctor, *mask_ctor, *result_ctor;
+
+  int count[GFC_MAX_DIMENSIONS], extent[GFC_MAX_DIMENSIONS],
+      sstride[GFC_MAX_DIMENSIONS], dstride[GFC_MAX_DIMENSIONS],
+      tmpstride[GFC_MAX_DIMENSIONS];
+
+  /* Shortcut for constant .FALSE. MASK.  */
+  if (mask
+      && mask->expr_type == EXPR_CONSTANT
+      && !mask->value.logical)
+    return result;
+
+  /* Build an indexed table for array element expressions to minimize
+     linked-list traversal. Masked elements are set to NULL.  */
+  gfc_array_size (array, &size);
+  arraysize = mpz_get_ui (size);
+
+  arrayvec = (gfc_expr**) gfc_getmem (sizeof (gfc_expr*) * arraysize);
+
+  array_ctor = array->value.constructor;
+  mask_ctor = NULL;
+  if (mask && mask->expr_type == EXPR_ARRAY)
+    mask_ctor = mask->value.constructor;
+
+  for (i = 0; i < arraysize; ++i)
+    {
+      arrayvec[i] = array_ctor->expr;
+      array_ctor = array_ctor->next;
+
+      if (mask_ctor)
+       {
+         if (!mask_ctor->expr->value.logical)
+           arrayvec[i] = NULL;
+
+         mask_ctor = mask_ctor->next;
+       }
+    }
+
+  /* Same for the result expression.  */
+  gfc_array_size (result, &size);
+  resultsize = mpz_get_ui (size);
+  mpz_clear (size);
+
+  resultvec = (gfc_expr**) gfc_getmem (sizeof (gfc_expr*) * resultsize);
+  result_ctor = result->value.constructor;
+  for (i = 0; i < resultsize; ++i)
+    {
+      resultvec[i] = result_ctor->expr;
+      result_ctor = result_ctor->next;
+    }
+
+  gfc_extract_int (dim, &dim_index);
+  dim_index -= 1;               /* zero-base index */
+  dim_extent = 0;
+  dim_stride = 0;
+
+  for (i = 0, n = 0; i < array->rank; ++i)
+    {
+      count[i] = 0;
+      tmpstride[i] = (i == 0) ? 1 : tmpstride[i-1] * mpz_get_si (array->shape[i-1]);
+      if (i == dim_index)
+       {
+         dim_extent = mpz_get_si (array->shape[i]);
+         dim_stride = tmpstride[i];
+         continue;
+       }
+
+      extent[n] = mpz_get_si (array->shape[i]);
+      sstride[n] = tmpstride[i];
+      dstride[n] = (n == 0) ? 1 : dstride[n-1] * extent[n-1];
+      n += 1;
+    }
+
+  done = false;
+  base = arrayvec;
+  dest = resultvec;
+  while (!done)
+    {
+      for (src = base, n = 0; n < dim_extent; src += dim_stride, ++n)
+       if (*src)
+         *dest = op (*dest, gfc_copy_expr (*src));
+
+      count[0]++;
+      base += sstride[0];
+      dest += dstride[0];
+
+      n = 0;
+      while (!done && count[n] == extent[n])
+       {
+         count[n] = 0;
+         base -= sstride[n] * extent[n];
+         dest -= dstride[n] * extent[n];
+
+         n++;
+         if (n < result->rank)
+           {
+             count [n]++;
+             base += sstride[n];
+             dest += dstride[n];
+           }
+         else
+           done = true;
+       }
+    }
+
+  /* Place updated expression in result constructor.  */
+  result_ctor = result->value.constructor;
+  for (i = 0; i < resultsize; ++i)
+    {
+      result_ctor->expr = resultvec[i];
+      result_ctor = result_ctor->next;
+    }
+
+  gfc_free (arrayvec);
+  gfc_free (resultvec);
+  return result;
+}
+
+
+
 /********************** Simplification functions *****************************/
 
 gfc_expr *
@@ -420,8 +639,7 @@ gfc_simplify_abs (gfc_expr *e)
 
       gfc_set_model_kind (e->ts.kind);
 
-      mpfr_hypot (result->value.real, e->value.complex.r, 
-                 e->value.complex.i, GFC_RND_MODE);
+      mpc_abs (result->value.real, e->value.complex, GFC_RND_MODE);
       result = range_check (result, "CABS");
       break;
 
@@ -507,17 +725,27 @@ gfc_simplify_acos (gfc_expr *x)
   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)
+  switch (x->ts.type)
     {
-      gfc_error ("Argument of ACOS at %L must be between -1 and 1",
-                &x->where);
-      return &gfc_bad_expr;
+      case BT_REAL:
+       if (mpfr_cmp_si (x->value.real, 1) > 0
+           || mpfr_cmp_si (x->value.real, -1) < 0)
+         {
+           gfc_error ("Argument of ACOS at %L must be between -1 and 1",
+                      &x->where);
+           return &gfc_bad_expr;
+         }
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpfr_acos (result->value.real, x->value.real, GFC_RND_MODE);
+       break;
+      case BT_COMPLEX:
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpc_acos (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+       break;
+      default:
+       gfc_internal_error ("in gfc_simplify_acos(): Bad type");
     }
 
-  result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
-
-  mpfr_acos (result->value.real, x->value.real, GFC_RND_MODE);
 
   return range_check (result, "ACOS");
 }
@@ -530,16 +758,26 @@ gfc_simplify_acosh (gfc_expr *x)
   if (x->expr_type != EXPR_CONSTANT)
     return NULL;
 
-  if (mpfr_cmp_si (x->value.real, 1) < 0)
+  switch (x->ts.type)
     {
-      gfc_error ("Argument of ACOSH at %L must not be less than 1",
-                &x->where);
-      return &gfc_bad_expr;
-    }
-
-  result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+      case BT_REAL:
+       if (mpfr_cmp_si (x->value.real, 1) < 0)
+         {
+           gfc_error ("Argument of ACOSH at %L must not be less than 1",
+                      &x->where);
+           return &gfc_bad_expr;
+         }
 
-  mpfr_acosh (result->value.real, x->value.real, GFC_RND_MODE);
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpfr_acosh (result->value.real, x->value.real, GFC_RND_MODE);
+       break;
+      case BT_COMPLEX:
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpc_acosh (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+       break;
+      default:
+       gfc_internal_error ("in gfc_simplify_acosh(): Bad type");
+    }
 
   return range_check (result, "ACOSH");
 }
@@ -627,7 +865,7 @@ gfc_simplify_aimag (gfc_expr *e)
     return NULL;
 
   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
-  mpfr_set (result->value.real, e->value.complex.i, GFC_RND_MODE);
+  mpfr_set (result->value.real, mpc_imagref (e->value.complex), GFC_RND_MODE);
 
   return range_check (result, "AIMAG");
 }
@@ -658,6 +896,25 @@ gfc_simplify_aint (gfc_expr *e, gfc_expr *k)
 
 
 gfc_expr *
+gfc_simplify_all (gfc_expr *mask, gfc_expr *dim)
+{
+  gfc_expr *result;
+
+  if (!is_constant_array_expr (mask)
+      || !gfc_is_constant_expr (dim))
+    return NULL;
+
+  result = transformational_result (mask, dim, mask->ts.type,
+                                   mask->ts.kind, &mask->where);
+  init_result_expr (result, true, NULL);
+
+  return !dim || mask->rank == 1 ?
+    simplify_transformation_to_scalar (result, mask, NULL, gfc_and) :
+    simplify_transformation_to_array (result, mask, dim, NULL, gfc_and);
+}
+
+
+gfc_expr *
 gfc_simplify_dint (gfc_expr *e)
 {
   gfc_expr *rtrunc, *result;
@@ -723,6 +980,25 @@ gfc_simplify_and (gfc_expr *x, gfc_expr *y)
 
 
 gfc_expr *
+gfc_simplify_any (gfc_expr *mask, gfc_expr *dim)
+{
+  gfc_expr *result;
+
+  if (!is_constant_array_expr (mask)
+      || !gfc_is_constant_expr (dim))
+    return NULL;
+
+  result = transformational_result (mask, dim, mask->ts.type,
+                                   mask->ts.kind, &mask->where);
+  init_result_expr (result, false, NULL);
+
+  return !dim || mask->rank == 1 ?
+    simplify_transformation_to_scalar (result, mask, NULL, gfc_or) :
+    simplify_transformation_to_array (result, mask, dim, NULL, gfc_or);
+}
+
+
+gfc_expr *
 gfc_simplify_dnint (gfc_expr *e)
 {
   gfc_expr *result;
@@ -746,18 +1022,27 @@ gfc_simplify_asin (gfc_expr *x)
   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)
+  switch (x->ts.type)
     {
-      gfc_error ("Argument of ASIN at %L must be between -1 and 1",
-                &x->where);
-      return &gfc_bad_expr;
+      case BT_REAL:
+       if (mpfr_cmp_si (x->value.real, 1) > 0
+           || mpfr_cmp_si (x->value.real, -1) < 0)
+         {
+           gfc_error ("Argument of ASIN at %L must be between -1 and 1",
+                      &x->where);
+           return &gfc_bad_expr;
+         }
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpfr_asin (result->value.real, x->value.real, GFC_RND_MODE);
+       break;
+      case BT_COMPLEX:
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpc_asin (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+       break;
+      default:
+       gfc_internal_error ("in gfc_simplify_asin(): Bad type");
     }
 
-  result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
-
-  mpfr_asin (result->value.real, x->value.real, GFC_RND_MODE);
-
   return range_check (result, "ASIN");
 }
 
@@ -770,9 +1055,19 @@ gfc_simplify_asinh (gfc_expr *x)
   if (x->expr_type != EXPR_CONSTANT)
     return NULL;
 
-  result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
-
-  mpfr_asinh (result->value.real, x->value.real, GFC_RND_MODE);
+  switch (x->ts.type)
+    {
+      case BT_REAL:
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpfr_asinh (result->value.real, x->value.real, GFC_RND_MODE);
+       break;
+      case BT_COMPLEX:
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpc_asinh (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+       break;
+      default:
+       gfc_internal_error ("in gfc_simplify_asinh(): Bad type");
+    }
 
   return range_check (result, "ASINH");
 }
@@ -786,9 +1081,19 @@ gfc_simplify_atan (gfc_expr *x)
   if (x->expr_type != EXPR_CONSTANT)
     return NULL;
     
-  result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
-
-  mpfr_atan (result->value.real, x->value.real, GFC_RND_MODE);
+  switch (x->ts.type)
+    {
+      case BT_REAL:
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpfr_atan (result->value.real, x->value.real, GFC_RND_MODE);
+       break;
+      case BT_COMPLEX:
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpc_atan (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+       break;
+      default:
+       gfc_internal_error ("in gfc_simplify_atan(): Bad type");
+    }
 
   return range_check (result, "ATAN");
 }
@@ -802,17 +1107,27 @@ gfc_simplify_atanh (gfc_expr *x)
   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)
+  switch (x->ts.type)
     {
-      gfc_error ("Argument of ATANH at %L must be inside the range -1 to 1",
-                &x->where);
-      return &gfc_bad_expr;
-    }
-
-  result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+      case BT_REAL:
+       if (mpfr_cmp_si (x->value.real, 1) >= 0
+           || mpfr_cmp_si (x->value.real, -1) <= 0)
+         {
+           gfc_error ("Argument of ATANH at %L must be inside the range -1 "
+                      "to 1", &x->where);
+           return &gfc_bad_expr;
+         }
 
-  mpfr_atanh (result->value.real, x->value.real, GFC_RND_MODE);
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpfr_atanh (result->value.real, x->value.real, GFC_RND_MODE);
+       break;
+      case BT_COMPLEX:
+       result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+       mpc_atanh (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+       break;
+      default:
+       gfc_internal_error ("in gfc_simplify_atanh(): Bad type");
+    }
 
   return range_check (result, "ATANH");
 }
@@ -1008,22 +1323,19 @@ simplify_cmplx (const char *name, gfc_expr *x, gfc_expr *y, int kind)
 
   result = gfc_constant_result (BT_COMPLEX, kind, &x->where);
 
-  mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
-
   switch (x->ts.type)
     {
     case BT_INTEGER:
       if (!x->is_boz)
-       mpfr_set_z (result->value.complex.r, x->value.integer, GFC_RND_MODE);
+       mpc_set_z (result->value.complex, x->value.integer, GFC_MPC_RND_MODE);
       break;
 
     case BT_REAL:
-      mpfr_set (result->value.complex.r, x->value.real, GFC_RND_MODE);
+      mpc_set_fr (result->value.complex, x->value.real, GFC_RND_MODE);
       break;
 
     case BT_COMPLEX:
-      mpfr_set (result->value.complex.r, x->value.complex.r, GFC_RND_MODE);
-      mpfr_set (result->value.complex.i, x->value.complex.i, GFC_RND_MODE);
+      mpc_set (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
       break;
 
     default:
@@ -1036,12 +1348,13 @@ simplify_cmplx (const char *name, gfc_expr *x, gfc_expr *y, int kind)
        {
        case BT_INTEGER:
          if (!y->is_boz)
-           mpfr_set_z (result->value.complex.i, y->value.integer,
-                       GFC_RND_MODE);
+           mpfr_set_z (mpc_imagref (result->value.complex),
+                       y->value.integer, GFC_RND_MODE);
          break;
 
        case BT_REAL:
-         mpfr_set (result->value.complex.i, y->value.real, GFC_RND_MODE);
+         mpfr_set (mpc_imagref (result->value.complex),
+                   y->value.real, GFC_RND_MODE);
          break;
 
        default:
@@ -1058,7 +1371,8 @@ simplify_cmplx (const char *name, gfc_expr *x, gfc_expr *y, int kind)
       ts.type = BT_REAL;
       if (!gfc_convert_boz (x, &ts))
        return &gfc_bad_expr;
-      mpfr_set (result->value.complex.r, x->value.real, GFC_RND_MODE);
+      mpfr_set (mpc_realref (result->value.complex),
+               x->value.real, GFC_RND_MODE);
     }
 
   if (y && y->is_boz)
@@ -1069,7 +1383,8 @@ simplify_cmplx (const char *name, gfc_expr *x, gfc_expr *y, int kind)
       ts.type = BT_REAL;
       if (!gfc_convert_boz (y, &ts))
        return &gfc_bad_expr;
-      mpfr_set (result->value.complex.i, y->value.real, GFC_RND_MODE);
+      mpfr_set (mpc_imagref (result->value.complex),
+               y->value.real, GFC_RND_MODE);
     }
 
   return range_check (result, name);
@@ -1151,8 +1466,7 @@ gfc_simplify_conjg (gfc_expr *e)
     return NULL;
 
   result = gfc_copy_expr (e);
-  mpfr_neg (result->value.complex.i, result->value.complex.i, GFC_RND_MODE);
-
+  mpc_conj (result->value.complex, result->value.complex, GFC_MPC_RND_MODE);
   return range_check (result, "CONJG");
 }
 
@@ -1174,27 +1488,7 @@ gfc_simplify_cos (gfc_expr *x)
       break;
     case BT_COMPLEX:
       gfc_set_model_kind (x->ts.kind);
-#ifdef HAVE_mpc
-      call_mpc_func (result->value.complex.r, result->value.complex.i,
-                    x->value.complex.r, x->value.complex.i, mpc_cos);
-#else
-    {
-      mpfr_t xp, xq;
-      mpfr_init (xp);
-      mpfr_init (xq);
-
-      mpfr_cos  (xp, x->value.complex.r, GFC_RND_MODE);
-      mpfr_cosh (xq, x->value.complex.i, GFC_RND_MODE);
-      mpfr_mul (result->value.complex.r, xp, xq, GFC_RND_MODE);
-
-      mpfr_sin  (xp, x->value.complex.r, GFC_RND_MODE);
-      mpfr_sinh (xq, x->value.complex.i, GFC_RND_MODE);
-      mpfr_mul (xp, xp, xq, GFC_RND_MODE);
-      mpfr_neg (result->value.complex.i, xp, GFC_RND_MODE );
-
-      mpfr_clears (xp, xq, NULL);
-    }
-#endif
+      mpc_cos (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
       break;
     default:
       gfc_internal_error ("in gfc_simplify_cos(): Bad type");
@@ -1215,13 +1509,44 @@ gfc_simplify_cosh (gfc_expr *x)
 
   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 
-  mpfr_cosh (result->value.real, x->value.real, GFC_RND_MODE);
+  if (x->ts.type == BT_REAL)
+    mpfr_cosh (result->value.real, x->value.real, GFC_RND_MODE);
+  else if (x->ts.type == BT_COMPLEX)
+    mpc_cosh (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+  else
+    gcc_unreachable ();
 
   return range_check (result, "COSH");
 }
 
 
 gfc_expr *
+gfc_simplify_count (gfc_expr *mask, gfc_expr *dim, gfc_expr *kind)
+{
+  gfc_expr *result;
+
+  if (!is_constant_array_expr (mask)
+      || !gfc_is_constant_expr (dim)
+      || !gfc_is_constant_expr (kind))
+    return NULL;
+
+  result = transformational_result (mask, dim,
+                                   BT_INTEGER,
+                                   get_kind (BT_INTEGER, kind, "COUNT",
+                                             gfc_default_integer_kind),
+                                   &mask->where);
+
+  init_result_expr (result, 0, NULL);
+
+  /* Passing MASK twice, once as data array, once as mask.
+     Whenever gfc_count is called, '1' is added to the result.  */
+  return !dim || mask->rank == 1 ?
+    simplify_transformation_to_scalar (result, mask, mask, gfc_count) :
+    simplify_transformation_to_array (result, mask, dim, mask, gfc_count);
+}
+
+
+gfc_expr *
 gfc_simplify_dcmplx (gfc_expr *x, gfc_expr *y)
 {
 
@@ -1593,22 +1918,7 @@ gfc_simplify_exp (gfc_expr *x)
 
     case BT_COMPLEX:
       gfc_set_model_kind (x->ts.kind);
-#ifdef HAVE_mpc
-      call_mpc_func (result->value.complex.r, result->value.complex.i,
-                    x->value.complex.r, x->value.complex.i, mpc_exp);
-#else
-    {
-      mpfr_t xp, xq;
-      mpfr_init (xp);
-      mpfr_init (xq);
-      mpfr_exp (xq, x->value.complex.r, GFC_RND_MODE);
-      mpfr_cos (xp, x->value.complex.i, GFC_RND_MODE);
-      mpfr_mul (result->value.complex.r, xq, xp, GFC_RND_MODE);
-      mpfr_sin (xp, x->value.complex.i, GFC_RND_MODE);
-      mpfr_mul (result->value.complex.i, xq, xp, GFC_RND_MODE);
-      mpfr_clears (xp, xq, NULL);
-    }
-#endif
+      mpc_exp (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
       break;
 
     default:
@@ -2322,6 +2632,54 @@ gfc_simplify_ior (gfc_expr *x, gfc_expr *y)
 
 
 gfc_expr *
+gfc_simplify_is_iostat_end (gfc_expr *x)
+{
+  gfc_expr *result;
+
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
+                               &x->where);
+  result->value.logical = (mpz_cmp_si (x->value.integer, LIBERROR_END) == 0);
+
+  return result;
+}
+
+
+gfc_expr *
+gfc_simplify_is_iostat_eor (gfc_expr *x)
+{
+  gfc_expr *result;
+
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
+                               &x->where);
+  result->value.logical = (mpz_cmp_si (x->value.integer, LIBERROR_EOR) == 0);
+
+  return result;
+}
+
+
+gfc_expr *
+gfc_simplify_isnan (gfc_expr *x)
+{
+  gfc_expr *result;
+
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
+                               &x->where);
+  result->value.logical = mpfr_nan_p (x->value.real);
+
+  return result;
+}
+
+
+gfc_expr *
 gfc_simplify_ishft (gfc_expr *e, gfc_expr *s)
 {
   gfc_expr *result;
@@ -2798,12 +3156,12 @@ gfc_simplify_len (gfc_expr *e, gfc_expr *kind)
        }
     }
 
-  if (e->ts.cl != NULL && e->ts.cl->length != NULL
-      && e->ts.cl->length->expr_type == EXPR_CONSTANT
-      && e->ts.cl->length->ts.type == BT_INTEGER)
+  if (e->ts.u.cl != NULL && e->ts.u.cl->length != NULL
+      && e->ts.u.cl->length->expr_type == EXPR_CONSTANT
+      && e->ts.u.cl->length->ts.type == BT_INTEGER)
     {
       result = gfc_constant_result (BT_INTEGER, k, &e->where);
-      mpz_set (result->value.integer, e->ts.cl->length->value.integer);
+      mpz_set (result->value.integer, e->ts.u.cl->length->value.integer);
       if (gfc_range_check (result) == ARITH_OK)
        return result;
       else
@@ -2929,8 +3287,8 @@ gfc_simplify_log (gfc_expr *x)
       break;
 
     case BT_COMPLEX:
-      if ((mpfr_sgn (x->value.complex.r) == 0)
-         && (mpfr_sgn (x->value.complex.i) == 0))
+      if ((mpfr_sgn (mpc_realref (x->value.complex)) == 0)
+         && (mpfr_sgn (mpc_imagref (x->value.complex)) == 0))
        {
          gfc_error ("Complex argument of LOG at %L cannot be zero",
                     &x->where);
@@ -2939,27 +3297,7 @@ gfc_simplify_log (gfc_expr *x)
        }
 
       gfc_set_model_kind (x->ts.kind);
-#ifdef HAVE_mpc
-      call_mpc_func (result->value.complex.r, result->value.complex.i,
-                    x->value.complex.r, x->value.complex.i, mpc_log);
-#else
-    {
-      mpfr_t xr, xi;
-      mpfr_init (xr);
-      mpfr_init (xi);
-
-      mpfr_atan2 (result->value.complex.i, x->value.complex.i,
-                 x->value.complex.r, GFC_RND_MODE);
-
-      mpfr_mul (xr, x->value.complex.r, x->value.complex.r, GFC_RND_MODE);
-      mpfr_mul (xi, x->value.complex.i, x->value.complex.i, GFC_RND_MODE);
-      mpfr_add (xr, xr, xi, GFC_RND_MODE);
-      mpfr_sqrt (xr, xr, GFC_RND_MODE);
-      mpfr_log (result->value.complex.r, xr, GFC_RND_MODE);
-
-      mpfr_clears (xr, xi, NULL);
-    }
-#endif
+      mpc_log (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
       break;
 
     default:
@@ -3684,7 +4022,7 @@ gfc_simplify_pack (gfc_expr *array, gfc_expr *mask, gfc_expr *vector)
   gfc_array_size (result, &result->shape[0]);
 
   if (array->ts.type == BT_CHARACTER)
-    result->ts.cl = array->ts.cl;
+    result->ts.u.cl = array->ts.u.cl;
 
   return result;
 }
@@ -3706,6 +4044,30 @@ gfc_simplify_precision (gfc_expr *e)
 
 
 gfc_expr *
+gfc_simplify_product (gfc_expr *array, gfc_expr *dim, gfc_expr *mask)
+{
+  gfc_expr *result;
+
+  if (!is_constant_array_expr (array)
+      || !gfc_is_constant_expr (dim))
+    return NULL;
+
+  if (mask
+      && !is_constant_array_expr (mask)
+      && mask->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  result = transformational_result (array, dim, array->ts.type,
+                                   array->ts.kind, &array->where);
+  init_result_expr (result, 1, NULL);
+
+  return !dim || array->rank == 1 ?
+    simplify_transformation_to_scalar (result, array, mask, gfc_multiply) :
+    simplify_transformation_to_array (result, array, dim, mask, gfc_multiply);
+}
+
+
+gfc_expr *
 gfc_simplify_radix (gfc_expr *e)
 {
   gfc_expr *result;
@@ -3828,8 +4190,7 @@ gfc_simplify_realpart (gfc_expr *e)
     return NULL;
 
   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
-  mpfr_set (result->value.real, e->value.complex.r, GFC_RND_MODE);
-
+  mpc_real (result->value.real, e->value.complex, GFC_RND_MODE);
   return range_check (result, "REALPART");
 }
 
@@ -3854,14 +4215,14 @@ gfc_simplify_repeat (gfc_expr *e, gfc_expr *n)
     }
 
   /* If we don't know the character length, we can do no more.  */
-  if (e->ts.cl && e->ts.cl->length
-       && e->ts.cl->length->expr_type == EXPR_CONSTANT)
+  if (e->ts.u.cl && e->ts.u.cl->length
+       && e->ts.u.cl->length->expr_type == EXPR_CONSTANT)
     {
-      len = mpz_get_si (e->ts.cl->length->value.integer);
+      len = mpz_get_si (e->ts.u.cl->length->value.integer);
       have_length = true;
     }
   else if (e->expr_type == EXPR_CONSTANT
-            && (e->ts.cl == NULL || e->ts.cl->length == NULL))
+            && (e->ts.u.cl == NULL || e->ts.u.cl->length == NULL))
     {
       len = e->value.character.length;
     }
@@ -3889,7 +4250,7 @@ gfc_simplify_repeat (gfc_expr *e, gfc_expr *n)
       if (have_length)
        {
          mpz_tdiv_q (max, gfc_integer_kinds[i].huge,
-                     e->ts.cl->length->value.integer);
+                     e->ts.u.cl->length->value.integer);
        }
       else
        {
@@ -3918,8 +4279,8 @@ gfc_simplify_repeat (gfc_expr *e, gfc_expr *n)
     return NULL;
 
   if (len || 
-      (e->ts.cl->length && 
-       mpz_sgn (e->ts.cl->length->value.integer)) != 0)
+      (e->ts.u.cl->length && 
+       mpz_sgn (e->ts.u.cl->length->value.integer)) != 0)
     {
       const char *res = gfc_extract_int (n, &ncop);
       gcc_assert (res == NULL);
@@ -4571,16 +4932,15 @@ gfc_simplify_sign (gfc_expr *x, gfc_expr *y)
       mpz_abs (result->value.integer, x->value.integer);
       if (mpz_sgn (y->value.integer) < 0)
        mpz_neg (result->value.integer, result->value.integer);
-
       break;
 
     case BT_REAL:
-      /* TODO: Handle -0.0 and +0.0 correctly on machines that support
-        it.  */
-      mpfr_abs (result->value.real, x->value.real, GFC_RND_MODE);
-      if (mpfr_sgn (y->value.real) < 0)
-       mpfr_neg (result->value.real, result->value.real, GFC_RND_MODE);
-
+      if (gfc_option.flag_sign_zero)
+       mpfr_copysign (result->value.real, x->value.real, y->value.real,
+                      GFC_RND_MODE);
+      else
+       mpfr_setsign (result->value.real, x->value.real,
+                     mpfr_sgn (y->value.real) < 0 ? 1 : 0, GFC_RND_MODE);
       break;
 
     default:
@@ -4609,26 +4969,7 @@ gfc_simplify_sin (gfc_expr *x)
 
     case BT_COMPLEX:
       gfc_set_model (x->value.real);
-#ifdef HAVE_mpc
-      call_mpc_func (result->value.complex.r, result->value.complex.i,
-                    x->value.complex.r, x->value.complex.i, mpc_sin);
-#else
-    {
-      mpfr_t xp, xq;
-      mpfr_init (xp);
-      mpfr_init (xq);
-
-      mpfr_sin  (xp, x->value.complex.r, GFC_RND_MODE);
-      mpfr_cosh (xq, x->value.complex.i, GFC_RND_MODE);
-      mpfr_mul (result->value.complex.r, xp, xq, GFC_RND_MODE);
-
-      mpfr_cos  (xp, x->value.complex.r, GFC_RND_MODE);
-      mpfr_sinh (xq, x->value.complex.i, GFC_RND_MODE);
-      mpfr_mul (result->value.complex.i, xp, xq, GFC_RND_MODE);
-
-      mpfr_clears (xp, xq, NULL);
-    }
-#endif
+      mpc_sin (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
       break;
 
     default:
@@ -4649,7 +4990,13 @@ gfc_simplify_sinh (gfc_expr *x)
 
   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 
-  mpfr_sinh (result->value.real, x->value.real, GFC_RND_MODE);
+  if (x->ts.type == BT_REAL)
+    mpfr_sinh (result->value.real, x->value.real, GFC_RND_MODE);
+  else if (x->ts.type == BT_COMPLEX)
+    mpc_sinh (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+  else
+    gcc_unreachable ();
+
 
   return range_check (result, "SINH");
 }
@@ -4710,6 +5057,113 @@ gfc_simplify_spacing (gfc_expr *x)
 
 
 gfc_expr *
+gfc_simplify_spread (gfc_expr *source, gfc_expr *dim_expr, gfc_expr *ncopies_expr)
+{
+  gfc_expr *result = 0L;
+  int i, j, dim, ncopies;
+  mpz_t size;
+
+  if ((!gfc_is_constant_expr (source)
+       && !is_constant_array_expr (source))
+      || !gfc_is_constant_expr (dim_expr)
+      || !gfc_is_constant_expr (ncopies_expr))
+    return NULL;
+
+  gcc_assert (dim_expr->ts.type == BT_INTEGER);
+  gfc_extract_int (dim_expr, &dim);
+  dim -= 1;   /* zero-base DIM */
+
+  gcc_assert (ncopies_expr->ts.type == BT_INTEGER);
+  gfc_extract_int (ncopies_expr, &ncopies);
+  ncopies = MAX (ncopies, 0);
+
+  /* Do not allow the array size to exceed the limit for an array
+     constructor.  */
+  if (source->expr_type == EXPR_ARRAY)
+    {
+      if (gfc_array_size (source, &size) == FAILURE)
+       gfc_internal_error ("Failure getting length of a constant array.");
+    }
+  else
+    mpz_init_set_ui (size, 1);
+
+  if (mpz_get_si (size)*ncopies > gfc_option.flag_max_array_constructor)
+    return NULL;
+
+  if (source->expr_type == EXPR_CONSTANT)
+    {
+      gcc_assert (dim == 0);
+
+      result = gfc_start_constructor (source->ts.type,
+                                     source->ts.kind,
+                                     &source->where);
+      result->rank = 1;
+      result->shape = gfc_get_shape (result->rank);
+      mpz_init_set_si (result->shape[0], ncopies);
+
+      for (i = 0; i < ncopies; ++i)
+        gfc_append_constructor (result, gfc_copy_expr (source));
+    }
+  else if (source->expr_type == EXPR_ARRAY)
+    {
+      int result_size, rstride[GFC_MAX_DIMENSIONS], extent[GFC_MAX_DIMENSIONS];
+      gfc_constructor *ctor, *source_ctor, *result_ctor;
+
+      gcc_assert (source->rank < GFC_MAX_DIMENSIONS);
+      gcc_assert (dim >= 0 && dim <= source->rank);
+
+      result = gfc_start_constructor (source->ts.type,
+                                     source->ts.kind,
+                                     &source->where);
+      result->rank = source->rank + 1;
+      result->shape = gfc_get_shape (result->rank);
+
+      result_size = 1;
+      for (i = 0, j = 0; i < result->rank; ++i)
+       {
+         if (i != dim)
+           mpz_init_set (result->shape[i], source->shape[j++]);
+         else
+           mpz_init_set_si (result->shape[i], ncopies);
+
+         extent[i] = mpz_get_si (result->shape[i]);
+         rstride[i] = (i == 0) ? 1 : rstride[i-1] * extent[i-1];
+         result_size *= extent[i];
+       }
+
+      for (i = 0; i < result_size; ++i)
+       gfc_append_constructor (result, NULL);
+
+      source_ctor = source->value.constructor;
+      result_ctor = result->value.constructor;
+      while (source_ctor)
+       {
+         ctor = result_ctor;
+
+         for (i = 0; i < ncopies; ++i)
+         {
+           ctor->expr = gfc_copy_expr (source_ctor->expr);
+           ADVANCE (ctor, rstride[dim]);
+         }
+
+         ADVANCE (result_ctor, (dim == 0 ? ncopies : 1));
+         ADVANCE (source_ctor, 1);
+       }
+    }
+  else
+    /* FIXME: Returning here avoids a regression in array_simplify_1.f90.
+       Replace NULL with gcc_unreachable() after implementing
+       gfc_simplify_cshift(). */
+    return NULL;
+
+  if (source->ts.type == BT_CHARACTER)
+    result->ts.u.cl = source->ts.u.cl;
+
+  return result;
+}
+
+
+gfc_expr *
 gfc_simplify_sqrt (gfc_expr *e)
 {
   gfc_expr *result;
@@ -4730,88 +5184,7 @@ gfc_simplify_sqrt (gfc_expr *e)
 
     case BT_COMPLEX:
       gfc_set_model (e->value.real);
-#ifdef HAVE_mpc
-      call_mpc_func (result->value.complex.r, result->value.complex.i,
-                    e->value.complex.r, e->value.complex.i, mpc_sqrt);
-#else
-    {
-      /* Formula taken from Numerical Recipes to avoid over- and
-        underflow.  */
-
-      mpfr_t ac, ad, s, t, w;
-      mpfr_init (ac);
-      mpfr_init (ad);
-      mpfr_init (s);
-      mpfr_init (t);
-      mpfr_init (w);
-
-      if (mpfr_cmp_ui (e->value.complex.r, 0) == 0
-         && mpfr_cmp_ui (e->value.complex.i, 0) == 0)
-       {
-         mpfr_set_ui (result->value.complex.r, 0, GFC_RND_MODE);
-         mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
-         break;
-       }
-
-      mpfr_abs (ac, e->value.complex.r, GFC_RND_MODE);
-      mpfr_abs (ad, e->value.complex.i, GFC_RND_MODE);
-
-      if (mpfr_cmp (ac, ad) >= 0)
-       {
-         mpfr_div (t, e->value.complex.i, e->value.complex.r, GFC_RND_MODE);
-         mpfr_mul (t, t, t, GFC_RND_MODE);
-         mpfr_add_ui (t, t, 1, GFC_RND_MODE);
-         mpfr_sqrt (t, t, GFC_RND_MODE);
-         mpfr_add_ui (t, t, 1, GFC_RND_MODE);
-         mpfr_div_ui (t, t, 2, GFC_RND_MODE);
-         mpfr_sqrt (t, t, GFC_RND_MODE);
-         mpfr_sqrt (s, ac, GFC_RND_MODE);
-         mpfr_mul (w, s, t, GFC_RND_MODE);
-       }
-      else
-       {
-         mpfr_div (s, e->value.complex.r, e->value.complex.i, GFC_RND_MODE);
-         mpfr_mul (t, s, s, GFC_RND_MODE);
-         mpfr_add_ui (t, t, 1, GFC_RND_MODE);
-         mpfr_sqrt (t, t, GFC_RND_MODE);
-         mpfr_abs (s, s, GFC_RND_MODE);
-         mpfr_add (t, t, s, GFC_RND_MODE);
-         mpfr_div_ui (t, t, 2, GFC_RND_MODE);
-         mpfr_sqrt (t, t, GFC_RND_MODE);
-         mpfr_sqrt (s, ad, GFC_RND_MODE);
-         mpfr_mul (w, s, t, GFC_RND_MODE);
-       }
-
-      if (mpfr_cmp_ui (w, 0) != 0 && mpfr_cmp_ui (e->value.complex.r, 0) >= 0)
-       {
-         mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
-         mpfr_div (result->value.complex.i, e->value.complex.i, t, GFC_RND_MODE);
-         mpfr_set (result->value.complex.r, w, GFC_RND_MODE);
-       }
-      else if (mpfr_cmp_ui (w, 0) != 0
-              && mpfr_cmp_ui (e->value.complex.r, 0) < 0
-              && mpfr_cmp_ui (e->value.complex.i, 0) >= 0)
-       {
-         mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
-         mpfr_div (result->value.complex.r, e->value.complex.i, t, GFC_RND_MODE);
-         mpfr_set (result->value.complex.i, w, GFC_RND_MODE);
-       }
-      else if (mpfr_cmp_ui (w, 0) != 0
-              && mpfr_cmp_ui (e->value.complex.r, 0) < 0
-              && mpfr_cmp_ui (e->value.complex.i, 0) < 0)
-       {
-         mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
-         mpfr_div (result->value.complex.r, ad, t, GFC_RND_MODE);
-         mpfr_neg (w, w, GFC_RND_MODE);
-         mpfr_set (result->value.complex.i, w, GFC_RND_MODE);
-       }
-      else
-       gfc_internal_error ("invalid complex argument of SQRT at %L",
-                           &e->where);
-
-      mpfr_clears (s, t, ac, ad, w, NULL);
-    }
-#endif
+      mpc_sqrt (result->value.complex, e->value.complex, GFC_MPC_RND_MODE);
       break;
 
     default:
@@ -4828,19 +5201,45 @@ negative_arg:
 
 
 gfc_expr *
+gfc_simplify_sum (gfc_expr *array, gfc_expr *dim, gfc_expr *mask)
+{
+  gfc_expr *result;
+
+  if (!is_constant_array_expr (array)
+      || !gfc_is_constant_expr (dim))
+    return NULL;
+
+  if (mask
+      && !is_constant_array_expr (mask)
+      && mask->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  result = transformational_result (array, dim, array->ts.type,
+                                   array->ts.kind, &array->where);
+  init_result_expr (result, 0, NULL);
+
+  return !dim || array->rank == 1 ?
+    simplify_transformation_to_scalar (result, array, mask, gfc_add) :
+    simplify_transformation_to_array (result, array, dim, mask, gfc_add);
+}
+
+
+gfc_expr *
 gfc_simplify_tan (gfc_expr *x)
 {
-  int i;
   gfc_expr *result;
 
   if (x->expr_type != EXPR_CONSTANT)
     return NULL;
 
-  i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
-
   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 
-  mpfr_tan (result->value.real, x->value.real, GFC_RND_MODE);
+  if (x->ts.type == BT_REAL)
+    mpfr_tan (result->value.real, x->value.real, GFC_RND_MODE);
+  else if (x->ts.type == BT_COMPLEX)
+    mpc_tan (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+  else
+    gcc_unreachable ();
 
   return range_check (result, "TAN");
 }
@@ -4856,7 +5255,12 @@ gfc_simplify_tanh (gfc_expr *x)
 
   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 
-  mpfr_tanh (result->value.real, x->value.real, GFC_RND_MODE);
+  if (x->ts.type == BT_REAL)
+    mpfr_tanh (result->value.real, x->value.real, GFC_RND_MODE);
+  else if (x->ts.type == BT_COMPLEX)
+    mpc_tanh (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
+  else
+    gcc_unreachable ();
 
   return range_check (result, "TANH");
 
@@ -5015,7 +5419,7 @@ gfc_simplify_transpose (gfc_expr *matrix)
   mpz_set (result->shape[1], matrix->shape[0]);
 
   if (matrix->ts.type == BT_CHARACTER)
-    result->ts.cl = matrix->ts.cl;
+    result->ts.u.cl = matrix->ts.u.cl;
 
   matrix_rows = mpz_get_si (matrix->shape[0]);
   matrix_ctor = matrix->value.constructor;
@@ -5080,6 +5484,54 @@ gfc_simplify_ubound (gfc_expr *array, gfc_expr *dim, gfc_expr *kind)
 
 
 gfc_expr *
+gfc_simplify_unpack (gfc_expr *vector, gfc_expr *mask, gfc_expr *field)
+{
+  gfc_expr *result, *e;
+  gfc_constructor *vector_ctor, *mask_ctor, *field_ctor;
+
+  if (!is_constant_array_expr (vector)
+      || !is_constant_array_expr (mask)
+      || (!gfc_is_constant_expr (field)
+         && !is_constant_array_expr(field)))
+    return NULL;
+
+  result = gfc_start_constructor (vector->ts.type,
+                                 vector->ts.kind,
+                                 &vector->where);
+  result->rank = mask->rank;
+  result->shape = gfc_copy_shape (mask->shape, mask->rank);
+
+  if (vector->ts.type == BT_CHARACTER)
+    result->ts.u.cl = vector->ts.u.cl;
+
+  vector_ctor = vector->value.constructor;
+  mask_ctor = mask->value.constructor;
+  field_ctor = field->expr_type == EXPR_ARRAY ? field->value.constructor : NULL;
+
+  while (mask_ctor)
+    {
+      if (mask_ctor->expr->value.logical)
+       {
+         gcc_assert (vector_ctor);
+         e = gfc_copy_expr (vector_ctor->expr);
+         ADVANCE (vector_ctor, 1);
+       }
+      else if (field->expr_type == EXPR_ARRAY)
+       e = gfc_copy_expr (field_ctor->expr);
+      else
+       e = gfc_copy_expr (field);
+
+      gfc_append_constructor (result, e);
+
+      ADVANCE (mask_ctor, 1);
+      ADVANCE (field_ctor, 1);
+    }
+
+  return result;
+}
+
+
+gfc_expr *
 gfc_simplify_verify (gfc_expr *s, gfc_expr *set, gfc_expr *b, gfc_expr *kind)
 {
   gfc_expr *result;
@@ -5431,7 +5883,7 @@ gfc_convert_char_constant (gfc_expr *e, bt type ATTRIBUTE_UNUSED, int kind)
       result->shape = gfc_copy_shape (e->shape, e->rank);
       result->where = e->where;
       result->rank = e->rank;
-      result->ts.cl = e->ts.cl;
+      result->ts.u.cl = e->ts.u.cl;
 
       return result;
     }