OSDN Git Service

2009-07-08 Tobias Burnus <burnus@net-b.de>
[pf3gnuchains/gcc-fork.git] / gcc / fortran / simplify.c
1 /* Simplify intrinsic functions at compile-time.
2    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009
3    Free Software Foundation, Inc.
4    Contributed by Andy Vaught & Katherine Holcomb
5
6 This file is part of GCC.
7
8 GCC is free software; you can redistribute it and/or modify it under
9 the terms of the GNU General Public License as published by the Free
10 Software Foundation; either version 3, or (at your option) any later
11 version.
12
13 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
14 WARRANTY; without even the implied warranty of MERCHANTABILITY or
15 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
16 for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with GCC; see the file COPYING3.  If not see
20 <http://www.gnu.org/licenses/>.  */
21
22 #include "config.h"
23 #include "system.h"
24 #include "flags.h"
25 #include "gfortran.h"
26 #include "arith.h"
27 #include "intrinsic.h"
28 #include "target-memory.h"
29
30 /* Savely advance an array constructor by 'n' elements.
31    Mainly used by simplifiers of transformational intrinsics.  */
32 #define ADVANCE(ctor, n) do { int i; for (i = 0; i < n && ctor; ++i) ctor = ctor->next; } while (0)
33
34 gfc_expr gfc_bad_expr;
35
36
37 /* Note that 'simplification' is not just transforming expressions.
38    For functions that are not simplified at compile time, range
39    checking is done if possible.
40
41    The return convention is that each simplification function returns:
42
43      A new expression node corresponding to the simplified arguments.
44      The original arguments are destroyed by the caller, and must not
45      be a part of the new expression.
46
47      NULL pointer indicating that no simplification was possible and
48      the original expression should remain intact.  If the
49      simplification function sets the type and/or the function name
50      via the pointer gfc_simple_expression, then this type is
51      retained.
52
53      An expression pointer to gfc_bad_expr (a static placeholder)
54      indicating that some error has prevented simplification.  For
55      example, sqrt(-1.0).  The error is generated within the function
56      and should be propagated upwards
57
58    By the time a simplification function gets control, it has been
59    decided that the function call is really supposed to be the
60    intrinsic.  No type checking is strictly necessary, since only
61    valid types will be passed on.  On the other hand, a simplification
62    subroutine may have to look at the type of an argument as part of
63    its processing.
64
65    Array arguments are never passed to these subroutines.
66
67    The functions in this file don't have much comment with them, but
68    everything is reasonably straight-forward.  The Standard, chapter 13
69    is the best comment you'll find for this file anyway.  */
70
71 /* Range checks an expression node.  If all goes well, returns the
72    node, otherwise returns &gfc_bad_expr and frees the node.  */
73
74 static gfc_expr *
75 range_check (gfc_expr *result, const char *name)
76 {
77   if (result == NULL)
78     return &gfc_bad_expr;
79
80   switch (gfc_range_check (result))
81     {
82       case ARITH_OK:
83         return result;
84  
85       case ARITH_OVERFLOW:
86         gfc_error ("Result of %s overflows its kind at %L", name,
87                    &result->where);
88         break;
89
90       case ARITH_UNDERFLOW:
91         gfc_error ("Result of %s underflows its kind at %L", name,
92                    &result->where);
93         break;
94
95       case ARITH_NAN:
96         gfc_error ("Result of %s is NaN at %L", name, &result->where);
97         break;
98
99       default:
100         gfc_error ("Result of %s gives range error for its kind at %L", name,
101                    &result->where);
102         break;
103     }
104
105   gfc_free_expr (result);
106   return &gfc_bad_expr;
107 }
108
109
110 /* A helper function that gets an optional and possibly missing
111    kind parameter.  Returns the kind, -1 if something went wrong.  */
112
113 static int
114 get_kind (bt type, gfc_expr *k, const char *name, int default_kind)
115 {
116   int kind;
117
118   if (k == NULL)
119     return default_kind;
120
121   if (k->expr_type != EXPR_CONSTANT)
122     {
123       gfc_error ("KIND parameter of %s at %L must be an initialization "
124                  "expression", name, &k->where);
125       return -1;
126     }
127
128   if (gfc_extract_int (k, &kind) != NULL
129       || gfc_validate_kind (type, kind, true) < 0)
130     {
131       gfc_error ("Invalid KIND parameter of %s at %L", name, &k->where);
132       return -1;
133     }
134
135   return kind;
136 }
137
138
139 /* Helper function to get an integer constant with a kind number given
140    by an integer constant expression.  */
141 static gfc_expr *
142 int_expr_with_kind (int i, gfc_expr *kind, const char *name)
143 {
144   gfc_expr *res = gfc_int_expr (i);
145   res->ts.kind = get_kind (BT_INTEGER, kind, name, gfc_default_integer_kind); 
146   if (res->ts.kind == -1)
147     return NULL;
148   else
149     return res;
150 }
151
152
153 /* Converts an mpz_t signed variable into an unsigned one, assuming
154    two's complement representations and a binary width of bitsize.
155    The conversion is a no-op unless x is negative; otherwise, it can
156    be accomplished by masking out the high bits.  */
157
158 static void
159 convert_mpz_to_unsigned (mpz_t x, int bitsize)
160 {
161   mpz_t mask;
162
163   if (mpz_sgn (x) < 0)
164     {
165       /* Confirm that no bits above the signed range are unset.  */
166       gcc_assert (mpz_scan0 (x, bitsize-1) == ULONG_MAX);
167
168       mpz_init_set_ui (mask, 1);
169       mpz_mul_2exp (mask, mask, bitsize);
170       mpz_sub_ui (mask, mask, 1);
171
172       mpz_and (x, x, mask);
173
174       mpz_clear (mask);
175     }
176   else
177     {
178       /* Confirm that no bits above the signed range are set.  */
179       gcc_assert (mpz_scan1 (x, bitsize-1) == ULONG_MAX);
180     }
181 }
182
183
184 /* Converts an mpz_t unsigned variable into a signed one, assuming
185    two's complement representations and a binary width of bitsize.
186    If the bitsize-1 bit is set, this is taken as a sign bit and
187    the number is converted to the corresponding negative number.  */
188
189 static void
190 convert_mpz_to_signed (mpz_t x, int bitsize)
191 {
192   mpz_t mask;
193
194   /* Confirm that no bits above the unsigned range are set.  */
195   gcc_assert (mpz_scan1 (x, bitsize) == ULONG_MAX);
196
197   if (mpz_tstbit (x, bitsize - 1) == 1)
198     {
199       mpz_init_set_ui (mask, 1);
200       mpz_mul_2exp (mask, mask, bitsize);
201       mpz_sub_ui (mask, mask, 1);
202
203       /* We negate the number by hand, zeroing the high bits, that is
204          make it the corresponding positive number, and then have it
205          negated by GMP, giving the correct representation of the
206          negative number.  */
207       mpz_com (x, x);
208       mpz_add_ui (x, x, 1);
209       mpz_and (x, x, mask);
210
211       mpz_neg (x, x);
212
213       mpz_clear (mask);
214     }
215 }
216
217 /* Test that the expression is an constant array.  */
218
219 static bool
220 is_constant_array_expr (gfc_expr *e)
221 {
222   gfc_constructor *c;
223
224   if (e == NULL)
225     return true;
226
227   if (e->expr_type != EXPR_ARRAY || !gfc_is_constant_expr (e))
228     return false;
229
230   for (c = e->value.constructor; c; c = c->next)
231     if (c->expr->expr_type != EXPR_CONSTANT)
232       return false;
233
234   return true;
235 }
236
237
238 /* Initialize a transformational result expression with a given value.  */
239
240 static void
241 init_result_expr (gfc_expr *e, int init, gfc_expr *array)
242 {
243   if (e && e->expr_type == EXPR_ARRAY)
244     {
245       gfc_constructor *ctor = e->value.constructor;
246       while (ctor)
247         {
248           init_result_expr (ctor->expr, init, array);
249           ctor = ctor->next;
250         }
251     }
252   else if (e && e->expr_type == EXPR_CONSTANT)
253     {
254       int i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
255       int length;
256       gfc_char_t *string;
257
258       switch (e->ts.type)
259         {
260           case BT_LOGICAL:
261             e->value.logical = (init ? 1 : 0);
262             break;
263
264           case BT_INTEGER:
265             if (init == INT_MIN)
266               mpz_set (e->value.integer, gfc_integer_kinds[i].min_int);
267             else if (init == INT_MAX)
268               mpz_set (e->value.integer, gfc_integer_kinds[i].huge);
269             else
270               mpz_set_si (e->value.integer, init);
271             break;
272
273           case BT_REAL:
274             if (init == INT_MIN)
275               {
276                 mpfr_set (e->value.real, gfc_real_kinds[i].huge, GFC_RND_MODE);
277                 mpfr_neg (e->value.real, e->value.real, GFC_RND_MODE);
278               }
279             else if (init == INT_MAX)
280               mpfr_set (e->value.real, gfc_real_kinds[i].huge, GFC_RND_MODE);
281             else
282               mpfr_set_si (e->value.real, init, GFC_RND_MODE);
283             break;
284
285           case BT_COMPLEX:
286 #ifdef HAVE_mpc
287             mpc_set_si (e->value.complex, init, GFC_MPC_RND_MODE);
288 #else
289             mpfr_set_si (e->value.complex.r, init, GFC_RND_MODE);
290             mpfr_set_si (e->value.complex.i, 0, GFC_RND_MODE);
291 #endif
292             break;
293
294           case BT_CHARACTER:
295             if (init == INT_MIN)
296               {
297                 gfc_expr *len = gfc_simplify_len (array, NULL);
298                 gfc_extract_int (len, &length);
299                 string = gfc_get_wide_string (length + 1);
300                 gfc_wide_memset (string, 0, length);
301               }
302             else if (init == INT_MAX)
303               {
304                 gfc_expr *len = gfc_simplify_len (array, NULL);
305                 gfc_extract_int (len, &length);
306                 string = gfc_get_wide_string (length + 1);
307                 gfc_wide_memset (string, 255, length);
308               }
309             else
310               {
311                 length = 0;
312                 string = gfc_get_wide_string (1);
313               }
314
315             string[length] = '\0';
316             e->value.character.length = length;
317             e->value.character.string = string;
318             break;
319
320           default:
321             gcc_unreachable();
322         }
323     }
324   else
325     gcc_unreachable();
326 }
327
328
329 /* Helper function for gfc_simplify_dot_product() and gfc_simplify_matmul.  */
330
331 static gfc_expr *
332 compute_dot_product (gfc_constructor *ctor_a, int stride_a,
333                      gfc_constructor *ctor_b, int stride_b)
334 {
335   gfc_expr *result;
336   gfc_expr *a = ctor_a->expr, *b = ctor_b->expr;
337
338   gcc_assert (gfc_compare_types (&a->ts, &b->ts));
339
340   result = gfc_constant_result (a->ts.type, a->ts.kind, &a->where);
341   init_result_expr (result, 0, NULL);
342
343   while (ctor_a && ctor_b)
344     {
345       /* Copying of expressions is required as operands are free'd
346          by the gfc_arith routines.  */
347       switch (result->ts.type)
348         {
349           case BT_LOGICAL:
350             result = gfc_or (result,
351                              gfc_and (gfc_copy_expr (ctor_a->expr),
352                                       gfc_copy_expr (ctor_b->expr)));
353             break;
354
355           case BT_INTEGER:
356           case BT_REAL:
357           case BT_COMPLEX:
358             result = gfc_add (result,
359                               gfc_multiply (gfc_copy_expr (ctor_a->expr),
360                                             gfc_copy_expr (ctor_b->expr)));
361             break;
362
363           default:
364             gcc_unreachable();
365         }
366
367       ADVANCE (ctor_a, stride_a);
368       ADVANCE (ctor_b, stride_b);
369     }
370
371   return result;
372 }
373
374
375 /* Build a result expression for transformational intrinsics, 
376    depending on DIM. */
377
378 static gfc_expr *
379 transformational_result (gfc_expr *array, gfc_expr *dim, bt type,
380                          int kind, locus* where)
381 {
382   gfc_expr *result;
383   int i, nelem;
384
385   if (!dim || array->rank == 1)
386     return gfc_constant_result (type, kind, where);
387
388   result = gfc_start_constructor (type, kind, where);
389   result->shape = gfc_copy_shape_excluding (array->shape, array->rank, dim);
390   result->rank = array->rank - 1;
391
392   /* gfc_array_size() would count the number of elements in the constructor,
393      we have not built those yet.  */
394   nelem = 1;
395   for  (i = 0; i < result->rank; ++i)
396     nelem *= mpz_get_ui (result->shape[i]);
397
398   for (i = 0; i < nelem; ++i)
399     {
400       gfc_expr *e = gfc_constant_result (type, kind, where);
401       gfc_append_constructor (result, e);
402     }
403
404   return result;
405 }
406
407
408 typedef gfc_expr* (*transformational_op)(gfc_expr*, gfc_expr*);
409
410 /* Wrapper function, implements 'op1 += 1'. Only called if MASK
411    of COUNT intrinsic is .TRUE..
412
413    Interface and implimentation mimics arith functions as
414    gfc_add, gfc_multiply, etc.  */
415
416 static gfc_expr* gfc_count (gfc_expr *op1, gfc_expr *op2)
417 {
418   gfc_expr *result;
419
420   gcc_assert (op1->ts.type == BT_INTEGER);
421   gcc_assert (op2->ts.type == BT_LOGICAL);
422   gcc_assert (op2->value.logical);
423
424   result = gfc_copy_expr (op1);
425   mpz_add_ui (result->value.integer, result->value.integer, 1);
426
427   gfc_free_expr (op1);
428   gfc_free_expr (op2);
429   return result;
430 }
431
432
433 /* Transforms an ARRAY with operation OP, according to MASK, to a
434    scalar RESULT. E.g. called if
435
436      REAL, PARAMETER :: array(n, m) = ...
437      REAL, PARAMETER :: s = SUM(array)
438
439   where OP == gfc_add().  */
440
441 static gfc_expr *
442 simplify_transformation_to_scalar (gfc_expr *result, gfc_expr *array, gfc_expr *mask,
443                                    transformational_op op)
444 {
445   gfc_expr *a, *m;
446   gfc_constructor *array_ctor, *mask_ctor;
447
448   /* Shortcut for constant .FALSE. MASK.  */
449   if (mask
450       && mask->expr_type == EXPR_CONSTANT
451       && !mask->value.logical)
452     return result;
453
454   array_ctor = array->value.constructor;
455   mask_ctor = NULL;
456   if (mask && mask->expr_type == EXPR_ARRAY)
457     mask_ctor = mask->value.constructor;
458
459   while (array_ctor)
460     {
461       a = array_ctor->expr;
462       array_ctor = array_ctor->next;
463
464       /* A constant MASK equals .TRUE. here and can be ignored.  */
465       if (mask_ctor)
466         {
467           m = mask_ctor->expr;
468           mask_ctor = mask_ctor->next;
469           if (!m->value.logical)
470             continue;
471         }
472
473       result = op (result, gfc_copy_expr (a));
474     }
475
476   return result;
477 }
478
479 /* Transforms an ARRAY with operation OP, according to MASK, to an
480    array RESULT. E.g. called if
481
482      REAL, PARAMETER :: array(n, m) = ...
483      REAL, PARAMETER :: s(n) = PROD(array, DIM=1)
484
485   where OP == gfc_multiply().  */
486
487 static gfc_expr *
488 simplify_transformation_to_array (gfc_expr *result, gfc_expr *array, gfc_expr *dim,
489                                   gfc_expr *mask, transformational_op op)
490 {
491   mpz_t size;
492   int done, i, n, arraysize, resultsize, dim_index, dim_extent, dim_stride;
493   gfc_expr **arrayvec, **resultvec, **base, **src, **dest;
494   gfc_constructor *array_ctor, *mask_ctor, *result_ctor;
495
496   int count[GFC_MAX_DIMENSIONS], extent[GFC_MAX_DIMENSIONS],
497       sstride[GFC_MAX_DIMENSIONS], dstride[GFC_MAX_DIMENSIONS],
498       tmpstride[GFC_MAX_DIMENSIONS];
499
500   /* Shortcut for constant .FALSE. MASK.  */
501   if (mask
502       && mask->expr_type == EXPR_CONSTANT
503       && !mask->value.logical)
504     return result;
505
506   /* Build an indexed table for array element expressions to minimize
507      linked-list traversal. Masked elements are set to NULL.  */
508   gfc_array_size (array, &size);
509   arraysize = mpz_get_ui (size);
510
511   arrayvec = (gfc_expr**) gfc_getmem (sizeof (gfc_expr*) * arraysize);
512
513   array_ctor = array->value.constructor;
514   mask_ctor = NULL;
515   if (mask && mask->expr_type == EXPR_ARRAY)
516     mask_ctor = mask->value.constructor;
517
518   for (i = 0; i < arraysize; ++i)
519     {
520       arrayvec[i] = array_ctor->expr;
521       array_ctor = array_ctor->next;
522
523       if (mask_ctor)
524         {
525           if (!mask_ctor->expr->value.logical)
526             arrayvec[i] = NULL;
527
528           mask_ctor = mask_ctor->next;
529         }
530     }
531
532   /* Same for the result expression.  */
533   gfc_array_size (result, &size);
534   resultsize = mpz_get_ui (size);
535   mpz_clear (size);
536
537   resultvec = (gfc_expr**) gfc_getmem (sizeof (gfc_expr*) * resultsize);
538   result_ctor = result->value.constructor;
539   for (i = 0; i < resultsize; ++i)
540     {
541       resultvec[i] = result_ctor->expr;
542       result_ctor = result_ctor->next;
543     }
544
545   gfc_extract_int (dim, &dim_index);
546   dim_index -= 1;               /* zero-base index */
547   dim_extent = 0;
548   dim_stride = 0;
549
550   for (i = 0, n = 0; i < array->rank; ++i)
551     {
552       count[i] = 0;
553       tmpstride[i] = (i == 0) ? 1 : tmpstride[i-1] * mpz_get_si (array->shape[i-1]);
554       if (i == dim_index)
555         {
556           dim_extent = mpz_get_si (array->shape[i]);
557           dim_stride = tmpstride[i];
558           continue;
559         }
560
561       extent[n] = mpz_get_si (array->shape[i]);
562       sstride[n] = tmpstride[i];
563       dstride[n] = (n == 0) ? 1 : dstride[n-1] * extent[n-1];
564       n += 1;
565     }
566
567   done = false;
568   base = arrayvec;
569   dest = resultvec;
570   while (!done)
571     {
572       for (src = base, n = 0; n < dim_extent; src += dim_stride, ++n)
573         if (*src)
574           *dest = op (*dest, gfc_copy_expr (*src));
575
576       count[0]++;
577       base += sstride[0];
578       dest += dstride[0];
579
580       n = 0;
581       while (!done && count[n] == extent[n])
582         {
583           count[n] = 0;
584           base -= sstride[n] * extent[n];
585           dest -= dstride[n] * extent[n];
586
587           n++;
588           if (n < result->rank)
589             {
590               count [n]++;
591               base += sstride[n];
592               dest += dstride[n];
593             }
594           else
595             done = true;
596        }
597     }
598
599   /* Place updated expression in result constructor.  */
600   result_ctor = result->value.constructor;
601   for (i = 0; i < resultsize; ++i)
602     {
603       result_ctor->expr = resultvec[i];
604       result_ctor = result_ctor->next;
605     }
606
607   gfc_free (arrayvec);
608   gfc_free (resultvec);
609   return result;
610 }
611
612
613
614 /********************** Simplification functions *****************************/
615
616 gfc_expr *
617 gfc_simplify_abs (gfc_expr *e)
618 {
619   gfc_expr *result;
620
621   if (e->expr_type != EXPR_CONSTANT)
622     return NULL;
623
624   switch (e->ts.type)
625     {
626     case BT_INTEGER:
627       result = gfc_constant_result (BT_INTEGER, e->ts.kind, &e->where);
628
629       mpz_abs (result->value.integer, e->value.integer);
630
631       result = range_check (result, "IABS");
632       break;
633
634     case BT_REAL:
635       result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
636
637       mpfr_abs (result->value.real, e->value.real, GFC_RND_MODE);
638
639       result = range_check (result, "ABS");
640       break;
641
642     case BT_COMPLEX:
643       result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
644
645       gfc_set_model_kind (e->ts.kind);
646
647 #ifdef HAVE_mpc
648       mpc_abs (result->value.real, e->value.complex, GFC_RND_MODE);
649 #else
650       mpfr_hypot (result->value.real, e->value.complex.r, 
651                   e->value.complex.i, GFC_RND_MODE);
652 #endif
653       result = range_check (result, "CABS");
654       break;
655
656     default:
657       gfc_internal_error ("gfc_simplify_abs(): Bad type");
658     }
659
660   return result;
661 }
662
663
664 static gfc_expr *
665 simplify_achar_char (gfc_expr *e, gfc_expr *k, const char *name, bool ascii)
666 {
667   gfc_expr *result;
668   int kind;
669   bool too_large = false;
670
671   if (e->expr_type != EXPR_CONSTANT)
672     return NULL;
673
674   kind = get_kind (BT_CHARACTER, k, name, gfc_default_character_kind);
675   if (kind == -1)
676     return &gfc_bad_expr;
677
678   if (mpz_cmp_si (e->value.integer, 0) < 0)
679     {
680       gfc_error ("Argument of %s function at %L is negative", name,
681                  &e->where);
682       return &gfc_bad_expr;
683     }
684
685   if (ascii && gfc_option.warn_surprising
686       && mpz_cmp_si (e->value.integer, 127) > 0)
687     gfc_warning ("Argument of %s function at %L outside of range [0,127]",
688                  name, &e->where);
689
690   if (kind == 1 && mpz_cmp_si (e->value.integer, 255) > 0)
691     too_large = true;
692   else if (kind == 4)
693     {
694       mpz_t t;
695       mpz_init_set_ui (t, 2);
696       mpz_pow_ui (t, t, 32);
697       mpz_sub_ui (t, t, 1);
698       if (mpz_cmp (e->value.integer, t) > 0)
699         too_large = true;
700       mpz_clear (t);
701     }
702
703   if (too_large)
704     {
705       gfc_error ("Argument of %s function at %L is too large for the "
706                  "collating sequence of kind %d", name, &e->where, kind);
707       return &gfc_bad_expr;
708     }
709
710   result = gfc_constant_result (BT_CHARACTER, kind, &e->where);
711   result->value.character.string = gfc_get_wide_string (2);
712   result->value.character.length = 1;
713   result->value.character.string[0] = mpz_get_ui (e->value.integer);
714   result->value.character.string[1] = '\0';     /* For debugger */
715   return result;
716 }
717
718
719
720 /* We use the processor's collating sequence, because all
721    systems that gfortran currently works on are ASCII.  */
722
723 gfc_expr *
724 gfc_simplify_achar (gfc_expr *e, gfc_expr *k)
725 {
726   return simplify_achar_char (e, k, "ACHAR", true);
727 }
728
729
730 gfc_expr *
731 gfc_simplify_acos (gfc_expr *x)
732 {
733   gfc_expr *result;
734
735   if (x->expr_type != EXPR_CONSTANT)
736     return NULL;
737
738   if (mpfr_cmp_si (x->value.real, 1) > 0
739       || mpfr_cmp_si (x->value.real, -1) < 0)
740     {
741       gfc_error ("Argument of ACOS at %L must be between -1 and 1",
742                  &x->where);
743       return &gfc_bad_expr;
744     }
745
746   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
747
748   mpfr_acos (result->value.real, x->value.real, GFC_RND_MODE);
749
750   return range_check (result, "ACOS");
751 }
752
753 gfc_expr *
754 gfc_simplify_acosh (gfc_expr *x)
755 {
756   gfc_expr *result;
757
758   if (x->expr_type != EXPR_CONSTANT)
759     return NULL;
760
761   if (mpfr_cmp_si (x->value.real, 1) < 0)
762     {
763       gfc_error ("Argument of ACOSH at %L must not be less than 1",
764                  &x->where);
765       return &gfc_bad_expr;
766     }
767
768   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
769
770   mpfr_acosh (result->value.real, x->value.real, GFC_RND_MODE);
771
772   return range_check (result, "ACOSH");
773 }
774
775 gfc_expr *
776 gfc_simplify_adjustl (gfc_expr *e)
777 {
778   gfc_expr *result;
779   int count, i, len;
780   gfc_char_t ch;
781
782   if (e->expr_type != EXPR_CONSTANT)
783     return NULL;
784
785   len = e->value.character.length;
786
787   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
788
789   result->value.character.length = len;
790   result->value.character.string = gfc_get_wide_string (len + 1);
791
792   for (count = 0, i = 0; i < len; ++i)
793     {
794       ch = e->value.character.string[i];
795       if (ch != ' ')
796         break;
797       ++count;
798     }
799
800   for (i = 0; i < len - count; ++i)
801     result->value.character.string[i] = e->value.character.string[count + i];
802
803   for (i = len - count; i < len; ++i)
804     result->value.character.string[i] = ' ';
805
806   result->value.character.string[len] = '\0';   /* For debugger */
807
808   return result;
809 }
810
811
812 gfc_expr *
813 gfc_simplify_adjustr (gfc_expr *e)
814 {
815   gfc_expr *result;
816   int count, i, len;
817   gfc_char_t ch;
818
819   if (e->expr_type != EXPR_CONSTANT)
820     return NULL;
821
822   len = e->value.character.length;
823
824   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
825
826   result->value.character.length = len;
827   result->value.character.string = gfc_get_wide_string (len + 1);
828
829   for (count = 0, i = len - 1; i >= 0; --i)
830     {
831       ch = e->value.character.string[i];
832       if (ch != ' ')
833         break;
834       ++count;
835     }
836
837   for (i = 0; i < count; ++i)
838     result->value.character.string[i] = ' ';
839
840   for (i = count; i < len; ++i)
841     result->value.character.string[i] = e->value.character.string[i - count];
842
843   result->value.character.string[len] = '\0';   /* For debugger */
844
845   return result;
846 }
847
848
849 gfc_expr *
850 gfc_simplify_aimag (gfc_expr *e)
851 {
852   gfc_expr *result;
853
854   if (e->expr_type != EXPR_CONSTANT)
855     return NULL;
856
857   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
858   mpfr_set (result->value.real, mpc_imagref (e->value.complex), GFC_RND_MODE);
859
860   return range_check (result, "AIMAG");
861 }
862
863
864 gfc_expr *
865 gfc_simplify_aint (gfc_expr *e, gfc_expr *k)
866 {
867   gfc_expr *rtrunc, *result;
868   int kind;
869
870   kind = get_kind (BT_REAL, k, "AINT", e->ts.kind);
871   if (kind == -1)
872     return &gfc_bad_expr;
873
874   if (e->expr_type != EXPR_CONSTANT)
875     return NULL;
876
877   rtrunc = gfc_copy_expr (e);
878
879   mpfr_trunc (rtrunc->value.real, e->value.real);
880
881   result = gfc_real2real (rtrunc, kind);
882   gfc_free_expr (rtrunc);
883
884   return range_check (result, "AINT");
885 }
886
887
888 gfc_expr *
889 gfc_simplify_all (gfc_expr *mask, gfc_expr *dim)
890 {
891   gfc_expr *result;
892
893   if (!is_constant_array_expr (mask)
894       || !gfc_is_constant_expr (dim))
895     return NULL;
896
897   result = transformational_result (mask, dim, mask->ts.type,
898                                     mask->ts.kind, &mask->where);
899   init_result_expr (result, true, NULL);
900
901   return !dim || mask->rank == 1 ?
902     simplify_transformation_to_scalar (result, mask, NULL, gfc_and) :
903     simplify_transformation_to_array (result, mask, dim, NULL, gfc_and);
904 }
905
906
907 gfc_expr *
908 gfc_simplify_dint (gfc_expr *e)
909 {
910   gfc_expr *rtrunc, *result;
911
912   if (e->expr_type != EXPR_CONSTANT)
913     return NULL;
914
915   rtrunc = gfc_copy_expr (e);
916
917   mpfr_trunc (rtrunc->value.real, e->value.real);
918
919   result = gfc_real2real (rtrunc, gfc_default_double_kind);
920   gfc_free_expr (rtrunc);
921
922   return range_check (result, "DINT");
923 }
924
925
926 gfc_expr *
927 gfc_simplify_anint (gfc_expr *e, gfc_expr *k)
928 {
929   gfc_expr *result;
930   int kind;
931
932   kind = get_kind (BT_REAL, k, "ANINT", e->ts.kind);
933   if (kind == -1)
934     return &gfc_bad_expr;
935
936   if (e->expr_type != EXPR_CONSTANT)
937     return NULL;
938
939   result = gfc_constant_result (e->ts.type, kind, &e->where);
940
941   mpfr_round (result->value.real, e->value.real);
942
943   return range_check (result, "ANINT");
944 }
945
946
947 gfc_expr *
948 gfc_simplify_and (gfc_expr *x, gfc_expr *y)
949 {
950   gfc_expr *result;
951   int kind;
952
953   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
954     return NULL;
955
956   kind = x->ts.kind > y->ts.kind ? x->ts.kind : y->ts.kind;
957   if (x->ts.type == BT_INTEGER)
958     {
959       result = gfc_constant_result (BT_INTEGER, kind, &x->where);
960       mpz_and (result->value.integer, x->value.integer, y->value.integer);
961       return range_check (result, "AND");
962     }
963   else /* BT_LOGICAL */
964     {
965       result = gfc_constant_result (BT_LOGICAL, kind, &x->where);
966       result->value.logical = x->value.logical && y->value.logical;
967       return result;
968     }
969 }
970
971
972 gfc_expr *
973 gfc_simplify_any (gfc_expr *mask, gfc_expr *dim)
974 {
975   gfc_expr *result;
976
977   if (!is_constant_array_expr (mask)
978       || !gfc_is_constant_expr (dim))
979     return NULL;
980
981   result = transformational_result (mask, dim, mask->ts.type,
982                                     mask->ts.kind, &mask->where);
983   init_result_expr (result, false, NULL);
984
985   return !dim || mask->rank == 1 ?
986     simplify_transformation_to_scalar (result, mask, NULL, gfc_or) :
987     simplify_transformation_to_array (result, mask, dim, NULL, gfc_or);
988 }
989
990
991 gfc_expr *
992 gfc_simplify_dnint (gfc_expr *e)
993 {
994   gfc_expr *result;
995
996   if (e->expr_type != EXPR_CONSTANT)
997     return NULL;
998
999   result = gfc_constant_result (BT_REAL, gfc_default_double_kind, &e->where);
1000
1001   mpfr_round (result->value.real, e->value.real);
1002
1003   return range_check (result, "DNINT");
1004 }
1005
1006
1007 gfc_expr *
1008 gfc_simplify_asin (gfc_expr *x)
1009 {
1010   gfc_expr *result;
1011
1012   if (x->expr_type != EXPR_CONSTANT)
1013     return NULL;
1014
1015   if (mpfr_cmp_si (x->value.real, 1) > 0
1016       || mpfr_cmp_si (x->value.real, -1) < 0)
1017     {
1018       gfc_error ("Argument of ASIN at %L must be between -1 and 1",
1019                  &x->where);
1020       return &gfc_bad_expr;
1021     }
1022
1023   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1024
1025   mpfr_asin (result->value.real, x->value.real, GFC_RND_MODE);
1026
1027   return range_check (result, "ASIN");
1028 }
1029
1030
1031 gfc_expr *
1032 gfc_simplify_asinh (gfc_expr *x)
1033 {
1034   gfc_expr *result;
1035
1036   if (x->expr_type != EXPR_CONSTANT)
1037     return NULL;
1038
1039   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1040
1041   mpfr_asinh (result->value.real, x->value.real, GFC_RND_MODE);
1042
1043   return range_check (result, "ASINH");
1044 }
1045
1046
1047 gfc_expr *
1048 gfc_simplify_atan (gfc_expr *x)
1049 {
1050   gfc_expr *result;
1051
1052   if (x->expr_type != EXPR_CONSTANT)
1053     return NULL;
1054     
1055   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1056
1057   mpfr_atan (result->value.real, x->value.real, GFC_RND_MODE);
1058
1059   return range_check (result, "ATAN");
1060 }
1061
1062
1063 gfc_expr *
1064 gfc_simplify_atanh (gfc_expr *x)
1065 {
1066   gfc_expr *result;
1067
1068   if (x->expr_type != EXPR_CONSTANT)
1069     return NULL;
1070
1071   if (mpfr_cmp_si (x->value.real, 1) >= 0
1072       || mpfr_cmp_si (x->value.real, -1) <= 0)
1073     {
1074       gfc_error ("Argument of ATANH at %L must be inside the range -1 to 1",
1075                  &x->where);
1076       return &gfc_bad_expr;
1077     }
1078
1079   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1080
1081   mpfr_atanh (result->value.real, x->value.real, GFC_RND_MODE);
1082
1083   return range_check (result, "ATANH");
1084 }
1085
1086
1087 gfc_expr *
1088 gfc_simplify_atan2 (gfc_expr *y, gfc_expr *x)
1089 {
1090   gfc_expr *result;
1091
1092   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1093     return NULL;
1094
1095   if (mpfr_sgn (y->value.real) == 0 && mpfr_sgn (x->value.real) == 0)
1096     {
1097       gfc_error ("If first argument of ATAN2 %L is zero, then the "
1098                  "second argument must not be zero", &x->where);
1099       return &gfc_bad_expr;
1100     }
1101
1102   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1103
1104   mpfr_atan2 (result->value.real, y->value.real, x->value.real, GFC_RND_MODE);
1105
1106   return range_check (result, "ATAN2");
1107 }
1108
1109
1110 gfc_expr *
1111 gfc_simplify_bessel_j0 (gfc_expr *x ATTRIBUTE_UNUSED)
1112 {
1113   gfc_expr *result;
1114
1115   if (x->expr_type != EXPR_CONSTANT)
1116     return NULL;
1117
1118   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1119   mpfr_j0 (result->value.real, x->value.real, GFC_RND_MODE);
1120
1121   return range_check (result, "BESSEL_J0");
1122 }
1123
1124
1125 gfc_expr *
1126 gfc_simplify_bessel_j1 (gfc_expr *x ATTRIBUTE_UNUSED)
1127 {
1128   gfc_expr *result;
1129
1130   if (x->expr_type != EXPR_CONSTANT)
1131     return NULL;
1132
1133   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1134   mpfr_j1 (result->value.real, x->value.real, GFC_RND_MODE);
1135
1136   return range_check (result, "BESSEL_J1");
1137 }
1138
1139
1140 gfc_expr *
1141 gfc_simplify_bessel_jn (gfc_expr *order ATTRIBUTE_UNUSED,
1142                         gfc_expr *x ATTRIBUTE_UNUSED)
1143 {
1144   gfc_expr *result;
1145   long n;
1146
1147   if (x->expr_type != EXPR_CONSTANT || order->expr_type != EXPR_CONSTANT)
1148     return NULL;
1149
1150   n = mpz_get_si (order->value.integer);
1151   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1152   mpfr_jn (result->value.real, n, x->value.real, GFC_RND_MODE);
1153
1154   return range_check (result, "BESSEL_JN");
1155 }
1156
1157
1158 gfc_expr *
1159 gfc_simplify_bessel_y0 (gfc_expr *x ATTRIBUTE_UNUSED)
1160 {
1161   gfc_expr *result;
1162
1163   if (x->expr_type != EXPR_CONSTANT)
1164     return NULL;
1165
1166   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1167   mpfr_y0 (result->value.real, x->value.real, GFC_RND_MODE);
1168
1169   return range_check (result, "BESSEL_Y0");
1170 }
1171
1172
1173 gfc_expr *
1174 gfc_simplify_bessel_y1 (gfc_expr *x ATTRIBUTE_UNUSED)
1175 {
1176   gfc_expr *result;
1177
1178   if (x->expr_type != EXPR_CONSTANT)
1179     return NULL;
1180
1181   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1182   mpfr_y1 (result->value.real, x->value.real, GFC_RND_MODE);
1183
1184   return range_check (result, "BESSEL_Y1");
1185 }
1186
1187
1188 gfc_expr *
1189 gfc_simplify_bessel_yn (gfc_expr *order ATTRIBUTE_UNUSED,
1190                         gfc_expr *x ATTRIBUTE_UNUSED)
1191 {
1192   gfc_expr *result;
1193   long n;
1194
1195   if (x->expr_type != EXPR_CONSTANT || order->expr_type != EXPR_CONSTANT)
1196     return NULL;
1197
1198   n = mpz_get_si (order->value.integer);
1199   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1200   mpfr_yn (result->value.real, n, x->value.real, GFC_RND_MODE);
1201
1202   return range_check (result, "BESSEL_YN");
1203 }
1204
1205
1206 gfc_expr *
1207 gfc_simplify_bit_size (gfc_expr *e)
1208 {
1209   gfc_expr *result;
1210   int i;
1211
1212   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
1213   result = gfc_constant_result (BT_INTEGER, e->ts.kind, &e->where);
1214   mpz_set_ui (result->value.integer, gfc_integer_kinds[i].bit_size);
1215
1216   return result;
1217 }
1218
1219
1220 gfc_expr *
1221 gfc_simplify_btest (gfc_expr *e, gfc_expr *bit)
1222 {
1223   int b;
1224
1225   if (e->expr_type != EXPR_CONSTANT || bit->expr_type != EXPR_CONSTANT)
1226     return NULL;
1227
1228   if (gfc_extract_int (bit, &b) != NULL || b < 0)
1229     return gfc_logical_expr (0, &e->where);
1230
1231   return gfc_logical_expr (mpz_tstbit (e->value.integer, b), &e->where);
1232 }
1233
1234
1235 gfc_expr *
1236 gfc_simplify_ceiling (gfc_expr *e, gfc_expr *k)
1237 {
1238   gfc_expr *ceil, *result;
1239   int kind;
1240
1241   kind = get_kind (BT_INTEGER, k, "CEILING", gfc_default_integer_kind);
1242   if (kind == -1)
1243     return &gfc_bad_expr;
1244
1245   if (e->expr_type != EXPR_CONSTANT)
1246     return NULL;
1247
1248   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
1249
1250   ceil = gfc_copy_expr (e);
1251
1252   mpfr_ceil (ceil->value.real, e->value.real);
1253   gfc_mpfr_to_mpz (result->value.integer, ceil->value.real, &e->where);
1254
1255   gfc_free_expr (ceil);
1256
1257   return range_check (result, "CEILING");
1258 }
1259
1260
1261 gfc_expr *
1262 gfc_simplify_char (gfc_expr *e, gfc_expr *k)
1263 {
1264   return simplify_achar_char (e, k, "CHAR", false);
1265 }
1266
1267
1268 /* Common subroutine for simplifying CMPLX and DCMPLX.  */
1269
1270 static gfc_expr *
1271 simplify_cmplx (const char *name, gfc_expr *x, gfc_expr *y, int kind)
1272 {
1273   gfc_expr *result;
1274
1275   result = gfc_constant_result (BT_COMPLEX, kind, &x->where);
1276
1277 #ifndef HAVE_mpc
1278   mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
1279 #endif
1280
1281   switch (x->ts.type)
1282     {
1283     case BT_INTEGER:
1284       if (!x->is_boz)
1285 #ifdef HAVE_mpc
1286         mpc_set_z (result->value.complex, x->value.integer, GFC_MPC_RND_MODE);
1287 #else
1288         mpfr_set_z (result->value.complex.r, x->value.integer, GFC_RND_MODE);
1289 #endif
1290       break;
1291
1292     case BT_REAL:
1293 #ifdef HAVE_mpc
1294       mpc_set_fr (result->value.complex, x->value.real, GFC_RND_MODE);
1295 #else
1296       mpfr_set (result->value.complex.r, x->value.real, GFC_RND_MODE);
1297 #endif
1298       break;
1299
1300     case BT_COMPLEX:
1301 #ifdef HAVE_mpc
1302       mpc_set (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
1303 #else
1304       mpfr_set (result->value.complex.r, x->value.complex.r, GFC_RND_MODE);
1305       mpfr_set (result->value.complex.i, x->value.complex.i, GFC_RND_MODE);
1306 #endif
1307       break;
1308
1309     default:
1310       gfc_internal_error ("gfc_simplify_dcmplx(): Bad type (x)");
1311     }
1312
1313   if (y != NULL)
1314     {
1315       switch (y->ts.type)
1316         {
1317         case BT_INTEGER:
1318           if (!y->is_boz)
1319             mpfr_set_z (mpc_imagref (result->value.complex),
1320                         y->value.integer, GFC_RND_MODE);
1321           break;
1322
1323         case BT_REAL:
1324           mpfr_set (mpc_imagref (result->value.complex),
1325                     y->value.real, GFC_RND_MODE);
1326           break;
1327
1328         default:
1329           gfc_internal_error ("gfc_simplify_dcmplx(): Bad type (y)");
1330         }
1331     }
1332
1333   /* Handle BOZ.  */
1334   if (x->is_boz)
1335     {
1336       gfc_typespec ts;
1337       gfc_clear_ts (&ts);
1338       ts.kind = result->ts.kind;
1339       ts.type = BT_REAL;
1340       if (!gfc_convert_boz (x, &ts))
1341         return &gfc_bad_expr;
1342       mpfr_set (mpc_realref (result->value.complex),
1343                 x->value.real, GFC_RND_MODE);
1344     }
1345
1346   if (y && y->is_boz)
1347     {
1348       gfc_typespec ts;
1349       gfc_clear_ts (&ts);
1350       ts.kind = result->ts.kind;
1351       ts.type = BT_REAL;
1352       if (!gfc_convert_boz (y, &ts))
1353         return &gfc_bad_expr;
1354       mpfr_set (mpc_imagref (result->value.complex),
1355                 y->value.real, GFC_RND_MODE);
1356     }
1357
1358   return range_check (result, name);
1359 }
1360
1361
1362 /* Function called when we won't simplify an expression like CMPLX (or
1363    COMPLEX or DCMPLX) but still want to convert BOZ arguments.  */
1364
1365 static gfc_expr *
1366 only_convert_cmplx_boz (gfc_expr *x, gfc_expr *y, int kind)
1367 {
1368   gfc_typespec ts;
1369   gfc_clear_ts (&ts);
1370   ts.type = BT_REAL;
1371   ts.kind = kind;
1372
1373   if (x->is_boz && !gfc_convert_boz (x, &ts))
1374     return &gfc_bad_expr;
1375
1376   if (y && y->is_boz && !gfc_convert_boz (y, &ts))
1377     return &gfc_bad_expr;
1378
1379   return NULL;
1380 }
1381
1382
1383 gfc_expr *
1384 gfc_simplify_cmplx (gfc_expr *x, gfc_expr *y, gfc_expr *k)
1385 {
1386   int kind;
1387
1388   kind = get_kind (BT_REAL, k, "CMPLX", gfc_default_real_kind);
1389   if (kind == -1)
1390     return &gfc_bad_expr;
1391
1392   if (x->expr_type != EXPR_CONSTANT
1393       || (y != NULL && y->expr_type != EXPR_CONSTANT))
1394     return only_convert_cmplx_boz (x, y, kind);
1395
1396   return simplify_cmplx ("CMPLX", x, y, kind);
1397 }
1398
1399
1400 gfc_expr *
1401 gfc_simplify_complex (gfc_expr *x, gfc_expr *y)
1402 {
1403   int kind;
1404
1405   if (x->ts.type == BT_INTEGER)
1406     {
1407       if (y->ts.type == BT_INTEGER)
1408         kind = gfc_default_real_kind;
1409       else
1410         kind = y->ts.kind;
1411     }
1412   else
1413     {
1414       if (y->ts.type == BT_REAL)
1415         kind = (x->ts.kind > y->ts.kind) ? x->ts.kind : y->ts.kind;
1416       else
1417         kind = x->ts.kind;
1418     }
1419
1420   if (x->expr_type != EXPR_CONSTANT
1421       || (y != NULL && y->expr_type != EXPR_CONSTANT))
1422     return only_convert_cmplx_boz (x, y, kind);
1423
1424   return simplify_cmplx ("COMPLEX", x, y, kind);
1425 }
1426
1427
1428 gfc_expr *
1429 gfc_simplify_conjg (gfc_expr *e)
1430 {
1431   gfc_expr *result;
1432
1433   if (e->expr_type != EXPR_CONSTANT)
1434     return NULL;
1435
1436   result = gfc_copy_expr (e);
1437 #ifdef HAVE_mpc
1438   mpc_conj (result->value.complex, result->value.complex, GFC_MPC_RND_MODE);
1439 #else
1440   mpfr_neg (result->value.complex.i, result->value.complex.i, GFC_RND_MODE);
1441 #endif
1442
1443   return range_check (result, "CONJG");
1444 }
1445
1446
1447 gfc_expr *
1448 gfc_simplify_cos (gfc_expr *x)
1449 {
1450   gfc_expr *result;
1451
1452   if (x->expr_type != EXPR_CONSTANT)
1453     return NULL;
1454
1455   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1456
1457   switch (x->ts.type)
1458     {
1459     case BT_REAL:
1460       mpfr_cos (result->value.real, x->value.real, GFC_RND_MODE);
1461       break;
1462     case BT_COMPLEX:
1463       gfc_set_model_kind (x->ts.kind);
1464 #ifdef HAVE_mpc
1465       mpc_cos (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
1466 #else
1467     {
1468       mpfr_t xp, xq;
1469       mpfr_init (xp);
1470       mpfr_init (xq);
1471
1472       mpfr_cos  (xp, x->value.complex.r, GFC_RND_MODE);
1473       mpfr_cosh (xq, x->value.complex.i, GFC_RND_MODE);
1474       mpfr_mul (result->value.complex.r, xp, xq, GFC_RND_MODE);
1475
1476       mpfr_sin  (xp, x->value.complex.r, GFC_RND_MODE);
1477       mpfr_sinh (xq, x->value.complex.i, GFC_RND_MODE);
1478       mpfr_mul (xp, xp, xq, GFC_RND_MODE);
1479       mpfr_neg (result->value.complex.i, xp, GFC_RND_MODE );
1480
1481       mpfr_clears (xp, xq, NULL);
1482     }
1483 #endif
1484       break;
1485     default:
1486       gfc_internal_error ("in gfc_simplify_cos(): Bad type");
1487     }
1488
1489   return range_check (result, "COS");
1490
1491 }
1492
1493
1494 gfc_expr *
1495 gfc_simplify_cosh (gfc_expr *x)
1496 {
1497   gfc_expr *result;
1498
1499   if (x->expr_type != EXPR_CONSTANT)
1500     return NULL;
1501
1502   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1503
1504   mpfr_cosh (result->value.real, x->value.real, GFC_RND_MODE);
1505
1506   return range_check (result, "COSH");
1507 }
1508
1509
1510 gfc_expr *
1511 gfc_simplify_count (gfc_expr *mask, gfc_expr *dim, gfc_expr *kind)
1512 {
1513   gfc_expr *result;
1514
1515   if (!is_constant_array_expr (mask)
1516       || !gfc_is_constant_expr (dim)
1517       || !gfc_is_constant_expr (kind))
1518     return NULL;
1519
1520   result = transformational_result (mask, dim,
1521                                     BT_INTEGER,
1522                                     get_kind (BT_INTEGER, kind, "COUNT",
1523                                               gfc_default_integer_kind),
1524                                     &mask->where);
1525
1526   init_result_expr (result, 0, NULL);
1527
1528   /* Passing MASK twice, once as data array, once as mask.
1529      Whenever gfc_count is called, '1' is added to the result.  */
1530   return !dim || mask->rank == 1 ?
1531     simplify_transformation_to_scalar (result, mask, mask, gfc_count) :
1532     simplify_transformation_to_array (result, mask, dim, mask, gfc_count);
1533 }
1534
1535
1536 gfc_expr *
1537 gfc_simplify_dcmplx (gfc_expr *x, gfc_expr *y)
1538 {
1539
1540   if (x->expr_type != EXPR_CONSTANT
1541       || (y != NULL && y->expr_type != EXPR_CONSTANT))
1542     return only_convert_cmplx_boz (x, y, gfc_default_double_kind);
1543
1544   return simplify_cmplx ("DCMPLX", x, y, gfc_default_double_kind);
1545 }
1546
1547
1548 gfc_expr *
1549 gfc_simplify_dble (gfc_expr *e)
1550 {
1551   gfc_expr *result = NULL;
1552
1553   if (e->expr_type != EXPR_CONSTANT)
1554     return NULL;
1555
1556   switch (e->ts.type)
1557     {
1558     case BT_INTEGER:
1559       if (!e->is_boz)
1560         result = gfc_int2real (e, gfc_default_double_kind);
1561       break;
1562
1563     case BT_REAL:
1564       result = gfc_real2real (e, gfc_default_double_kind);
1565       break;
1566
1567     case BT_COMPLEX:
1568       result = gfc_complex2real (e, gfc_default_double_kind);
1569       break;
1570
1571     default:
1572       gfc_internal_error ("gfc_simplify_dble(): bad type at %L", &e->where);
1573     }
1574
1575   if (e->ts.type == BT_INTEGER && e->is_boz)
1576     {
1577       gfc_typespec ts;
1578       gfc_clear_ts (&ts);
1579       ts.type = BT_REAL;
1580       ts.kind = gfc_default_double_kind;
1581       result = gfc_copy_expr (e);
1582       if (!gfc_convert_boz (result, &ts))
1583         {
1584           gfc_free_expr (result);
1585           return &gfc_bad_expr;
1586         }
1587     }
1588
1589   return range_check (result, "DBLE");
1590 }
1591
1592
1593 gfc_expr *
1594 gfc_simplify_digits (gfc_expr *x)
1595 {
1596   int i, digits;
1597
1598   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
1599   switch (x->ts.type)
1600     {
1601     case BT_INTEGER:
1602       digits = gfc_integer_kinds[i].digits;
1603       break;
1604
1605     case BT_REAL:
1606     case BT_COMPLEX:
1607       digits = gfc_real_kinds[i].digits;
1608       break;
1609
1610     default:
1611       gcc_unreachable ();
1612     }
1613
1614   return gfc_int_expr (digits);
1615 }
1616
1617
1618 gfc_expr *
1619 gfc_simplify_dim (gfc_expr *x, gfc_expr *y)
1620 {
1621   gfc_expr *result;
1622   int kind;
1623
1624   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1625     return NULL;
1626
1627   kind = x->ts.kind > y->ts.kind ? x->ts.kind : y->ts.kind;
1628   result = gfc_constant_result (x->ts.type, kind, &x->where);
1629
1630   switch (x->ts.type)
1631     {
1632     case BT_INTEGER:
1633       if (mpz_cmp (x->value.integer, y->value.integer) > 0)
1634         mpz_sub (result->value.integer, x->value.integer, y->value.integer);
1635       else
1636         mpz_set_ui (result->value.integer, 0);
1637
1638       break;
1639
1640     case BT_REAL:
1641       if (mpfr_cmp (x->value.real, y->value.real) > 0)
1642         mpfr_sub (result->value.real, x->value.real, y->value.real,
1643                   GFC_RND_MODE);
1644       else
1645         mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
1646
1647       break;
1648
1649     default:
1650       gfc_internal_error ("gfc_simplify_dim(): Bad type");
1651     }
1652
1653   return range_check (result, "DIM");
1654 }
1655
1656
1657 gfc_expr*
1658 gfc_simplify_dot_product (gfc_expr *vector_a, gfc_expr *vector_b)
1659 {
1660   gfc_expr *result;
1661
1662   if (!is_constant_array_expr (vector_a)
1663       || !is_constant_array_expr (vector_b))
1664     return NULL;
1665
1666   gcc_assert (vector_a->rank == 1);
1667   gcc_assert (vector_b->rank == 1);
1668   gcc_assert (gfc_compare_types (&vector_a->ts, &vector_b->ts));
1669
1670   if (vector_a->value.constructor && vector_b->value.constructor)
1671     return compute_dot_product (vector_a->value.constructor, 1,
1672                                 vector_b->value.constructor, 1);
1673
1674   /* Zero sized array ...  */
1675   result = gfc_constant_result (vector_a->ts.type,
1676                                 vector_a->ts.kind,
1677                                 &vector_a->where);
1678   init_result_expr (result, 0, NULL);
1679   return result;
1680 }
1681
1682
1683 gfc_expr *
1684 gfc_simplify_dprod (gfc_expr *x, gfc_expr *y)
1685 {
1686   gfc_expr *a1, *a2, *result;
1687
1688   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1689     return NULL;
1690
1691   result = gfc_constant_result (BT_REAL, gfc_default_double_kind, &x->where);
1692
1693   a1 = gfc_real2real (x, gfc_default_double_kind);
1694   a2 = gfc_real2real (y, gfc_default_double_kind);
1695
1696   mpfr_mul (result->value.real, a1->value.real, a2->value.real, GFC_RND_MODE);
1697
1698   gfc_free_expr (a1);
1699   gfc_free_expr (a2);
1700
1701   return range_check (result, "DPROD");
1702 }
1703
1704
1705 gfc_expr *
1706 gfc_simplify_erf (gfc_expr *x)
1707 {
1708   gfc_expr *result;
1709
1710   if (x->expr_type != EXPR_CONSTANT)
1711     return NULL;
1712
1713   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1714
1715   mpfr_erf (result->value.real, x->value.real, GFC_RND_MODE);
1716
1717   return range_check (result, "ERF");
1718 }
1719
1720
1721 gfc_expr *
1722 gfc_simplify_erfc (gfc_expr *x)
1723 {
1724   gfc_expr *result;
1725
1726   if (x->expr_type != EXPR_CONSTANT)
1727     return NULL;
1728
1729   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1730
1731   mpfr_erfc (result->value.real, x->value.real, GFC_RND_MODE);
1732
1733   return range_check (result, "ERFC");
1734 }
1735
1736
1737 /* Helper functions to simplify ERFC_SCALED(x) = ERFC(x) * EXP(X**2).  */
1738
1739 #define MAX_ITER 200
1740 #define ARG_LIMIT 12
1741
1742 /* Calculate ERFC_SCALED directly by its definition:
1743
1744      ERFC_SCALED(x) = ERFC(x) * EXP(X**2)
1745
1746    using a large precision for intermediate results.  This is used for all
1747    but large values of the argument.  */
1748 static void
1749 fullprec_erfc_scaled (mpfr_t res, mpfr_t arg)
1750 {
1751   mp_prec_t prec;
1752   mpfr_t a, b;
1753
1754   prec = mpfr_get_default_prec ();
1755   mpfr_set_default_prec (10 * prec);
1756
1757   mpfr_init (a);
1758   mpfr_init (b);
1759
1760   mpfr_set (a, arg, GFC_RND_MODE);
1761   mpfr_sqr (b, a, GFC_RND_MODE);
1762   mpfr_exp (b, b, GFC_RND_MODE);
1763   mpfr_erfc (a, a, GFC_RND_MODE);
1764   mpfr_mul (a, a, b, GFC_RND_MODE);
1765
1766   mpfr_set (res, a, GFC_RND_MODE);
1767   mpfr_set_default_prec (prec);
1768
1769   mpfr_clear (a);
1770   mpfr_clear (b);
1771 }
1772
1773 /* Calculate ERFC_SCALED using a power series expansion in 1/arg:
1774
1775     ERFC_SCALED(x) = 1 / (x * sqrt(pi))
1776                      * (1 + Sum_n (-1)**n * (1 * 3 * 5 * ... * (2n-1))
1777                                           / (2 * x**2)**n)
1778
1779   This is used for large values of the argument.  Intermediate calculations
1780   are performed with twice the precision.  We don't do a fixed number of
1781   iterations of the sum, but stop when it has converged to the required
1782   precision.  */
1783 static void
1784 asympt_erfc_scaled (mpfr_t res, mpfr_t arg)
1785 {
1786   mpfr_t sum, x, u, v, w, oldsum, sumtrunc;
1787   mpz_t num;
1788   mp_prec_t prec;
1789   unsigned i;
1790
1791   prec = mpfr_get_default_prec ();
1792   mpfr_set_default_prec (2 * prec);
1793
1794   mpfr_init (sum);
1795   mpfr_init (x);
1796   mpfr_init (u);
1797   mpfr_init (v);
1798   mpfr_init (w);
1799   mpz_init (num);
1800
1801   mpfr_init (oldsum);
1802   mpfr_init (sumtrunc);
1803   mpfr_set_prec (oldsum, prec);
1804   mpfr_set_prec (sumtrunc, prec);
1805
1806   mpfr_set (x, arg, GFC_RND_MODE);
1807   mpfr_set_ui (sum, 1, GFC_RND_MODE);
1808   mpz_set_ui (num, 1);
1809
1810   mpfr_set (u, x, GFC_RND_MODE);
1811   mpfr_sqr (u, u, GFC_RND_MODE);
1812   mpfr_mul_ui (u, u, 2, GFC_RND_MODE);
1813   mpfr_pow_si (u, u, -1, GFC_RND_MODE);
1814
1815   for (i = 1; i < MAX_ITER; i++)
1816   {
1817     mpfr_set (oldsum, sum, GFC_RND_MODE);
1818
1819     mpz_mul_ui (num, num, 2 * i - 1);
1820     mpz_neg (num, num);
1821
1822     mpfr_set (w, u, GFC_RND_MODE);
1823     mpfr_pow_ui (w, w, i, GFC_RND_MODE);
1824
1825     mpfr_set_z (v, num, GFC_RND_MODE);
1826     mpfr_mul (v, v, w, GFC_RND_MODE);
1827
1828     mpfr_add (sum, sum, v, GFC_RND_MODE);
1829
1830     mpfr_set (sumtrunc, sum, GFC_RND_MODE);
1831     if (mpfr_cmp (sumtrunc, oldsum) == 0)
1832       break;
1833   }
1834
1835   /* We should have converged by now; otherwise, ARG_LIMIT is probably
1836      set too low.  */
1837   gcc_assert (i < MAX_ITER);
1838
1839   /* Divide by x * sqrt(Pi).  */
1840   mpfr_const_pi (u, GFC_RND_MODE);
1841   mpfr_sqrt (u, u, GFC_RND_MODE);
1842   mpfr_mul (u, u, x, GFC_RND_MODE);
1843   mpfr_div (sum, sum, u, GFC_RND_MODE);
1844
1845   mpfr_set (res, sum, GFC_RND_MODE);
1846   mpfr_set_default_prec (prec);
1847
1848   mpfr_clears (sum, x, u, v, w, oldsum, sumtrunc, NULL);
1849   mpz_clear (num);
1850 }
1851
1852
1853 gfc_expr *
1854 gfc_simplify_erfc_scaled (gfc_expr *x)
1855 {
1856   gfc_expr *result;
1857
1858   if (x->expr_type != EXPR_CONSTANT)
1859     return NULL;
1860
1861   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1862   if (mpfr_cmp_d (x->value.real, ARG_LIMIT) >= 0)
1863     asympt_erfc_scaled (result->value.real, x->value.real);
1864   else
1865     fullprec_erfc_scaled (result->value.real, x->value.real);
1866
1867   return range_check (result, "ERFC_SCALED");
1868 }
1869
1870 #undef MAX_ITER
1871 #undef ARG_LIMIT
1872
1873
1874 gfc_expr *
1875 gfc_simplify_epsilon (gfc_expr *e)
1876 {
1877   gfc_expr *result;
1878   int i;
1879
1880   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
1881
1882   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
1883
1884   mpfr_set (result->value.real, gfc_real_kinds[i].epsilon, GFC_RND_MODE);
1885
1886   return range_check (result, "EPSILON");
1887 }
1888
1889
1890 gfc_expr *
1891 gfc_simplify_exp (gfc_expr *x)
1892 {
1893   gfc_expr *result;
1894
1895   if (x->expr_type != EXPR_CONSTANT)
1896     return NULL;
1897
1898   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1899
1900   switch (x->ts.type)
1901     {
1902     case BT_REAL:
1903       mpfr_exp (result->value.real, x->value.real, GFC_RND_MODE);
1904       break;
1905
1906     case BT_COMPLEX:
1907       gfc_set_model_kind (x->ts.kind);
1908 #ifdef HAVE_mpc
1909       mpc_exp (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
1910 #else
1911     {
1912       mpfr_t xp, xq;
1913       mpfr_init (xp);
1914       mpfr_init (xq);
1915       mpfr_exp (xq, x->value.complex.r, GFC_RND_MODE);
1916       mpfr_cos (xp, x->value.complex.i, GFC_RND_MODE);
1917       mpfr_mul (result->value.complex.r, xq, xp, GFC_RND_MODE);
1918       mpfr_sin (xp, x->value.complex.i, GFC_RND_MODE);
1919       mpfr_mul (result->value.complex.i, xq, xp, GFC_RND_MODE);
1920       mpfr_clears (xp, xq, NULL);
1921     }
1922 #endif
1923       break;
1924
1925     default:
1926       gfc_internal_error ("in gfc_simplify_exp(): Bad type");
1927     }
1928
1929   return range_check (result, "EXP");
1930 }
1931
1932 gfc_expr *
1933 gfc_simplify_exponent (gfc_expr *x)
1934 {
1935   int i;
1936   gfc_expr *result;
1937
1938   if (x->expr_type != EXPR_CONSTANT)
1939     return NULL;
1940
1941   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1942                                 &x->where);
1943
1944   gfc_set_model (x->value.real);
1945
1946   if (mpfr_sgn (x->value.real) == 0)
1947     {
1948       mpz_set_ui (result->value.integer, 0);
1949       return result;
1950     }
1951
1952   i = (int) mpfr_get_exp (x->value.real);
1953   mpz_set_si (result->value.integer, i);
1954
1955   return range_check (result, "EXPONENT");
1956 }
1957
1958
1959 gfc_expr *
1960 gfc_simplify_float (gfc_expr *a)
1961 {
1962   gfc_expr *result;
1963
1964   if (a->expr_type != EXPR_CONSTANT)
1965     return NULL;
1966
1967   if (a->is_boz)
1968     {
1969       gfc_typespec ts;
1970       gfc_clear_ts (&ts);
1971
1972       ts.type = BT_REAL;
1973       ts.kind = gfc_default_real_kind;
1974
1975       result = gfc_copy_expr (a);
1976       if (!gfc_convert_boz (result, &ts))
1977         {
1978           gfc_free_expr (result);
1979           return &gfc_bad_expr;
1980         }
1981     }
1982   else
1983     result = gfc_int2real (a, gfc_default_real_kind);
1984   return range_check (result, "FLOAT");
1985 }
1986
1987
1988 gfc_expr *
1989 gfc_simplify_floor (gfc_expr *e, gfc_expr *k)
1990 {
1991   gfc_expr *result;
1992   mpfr_t floor;
1993   int kind;
1994
1995   kind = get_kind (BT_INTEGER, k, "FLOOR", gfc_default_integer_kind);
1996   if (kind == -1)
1997     gfc_internal_error ("gfc_simplify_floor(): Bad kind");
1998
1999   if (e->expr_type != EXPR_CONSTANT)
2000     return NULL;
2001
2002   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
2003
2004   gfc_set_model_kind (kind);
2005   mpfr_init (floor);
2006   mpfr_floor (floor, e->value.real);
2007
2008   gfc_mpfr_to_mpz (result->value.integer, floor, &e->where);
2009
2010   mpfr_clear (floor);
2011
2012   return range_check (result, "FLOOR");
2013 }
2014
2015
2016 gfc_expr *
2017 gfc_simplify_fraction (gfc_expr *x)
2018 {
2019   gfc_expr *result;
2020   mpfr_t absv, exp, pow2;
2021
2022   if (x->expr_type != EXPR_CONSTANT)
2023     return NULL;
2024
2025   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
2026
2027   if (mpfr_sgn (x->value.real) == 0)
2028     {
2029       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
2030       return result;
2031     }
2032
2033   gfc_set_model_kind (x->ts.kind);
2034   mpfr_init (exp);
2035   mpfr_init (absv);
2036   mpfr_init (pow2);
2037
2038   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
2039   mpfr_log2 (exp, absv, GFC_RND_MODE);
2040
2041   mpfr_trunc (exp, exp);
2042   mpfr_add_ui (exp, exp, 1, GFC_RND_MODE);
2043
2044   mpfr_ui_pow (pow2, 2, exp, GFC_RND_MODE);
2045
2046   mpfr_div (result->value.real, absv, pow2, GFC_RND_MODE);
2047
2048   mpfr_clears (exp, absv, pow2, NULL);
2049
2050   return range_check (result, "FRACTION");
2051 }
2052
2053
2054 gfc_expr *
2055 gfc_simplify_gamma (gfc_expr *x)
2056 {
2057   gfc_expr *result;
2058
2059   if (x->expr_type != EXPR_CONSTANT)
2060     return NULL;
2061
2062   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
2063
2064   mpfr_gamma (result->value.real, x->value.real, GFC_RND_MODE);
2065
2066   return range_check (result, "GAMMA");
2067 }
2068
2069
2070 gfc_expr *
2071 gfc_simplify_huge (gfc_expr *e)
2072 {
2073   gfc_expr *result;
2074   int i;
2075
2076   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
2077
2078   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
2079
2080   switch (e->ts.type)
2081     {
2082     case BT_INTEGER:
2083       mpz_set (result->value.integer, gfc_integer_kinds[i].huge);
2084       break;
2085
2086     case BT_REAL:
2087       mpfr_set (result->value.real, gfc_real_kinds[i].huge, GFC_RND_MODE);
2088       break;
2089
2090     default:
2091       gcc_unreachable ();
2092     }
2093
2094   return result;
2095 }
2096
2097
2098 gfc_expr *
2099 gfc_simplify_hypot (gfc_expr *x, gfc_expr *y)
2100 {
2101   gfc_expr *result;
2102
2103   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
2104     return NULL;
2105
2106   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
2107   mpfr_hypot (result->value.real, x->value.real, y->value.real, GFC_RND_MODE);
2108   return range_check (result, "HYPOT");
2109 }
2110
2111
2112 /* We use the processor's collating sequence, because all
2113    systems that gfortran currently works on are ASCII.  */
2114
2115 gfc_expr *
2116 gfc_simplify_iachar (gfc_expr *e, gfc_expr *kind)
2117 {
2118   gfc_expr *result;
2119   gfc_char_t index;
2120
2121   if (e->expr_type != EXPR_CONSTANT)
2122     return NULL;
2123
2124   if (e->value.character.length != 1)
2125     {
2126       gfc_error ("Argument of IACHAR at %L must be of length one", &e->where);
2127       return &gfc_bad_expr;
2128     }
2129
2130   index = e->value.character.string[0];
2131
2132   if (gfc_option.warn_surprising && index > 127)
2133     gfc_warning ("Argument of IACHAR function at %L outside of range 0..127",
2134                  &e->where);
2135
2136   if ((result = int_expr_with_kind (index, kind, "IACHAR")) == NULL)
2137     return &gfc_bad_expr;
2138
2139   result->where = e->where;
2140
2141   return range_check (result, "IACHAR");
2142 }
2143
2144
2145 gfc_expr *
2146 gfc_simplify_iand (gfc_expr *x, gfc_expr *y)
2147 {
2148   gfc_expr *result;
2149
2150   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
2151     return NULL;
2152
2153   result = gfc_constant_result (BT_INTEGER, x->ts.kind, &x->where);
2154
2155   mpz_and (result->value.integer, x->value.integer, y->value.integer);
2156
2157   return range_check (result, "IAND");
2158 }
2159
2160
2161 gfc_expr *
2162 gfc_simplify_ibclr (gfc_expr *x, gfc_expr *y)
2163 {
2164   gfc_expr *result;
2165   int k, pos;
2166
2167   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
2168     return NULL;
2169
2170   if (gfc_extract_int (y, &pos) != NULL || pos < 0)
2171     {
2172       gfc_error ("Invalid second argument of IBCLR at %L", &y->where);
2173       return &gfc_bad_expr;
2174     }
2175
2176   k = gfc_validate_kind (x->ts.type, x->ts.kind, false);
2177
2178   if (pos >= gfc_integer_kinds[k].bit_size)
2179     {
2180       gfc_error ("Second argument of IBCLR exceeds bit size at %L",
2181                  &y->where);
2182       return &gfc_bad_expr;
2183     }
2184
2185   result = gfc_copy_expr (x);
2186
2187   convert_mpz_to_unsigned (result->value.integer,
2188                            gfc_integer_kinds[k].bit_size);
2189
2190   mpz_clrbit (result->value.integer, pos);
2191
2192   convert_mpz_to_signed (result->value.integer,
2193                          gfc_integer_kinds[k].bit_size);
2194
2195   return result;
2196 }
2197
2198
2199 gfc_expr *
2200 gfc_simplify_ibits (gfc_expr *x, gfc_expr *y, gfc_expr *z)
2201 {
2202   gfc_expr *result;
2203   int pos, len;
2204   int i, k, bitsize;
2205   int *bits;
2206
2207   if (x->expr_type != EXPR_CONSTANT
2208       || y->expr_type != EXPR_CONSTANT
2209       || z->expr_type != EXPR_CONSTANT)
2210     return NULL;
2211
2212   if (gfc_extract_int (y, &pos) != NULL || pos < 0)
2213     {
2214       gfc_error ("Invalid second argument of IBITS at %L", &y->where);
2215       return &gfc_bad_expr;
2216     }
2217
2218   if (gfc_extract_int (z, &len) != NULL || len < 0)
2219     {
2220       gfc_error ("Invalid third argument of IBITS at %L", &z->where);
2221       return &gfc_bad_expr;
2222     }
2223
2224   k = gfc_validate_kind (BT_INTEGER, x->ts.kind, false);
2225
2226   bitsize = gfc_integer_kinds[k].bit_size;
2227
2228   if (pos + len > bitsize)
2229     {
2230       gfc_error ("Sum of second and third arguments of IBITS exceeds "
2231                  "bit size at %L", &y->where);
2232       return &gfc_bad_expr;
2233     }
2234
2235   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
2236   convert_mpz_to_unsigned (result->value.integer,
2237                            gfc_integer_kinds[k].bit_size);
2238
2239   bits = XCNEWVEC (int, bitsize);
2240
2241   for (i = 0; i < bitsize; i++)
2242     bits[i] = 0;
2243
2244   for (i = 0; i < len; i++)
2245     bits[i] = mpz_tstbit (x->value.integer, i + pos);
2246
2247   for (i = 0; i < bitsize; i++)
2248     {
2249       if (bits[i] == 0)
2250         mpz_clrbit (result->value.integer, i);
2251       else if (bits[i] == 1)
2252         mpz_setbit (result->value.integer, i);
2253       else
2254         gfc_internal_error ("IBITS: Bad bit");
2255     }
2256
2257   gfc_free (bits);
2258
2259   convert_mpz_to_signed (result->value.integer,
2260                          gfc_integer_kinds[k].bit_size);
2261
2262   return result;
2263 }
2264
2265
2266 gfc_expr *
2267 gfc_simplify_ibset (gfc_expr *x, gfc_expr *y)
2268 {
2269   gfc_expr *result;
2270   int k, pos;
2271
2272   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
2273     return NULL;
2274
2275   if (gfc_extract_int (y, &pos) != NULL || pos < 0)
2276     {
2277       gfc_error ("Invalid second argument of IBSET at %L", &y->where);
2278       return &gfc_bad_expr;
2279     }
2280
2281   k = gfc_validate_kind (x->ts.type, x->ts.kind, false);
2282
2283   if (pos >= gfc_integer_kinds[k].bit_size)
2284     {
2285       gfc_error ("Second argument of IBSET exceeds bit size at %L",
2286                  &y->where);
2287       return &gfc_bad_expr;
2288     }
2289
2290   result = gfc_copy_expr (x);
2291
2292   convert_mpz_to_unsigned (result->value.integer,
2293                            gfc_integer_kinds[k].bit_size);
2294
2295   mpz_setbit (result->value.integer, pos);
2296
2297   convert_mpz_to_signed (result->value.integer,
2298                          gfc_integer_kinds[k].bit_size);
2299
2300   return result;
2301 }
2302
2303
2304 gfc_expr *
2305 gfc_simplify_ichar (gfc_expr *e, gfc_expr *kind)
2306 {
2307   gfc_expr *result;
2308   gfc_char_t index;
2309
2310   if (e->expr_type != EXPR_CONSTANT)
2311     return NULL;
2312
2313   if (e->value.character.length != 1)
2314     {
2315       gfc_error ("Argument of ICHAR at %L must be of length one", &e->where);
2316       return &gfc_bad_expr;
2317     }
2318
2319   index = e->value.character.string[0];
2320
2321   if ((result = int_expr_with_kind (index, kind, "ICHAR")) == NULL)
2322     return &gfc_bad_expr;
2323
2324   result->where = e->where;
2325   return range_check (result, "ICHAR");
2326 }
2327
2328
2329 gfc_expr *
2330 gfc_simplify_ieor (gfc_expr *x, gfc_expr *y)
2331 {
2332   gfc_expr *result;
2333
2334   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
2335     return NULL;
2336
2337   result = gfc_constant_result (BT_INTEGER, x->ts.kind, &x->where);
2338
2339   mpz_xor (result->value.integer, x->value.integer, y->value.integer);
2340
2341   return range_check (result, "IEOR");
2342 }
2343
2344
2345 gfc_expr *
2346 gfc_simplify_index (gfc_expr *x, gfc_expr *y, gfc_expr *b, gfc_expr *kind)
2347 {
2348   gfc_expr *result;
2349   int back, len, lensub;
2350   int i, j, k, count, index = 0, start;
2351
2352   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT 
2353       || ( b != NULL && b->expr_type !=  EXPR_CONSTANT))
2354     return NULL;
2355
2356   if (b != NULL && b->value.logical != 0)
2357     back = 1;
2358   else
2359     back = 0;
2360
2361   k = get_kind (BT_INTEGER, kind, "INDEX", gfc_default_integer_kind); 
2362   if (k == -1)
2363     return &gfc_bad_expr;
2364
2365   result = gfc_constant_result (BT_INTEGER, k, &x->where);
2366
2367   len = x->value.character.length;
2368   lensub = y->value.character.length;
2369
2370   if (len < lensub)
2371     {
2372       mpz_set_si (result->value.integer, 0);
2373       return result;
2374     }
2375
2376   if (back == 0)
2377     {
2378       if (lensub == 0)
2379         {
2380           mpz_set_si (result->value.integer, 1);
2381           return result;
2382         }
2383       else if (lensub == 1)
2384         {
2385           for (i = 0; i < len; i++)
2386             {
2387               for (j = 0; j < lensub; j++)
2388                 {
2389                   if (y->value.character.string[j]
2390                       == x->value.character.string[i])
2391                     {
2392                       index = i + 1;
2393                       goto done;
2394                     }
2395                 }
2396             }
2397         }
2398       else
2399         {
2400           for (i = 0; i < len; i++)
2401             {
2402               for (j = 0; j < lensub; j++)
2403                 {
2404                   if (y->value.character.string[j]
2405                       == x->value.character.string[i])
2406                     {
2407                       start = i;
2408                       count = 0;
2409
2410                       for (k = 0; k < lensub; k++)
2411                         {
2412                           if (y->value.character.string[k]
2413                               == x->value.character.string[k + start])
2414                             count++;
2415                         }
2416
2417                       if (count == lensub)
2418                         {
2419                           index = start + 1;
2420                           goto done;
2421                         }
2422                     }
2423                 }
2424             }
2425         }
2426
2427     }
2428   else
2429     {
2430       if (lensub == 0)
2431         {
2432           mpz_set_si (result->value.integer, len + 1);
2433           return result;
2434         }
2435       else if (lensub == 1)
2436         {
2437           for (i = 0; i < len; i++)
2438             {
2439               for (j = 0; j < lensub; j++)
2440                 {
2441                   if (y->value.character.string[j]
2442                       == x->value.character.string[len - i])
2443                     {
2444                       index = len - i + 1;
2445                       goto done;
2446                     }
2447                 }
2448             }
2449         }
2450       else
2451         {
2452           for (i = 0; i < len; i++)
2453             {
2454               for (j = 0; j < lensub; j++)
2455                 {
2456                   if (y->value.character.string[j]
2457                       == x->value.character.string[len - i])
2458                     {
2459                       start = len - i;
2460                       if (start <= len - lensub)
2461                         {
2462                           count = 0;
2463                           for (k = 0; k < lensub; k++)
2464                             if (y->value.character.string[k]
2465                                 == x->value.character.string[k + start])
2466                               count++;
2467
2468                           if (count == lensub)
2469                             {
2470                               index = start + 1;
2471                               goto done;
2472                             }
2473                         }
2474                       else
2475                         {
2476                           continue;
2477                         }
2478                     }
2479                 }
2480             }
2481         }
2482     }
2483
2484 done:
2485   mpz_set_si (result->value.integer, index);
2486   return range_check (result, "INDEX");
2487 }
2488
2489
2490 gfc_expr *
2491 gfc_simplify_int (gfc_expr *e, gfc_expr *k)
2492 {
2493   gfc_expr *result = NULL;
2494   int kind;
2495
2496   kind = get_kind (BT_INTEGER, k, "INT", gfc_default_integer_kind);
2497   if (kind == -1)
2498     return &gfc_bad_expr;
2499
2500   if (e->expr_type != EXPR_CONSTANT)
2501     return NULL;
2502
2503   switch (e->ts.type)
2504     {
2505     case BT_INTEGER:
2506       result = gfc_int2int (e, kind);
2507       break;
2508
2509     case BT_REAL:
2510       result = gfc_real2int (e, kind);
2511       break;
2512
2513     case BT_COMPLEX:
2514       result = gfc_complex2int (e, kind);
2515       break;
2516
2517     default:
2518       gfc_error ("Argument of INT at %L is not a valid type", &e->where);
2519       return &gfc_bad_expr;
2520     }
2521
2522   return range_check (result, "INT");
2523 }
2524
2525
2526 static gfc_expr *
2527 simplify_intconv (gfc_expr *e, int kind, const char *name)
2528 {
2529   gfc_expr *result = NULL;
2530
2531   if (e->expr_type != EXPR_CONSTANT)
2532     return NULL;
2533
2534   switch (e->ts.type)
2535     {
2536     case BT_INTEGER:
2537       result = gfc_int2int (e, kind);
2538       break;
2539
2540     case BT_REAL:
2541       result = gfc_real2int (e, kind);
2542       break;
2543
2544     case BT_COMPLEX:
2545       result = gfc_complex2int (e, kind);
2546       break;
2547
2548     default:
2549       gfc_error ("Argument of %s at %L is not a valid type", name, &e->where);
2550       return &gfc_bad_expr;
2551     }
2552
2553   return range_check (result, name);
2554 }
2555
2556
2557 gfc_expr *
2558 gfc_simplify_int2 (gfc_expr *e)
2559 {
2560   return simplify_intconv (e, 2, "INT2");
2561 }
2562
2563
2564 gfc_expr *
2565 gfc_simplify_int8 (gfc_expr *e)
2566 {
2567   return simplify_intconv (e, 8, "INT8");
2568 }
2569
2570
2571 gfc_expr *
2572 gfc_simplify_long (gfc_expr *e)
2573 {
2574   return simplify_intconv (e, 4, "LONG");
2575 }
2576
2577
2578 gfc_expr *
2579 gfc_simplify_ifix (gfc_expr *e)
2580 {
2581   gfc_expr *rtrunc, *result;
2582
2583   if (e->expr_type != EXPR_CONSTANT)
2584     return NULL;
2585
2586   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
2587                                 &e->where);
2588
2589   rtrunc = gfc_copy_expr (e);
2590
2591   mpfr_trunc (rtrunc->value.real, e->value.real);
2592   gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real, &e->where);
2593
2594   gfc_free_expr (rtrunc);
2595   return range_check (result, "IFIX");
2596 }
2597
2598
2599 gfc_expr *
2600 gfc_simplify_idint (gfc_expr *e)
2601 {
2602   gfc_expr *rtrunc, *result;
2603
2604   if (e->expr_type != EXPR_CONSTANT)
2605     return NULL;
2606
2607   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
2608                                 &e->where);
2609
2610   rtrunc = gfc_copy_expr (e);
2611
2612   mpfr_trunc (rtrunc->value.real, e->value.real);
2613   gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real, &e->where);
2614
2615   gfc_free_expr (rtrunc);
2616   return range_check (result, "IDINT");
2617 }
2618
2619
2620 gfc_expr *
2621 gfc_simplify_ior (gfc_expr *x, gfc_expr *y)
2622 {
2623   gfc_expr *result;
2624
2625   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
2626     return NULL;
2627
2628   result = gfc_constant_result (BT_INTEGER, x->ts.kind, &x->where);
2629
2630   mpz_ior (result->value.integer, x->value.integer, y->value.integer);
2631   return range_check (result, "IOR");
2632 }
2633
2634
2635 gfc_expr *
2636 gfc_simplify_is_iostat_end (gfc_expr *x)
2637 {
2638   gfc_expr *result;
2639
2640   if (x->expr_type != EXPR_CONSTANT)
2641     return NULL;
2642
2643   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
2644                                 &x->where);
2645   result->value.logical = (mpz_cmp_si (x->value.integer, LIBERROR_END) == 0);
2646
2647   return result;
2648 }
2649
2650
2651 gfc_expr *
2652 gfc_simplify_is_iostat_eor (gfc_expr *x)
2653 {
2654   gfc_expr *result;
2655
2656   if (x->expr_type != EXPR_CONSTANT)
2657     return NULL;
2658
2659   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
2660                                 &x->where);
2661   result->value.logical = (mpz_cmp_si (x->value.integer, LIBERROR_EOR) == 0);
2662
2663   return result;
2664 }
2665
2666
2667 gfc_expr *
2668 gfc_simplify_isnan (gfc_expr *x)
2669 {
2670   gfc_expr *result;
2671
2672   if (x->expr_type != EXPR_CONSTANT)
2673     return NULL;
2674
2675   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
2676                                 &x->where);
2677   result->value.logical = mpfr_nan_p (x->value.real);
2678
2679   return result;
2680 }
2681
2682
2683 gfc_expr *
2684 gfc_simplify_ishft (gfc_expr *e, gfc_expr *s)
2685 {
2686   gfc_expr *result;
2687   int shift, ashift, isize, k, *bits, i;
2688
2689   if (e->expr_type != EXPR_CONSTANT || s->expr_type != EXPR_CONSTANT)
2690     return NULL;
2691
2692   if (gfc_extract_int (s, &shift) != NULL)
2693     {
2694       gfc_error ("Invalid second argument of ISHFT at %L", &s->where);
2695       return &gfc_bad_expr;
2696     }
2697
2698   k = gfc_validate_kind (BT_INTEGER, e->ts.kind, false);
2699
2700   isize = gfc_integer_kinds[k].bit_size;
2701
2702   if (shift >= 0)
2703     ashift = shift;
2704   else
2705     ashift = -shift;
2706
2707   if (ashift > isize)
2708     {
2709       gfc_error ("Magnitude of second argument of ISHFT exceeds bit size "
2710                  "at %L", &s->where);
2711       return &gfc_bad_expr;
2712     }
2713
2714   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
2715
2716   if (shift == 0)
2717     {
2718       mpz_set (result->value.integer, e->value.integer);
2719       return range_check (result, "ISHFT");
2720     }
2721   
2722   bits = XCNEWVEC (int, isize);
2723
2724   for (i = 0; i < isize; i++)
2725     bits[i] = mpz_tstbit (e->value.integer, i);
2726
2727   if (shift > 0)
2728     {
2729       for (i = 0; i < shift; i++)
2730         mpz_clrbit (result->value.integer, i);
2731
2732       for (i = 0; i < isize - shift; i++)
2733         {
2734           if (bits[i] == 0)
2735             mpz_clrbit (result->value.integer, i + shift);
2736           else
2737             mpz_setbit (result->value.integer, i + shift);
2738         }
2739     }
2740   else
2741     {
2742       for (i = isize - 1; i >= isize - ashift; i--)
2743         mpz_clrbit (result->value.integer, i);
2744
2745       for (i = isize - 1; i >= ashift; i--)
2746         {
2747           if (bits[i] == 0)
2748             mpz_clrbit (result->value.integer, i - ashift);
2749           else
2750             mpz_setbit (result->value.integer, i - ashift);
2751         }
2752     }
2753
2754   convert_mpz_to_signed (result->value.integer, isize);
2755
2756   gfc_free (bits);
2757   return result;
2758 }
2759
2760
2761 gfc_expr *
2762 gfc_simplify_ishftc (gfc_expr *e, gfc_expr *s, gfc_expr *sz)
2763 {
2764   gfc_expr *result;
2765   int shift, ashift, isize, ssize, delta, k;
2766   int i, *bits;
2767
2768   if (e->expr_type != EXPR_CONSTANT || s->expr_type != EXPR_CONSTANT)
2769     return NULL;
2770
2771   if (gfc_extract_int (s, &shift) != NULL)
2772     {
2773       gfc_error ("Invalid second argument of ISHFTC at %L", &s->where);
2774       return &gfc_bad_expr;
2775     }
2776
2777   k = gfc_validate_kind (e->ts.type, e->ts.kind, false);
2778   isize = gfc_integer_kinds[k].bit_size;
2779
2780   if (sz != NULL)
2781     {
2782       if (sz->expr_type != EXPR_CONSTANT)
2783         return NULL;
2784
2785       if (gfc_extract_int (sz, &ssize) != NULL || ssize <= 0)
2786         {
2787           gfc_error ("Invalid third argument of ISHFTC at %L", &sz->where);
2788           return &gfc_bad_expr;
2789         }
2790
2791       if (ssize > isize)
2792         {
2793           gfc_error ("Magnitude of third argument of ISHFTC exceeds "
2794                      "BIT_SIZE of first argument at %L", &s->where);
2795           return &gfc_bad_expr;
2796         }
2797     }
2798   else
2799     ssize = isize;
2800
2801   if (shift >= 0)
2802     ashift = shift;
2803   else
2804     ashift = -shift;
2805
2806   if (ashift > ssize)
2807     {
2808       if (sz != NULL)
2809         gfc_error ("Magnitude of second argument of ISHFTC exceeds "
2810                    "third argument at %L", &s->where);
2811       else
2812         gfc_error ("Magnitude of second argument of ISHFTC exceeds "
2813                    "BIT_SIZE of first argument at %L", &s->where);
2814       return &gfc_bad_expr;
2815     }
2816
2817   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
2818
2819   mpz_set (result->value.integer, e->value.integer);
2820
2821   if (shift == 0)
2822     return result;
2823
2824   convert_mpz_to_unsigned (result->value.integer, isize);
2825
2826   bits = XCNEWVEC (int, ssize);
2827
2828   for (i = 0; i < ssize; i++)
2829     bits[i] = mpz_tstbit (e->value.integer, i);
2830
2831   delta = ssize - ashift;
2832
2833   if (shift > 0)
2834     {
2835       for (i = 0; i < delta; i++)
2836         {
2837           if (bits[i] == 0)
2838             mpz_clrbit (result->value.integer, i + shift);
2839           else
2840             mpz_setbit (result->value.integer, i + shift);
2841         }
2842
2843       for (i = delta; i < ssize; i++)
2844         {
2845           if (bits[i] == 0)
2846             mpz_clrbit (result->value.integer, i - delta);
2847           else
2848             mpz_setbit (result->value.integer, i - delta);
2849         }
2850     }
2851   else
2852     {
2853       for (i = 0; i < ashift; i++)
2854         {
2855           if (bits[i] == 0)
2856             mpz_clrbit (result->value.integer, i + delta);
2857           else
2858             mpz_setbit (result->value.integer, i + delta);
2859         }
2860
2861       for (i = ashift; i < ssize; i++)
2862         {
2863           if (bits[i] == 0)
2864             mpz_clrbit (result->value.integer, i + shift);
2865           else
2866             mpz_setbit (result->value.integer, i + shift);
2867         }
2868     }
2869
2870   convert_mpz_to_signed (result->value.integer, isize);
2871
2872   gfc_free (bits);
2873   return result;
2874 }
2875
2876
2877 gfc_expr *
2878 gfc_simplify_kind (gfc_expr *e)
2879 {
2880
2881   if (e->ts.type == BT_DERIVED)
2882     {
2883       gfc_error ("Argument of KIND at %L is a DERIVED type", &e->where);
2884       return &gfc_bad_expr;
2885     }
2886
2887   return gfc_int_expr (e->ts.kind);
2888 }
2889
2890
2891 static gfc_expr *
2892 simplify_bound_dim (gfc_expr *array, gfc_expr *kind, int d, int upper,
2893                     gfc_array_spec *as, gfc_ref *ref)
2894 {
2895   gfc_expr *l, *u, *result;
2896   int k;
2897
2898   /* The last dimension of an assumed-size array is special.  */
2899   if (d == as->rank && as->type == AS_ASSUMED_SIZE && !upper)
2900     {
2901       if (as->lower[d-1]->expr_type == EXPR_CONSTANT)
2902         return gfc_copy_expr (as->lower[d-1]);
2903       else
2904         return NULL;
2905     }
2906
2907   k = get_kind (BT_INTEGER, kind, upper ? "UBOUND" : "LBOUND",
2908                 gfc_default_integer_kind); 
2909   if (k == -1)
2910     return &gfc_bad_expr;
2911
2912   result = gfc_constant_result (BT_INTEGER, k, &array->where);
2913
2914
2915   /* Then, we need to know the extent of the given dimension.  */
2916   if (ref->u.ar.type == AR_FULL)
2917     {
2918       l = as->lower[d-1];
2919       u = as->upper[d-1];
2920
2921       if (l->expr_type != EXPR_CONSTANT || u->expr_type != EXPR_CONSTANT)
2922         return NULL;
2923
2924       if (mpz_cmp (l->value.integer, u->value.integer) > 0)
2925         {
2926           /* Zero extent.  */
2927           if (upper)
2928             mpz_set_si (result->value.integer, 0);
2929           else
2930             mpz_set_si (result->value.integer, 1);
2931         }
2932       else
2933         {
2934           /* Nonzero extent.  */
2935           if (upper)
2936             mpz_set (result->value.integer, u->value.integer);
2937           else
2938             mpz_set (result->value.integer, l->value.integer);
2939         }
2940     }
2941   else
2942     {
2943       if (upper)
2944         {
2945           if (gfc_ref_dimen_size (&ref->u.ar, d-1, &result->value.integer)
2946               != SUCCESS)
2947             return NULL;
2948         }
2949       else
2950         mpz_set_si (result->value.integer, (long int) 1);
2951     }
2952
2953   return range_check (result, upper ? "UBOUND" : "LBOUND");
2954 }
2955
2956
2957 static gfc_expr *
2958 simplify_bound (gfc_expr *array, gfc_expr *dim, gfc_expr *kind, int upper)
2959 {
2960   gfc_ref *ref;
2961   gfc_array_spec *as;
2962   int d;
2963
2964   if (array->expr_type != EXPR_VARIABLE)
2965     return NULL;
2966
2967   /* Follow any component references.  */
2968   as = array->symtree->n.sym->as;
2969   for (ref = array->ref; ref; ref = ref->next)
2970     {
2971       switch (ref->type)
2972         {
2973         case REF_ARRAY:
2974           switch (ref->u.ar.type)
2975             {
2976             case AR_ELEMENT:
2977               as = NULL;
2978               continue;
2979
2980             case AR_FULL:
2981               /* We're done because 'as' has already been set in the
2982                  previous iteration.  */
2983               if (!ref->next)
2984                 goto done;
2985
2986             /* Fall through.  */
2987
2988             case AR_UNKNOWN:
2989               return NULL;
2990
2991             case AR_SECTION:
2992               as = ref->u.ar.as;
2993               goto done;
2994             }
2995
2996           gcc_unreachable ();
2997
2998         case REF_COMPONENT:
2999           as = ref->u.c.component->as;
3000           continue;
3001
3002         case REF_SUBSTRING:
3003           continue;
3004         }
3005     }
3006
3007   gcc_unreachable ();
3008
3009  done:
3010
3011   if (as->type == AS_DEFERRED || as->type == AS_ASSUMED_SHAPE)
3012     return NULL;
3013
3014   if (dim == NULL)
3015     {
3016       /* Multi-dimensional bounds.  */
3017       gfc_expr *bounds[GFC_MAX_DIMENSIONS];
3018       gfc_expr *e;
3019       gfc_constructor *head, *tail;
3020       int k;
3021
3022       /* UBOUND(ARRAY) is not valid for an assumed-size array.  */
3023       if (upper && as->type == AS_ASSUMED_SIZE)
3024         {
3025           /* An error message will be emitted in
3026              check_assumed_size_reference (resolve.c).  */
3027           return &gfc_bad_expr;
3028         }
3029
3030       /* Simplify the bounds for each dimension.  */
3031       for (d = 0; d < array->rank; d++)
3032         {
3033           bounds[d] = simplify_bound_dim (array, kind, d + 1, upper, as, ref);
3034           if (bounds[d] == NULL || bounds[d] == &gfc_bad_expr)
3035             {
3036               int j;
3037
3038               for (j = 0; j < d; j++)
3039                 gfc_free_expr (bounds[j]);
3040               return bounds[d];
3041             }
3042         }
3043
3044       /* Allocate the result expression.  */
3045       e = gfc_get_expr ();
3046       e->where = array->where;
3047       e->expr_type = EXPR_ARRAY;
3048       e->ts.type = BT_INTEGER;
3049       k = get_kind (BT_INTEGER, kind, upper ? "UBOUND" : "LBOUND",
3050                     gfc_default_integer_kind); 
3051       if (k == -1)
3052         {
3053           gfc_free_expr (e);
3054           return &gfc_bad_expr;
3055         }
3056       e->ts.kind = k;
3057
3058       /* The result is a rank 1 array; its size is the rank of the first
3059          argument to {L,U}BOUND.  */
3060       e->rank = 1;
3061       e->shape = gfc_get_shape (1);
3062       mpz_init_set_ui (e->shape[0], array->rank);
3063
3064       /* Create the constructor for this array.  */
3065       head = tail = NULL;
3066       for (d = 0; d < array->rank; d++)
3067         {
3068           /* Get a new constructor element.  */
3069           if (head == NULL)
3070             head = tail = gfc_get_constructor ();
3071           else
3072             {
3073               tail->next = gfc_get_constructor ();
3074               tail = tail->next;
3075             }
3076
3077           tail->where = e->where;
3078           tail->expr = bounds[d];
3079         }
3080       e->value.constructor = head;
3081
3082       return e;
3083     }
3084   else
3085     {
3086       /* A DIM argument is specified.  */
3087       if (dim->expr_type != EXPR_CONSTANT)
3088         return NULL;
3089
3090       d = mpz_get_si (dim->value.integer);
3091
3092       if (d < 1 || d > as->rank
3093           || (d == as->rank && as->type == AS_ASSUMED_SIZE && upper))
3094         {
3095           gfc_error ("DIM argument at %L is out of bounds", &dim->where);
3096           return &gfc_bad_expr;
3097         }
3098
3099       return simplify_bound_dim (array, kind, d, upper, as, ref);
3100     }
3101 }
3102
3103
3104 gfc_expr *
3105 gfc_simplify_lbound (gfc_expr *array, gfc_expr *dim, gfc_expr *kind)
3106 {
3107   return simplify_bound (array, dim, kind, 0);
3108 }
3109
3110
3111 gfc_expr *
3112 gfc_simplify_leadz (gfc_expr *e)
3113 {
3114   gfc_expr *result;
3115   unsigned long lz, bs;
3116   int i;
3117
3118   if (e->expr_type != EXPR_CONSTANT)
3119     return NULL;
3120
3121   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
3122   bs = gfc_integer_kinds[i].bit_size;
3123   if (mpz_cmp_si (e->value.integer, 0) == 0)
3124     lz = bs;
3125   else if (mpz_cmp_si (e->value.integer, 0) < 0)
3126     lz = 0;
3127   else
3128     lz = bs - mpz_sizeinbase (e->value.integer, 2);
3129
3130   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
3131                                 &e->where);
3132   mpz_set_ui (result->value.integer, lz);
3133
3134   return result;
3135 }
3136
3137
3138 gfc_expr *
3139 gfc_simplify_len (gfc_expr *e, gfc_expr *kind)
3140 {
3141   gfc_expr *result;
3142   int k = get_kind (BT_INTEGER, kind, "LEN", gfc_default_integer_kind);
3143
3144   if (k == -1)
3145     return &gfc_bad_expr;
3146
3147   if (e->expr_type == EXPR_CONSTANT)
3148     {
3149       result = gfc_constant_result (BT_INTEGER, k, &e->where);
3150       mpz_set_si (result->value.integer, e->value.character.length);
3151       if (gfc_range_check (result) == ARITH_OK)
3152         return result;
3153       else
3154         {
3155           gfc_free_expr (result);
3156           return NULL;
3157         }
3158     }
3159
3160   if (e->ts.cl != NULL && e->ts.cl->length != NULL
3161       && e->ts.cl->length->expr_type == EXPR_CONSTANT
3162       && e->ts.cl->length->ts.type == BT_INTEGER)
3163     {
3164       result = gfc_constant_result (BT_INTEGER, k, &e->where);
3165       mpz_set (result->value.integer, e->ts.cl->length->value.integer);
3166       if (gfc_range_check (result) == ARITH_OK)
3167         return result;
3168       else
3169         {
3170           gfc_free_expr (result);
3171           return NULL;
3172         }
3173     }
3174
3175   return NULL;
3176 }
3177
3178
3179 gfc_expr *
3180 gfc_simplify_len_trim (gfc_expr *e, gfc_expr *kind)
3181 {
3182   gfc_expr *result;
3183   int count, len, lentrim, i;
3184   int k = get_kind (BT_INTEGER, kind, "LEN_TRIM", gfc_default_integer_kind);
3185
3186   if (k == -1)
3187     return &gfc_bad_expr;
3188
3189   if (e->expr_type != EXPR_CONSTANT)
3190     return NULL;
3191
3192   result = gfc_constant_result (BT_INTEGER, k, &e->where);
3193   len = e->value.character.length;
3194
3195   for (count = 0, i = 1; i <= len; i++)
3196     if (e->value.character.string[len - i] == ' ')
3197       count++;
3198     else
3199       break;
3200
3201   lentrim = len - count;
3202
3203   mpz_set_si (result->value.integer, lentrim);
3204   return range_check (result, "LEN_TRIM");
3205 }
3206
3207 gfc_expr *
3208 gfc_simplify_lgamma (gfc_expr *x ATTRIBUTE_UNUSED)
3209 {
3210   gfc_expr *result;
3211   int sg;
3212
3213   if (x->expr_type != EXPR_CONSTANT)
3214     return NULL;
3215
3216   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3217
3218   mpfr_lgamma (result->value.real, &sg, x->value.real, GFC_RND_MODE);
3219
3220   return range_check (result, "LGAMMA");
3221 }
3222
3223
3224 gfc_expr *
3225 gfc_simplify_lge (gfc_expr *a, gfc_expr *b)
3226 {
3227   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
3228     return NULL;
3229
3230   return gfc_logical_expr (gfc_compare_string (a, b) >= 0, &a->where);
3231 }
3232
3233
3234 gfc_expr *
3235 gfc_simplify_lgt (gfc_expr *a, gfc_expr *b)
3236 {
3237   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
3238     return NULL;
3239
3240   return gfc_logical_expr (gfc_compare_string (a, b) > 0,
3241                            &a->where);
3242 }
3243
3244
3245 gfc_expr *
3246 gfc_simplify_lle (gfc_expr *a, gfc_expr *b)
3247 {
3248   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
3249     return NULL;
3250
3251   return gfc_logical_expr (gfc_compare_string (a, b) <= 0, &a->where);
3252 }
3253
3254
3255 gfc_expr *
3256 gfc_simplify_llt (gfc_expr *a, gfc_expr *b)
3257 {
3258   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
3259     return NULL;
3260
3261   return gfc_logical_expr (gfc_compare_string (a, b) < 0, &a->where);
3262 }
3263
3264
3265 gfc_expr *
3266 gfc_simplify_log (gfc_expr *x)
3267 {
3268   gfc_expr *result;
3269
3270   if (x->expr_type != EXPR_CONSTANT)
3271     return NULL;
3272
3273   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3274
3275
3276   switch (x->ts.type)
3277     {
3278     case BT_REAL:
3279       if (mpfr_sgn (x->value.real) <= 0)
3280         {
3281           gfc_error ("Argument of LOG at %L cannot be less than or equal "
3282                      "to zero", &x->where);
3283           gfc_free_expr (result);
3284           return &gfc_bad_expr;
3285         }
3286
3287       mpfr_log (result->value.real, x->value.real, GFC_RND_MODE);
3288       break;
3289
3290     case BT_COMPLEX:
3291       if ((mpfr_sgn (mpc_realref (x->value.complex)) == 0)
3292           && (mpfr_sgn (mpc_imagref (x->value.complex)) == 0))
3293         {
3294           gfc_error ("Complex argument of LOG at %L cannot be zero",
3295                      &x->where);
3296           gfc_free_expr (result);
3297           return &gfc_bad_expr;
3298         }
3299
3300       gfc_set_model_kind (x->ts.kind);
3301 #ifdef HAVE_mpc
3302       mpc_log (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
3303 #else
3304     {
3305       mpfr_t xr, xi;
3306       mpfr_init (xr);
3307       mpfr_init (xi);
3308
3309       mpfr_atan2 (result->value.complex.i, x->value.complex.i,
3310                   x->value.complex.r, GFC_RND_MODE);
3311
3312       mpfr_mul (xr, x->value.complex.r, x->value.complex.r, GFC_RND_MODE);
3313       mpfr_mul (xi, x->value.complex.i, x->value.complex.i, GFC_RND_MODE);
3314       mpfr_add (xr, xr, xi, GFC_RND_MODE);
3315       mpfr_sqrt (xr, xr, GFC_RND_MODE);
3316       mpfr_log (result->value.complex.r, xr, GFC_RND_MODE);
3317
3318       mpfr_clears (xr, xi, NULL);
3319     }
3320 #endif
3321       break;
3322
3323     default:
3324       gfc_internal_error ("gfc_simplify_log: bad type");
3325     }
3326
3327   return range_check (result, "LOG");
3328 }
3329
3330
3331 gfc_expr *
3332 gfc_simplify_log10 (gfc_expr *x)
3333 {
3334   gfc_expr *result;
3335
3336   if (x->expr_type != EXPR_CONSTANT)
3337     return NULL;
3338
3339   if (mpfr_sgn (x->value.real) <= 0)
3340     {
3341       gfc_error ("Argument of LOG10 at %L cannot be less than or equal "
3342                  "to zero", &x->where);
3343       return &gfc_bad_expr;
3344     }
3345
3346   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3347
3348   mpfr_log10 (result->value.real, x->value.real, GFC_RND_MODE);
3349
3350   return range_check (result, "LOG10");
3351 }
3352
3353
3354 gfc_expr *
3355 gfc_simplify_logical (gfc_expr *e, gfc_expr *k)
3356 {
3357   gfc_expr *result;
3358   int kind;
3359
3360   kind = get_kind (BT_LOGICAL, k, "LOGICAL", gfc_default_logical_kind);
3361   if (kind < 0)
3362     return &gfc_bad_expr;
3363
3364   if (e->expr_type != EXPR_CONSTANT)
3365     return NULL;
3366
3367   result = gfc_constant_result (BT_LOGICAL, kind, &e->where);
3368
3369   result->value.logical = e->value.logical;
3370
3371   return result;
3372 }
3373
3374
3375 gfc_expr*
3376 gfc_simplify_matmul (gfc_expr *matrix_a, gfc_expr *matrix_b)
3377 {
3378   gfc_expr *result;
3379   gfc_constructor *ma_ctor, *mb_ctor;
3380   int row, result_rows, col, result_columns, stride_a, stride_b;
3381
3382   if (!is_constant_array_expr (matrix_a)
3383       || !is_constant_array_expr (matrix_b))
3384     return NULL;
3385
3386   gcc_assert (gfc_compare_types (&matrix_a->ts, &matrix_b->ts));
3387   result = gfc_start_constructor (matrix_a->ts.type,
3388                                   matrix_a->ts.kind,
3389                                   &matrix_a->where);
3390
3391   if (matrix_a->rank == 1 && matrix_b->rank == 2)
3392     {
3393       result_rows = 1;
3394       result_columns = mpz_get_si (matrix_b->shape[0]);
3395       stride_a = 1;
3396       stride_b = mpz_get_si (matrix_b->shape[0]);
3397
3398       result->rank = 1;
3399       result->shape = gfc_get_shape (result->rank);
3400       mpz_init_set_si (result->shape[0], result_columns);
3401     }
3402   else if (matrix_a->rank == 2 && matrix_b->rank == 1)
3403     {
3404       result_rows = mpz_get_si (matrix_b->shape[0]);
3405       result_columns = 1;
3406       stride_a = mpz_get_si (matrix_a->shape[0]);
3407       stride_b = 1;
3408
3409       result->rank = 1;
3410       result->shape = gfc_get_shape (result->rank);
3411       mpz_init_set_si (result->shape[0], result_rows);
3412     }
3413   else if (matrix_a->rank == 2 && matrix_b->rank == 2)
3414     {
3415       result_rows = mpz_get_si (matrix_a->shape[0]);
3416       result_columns = mpz_get_si (matrix_b->shape[1]);
3417       stride_a = mpz_get_si (matrix_a->shape[1]);
3418       stride_b = mpz_get_si (matrix_b->shape[0]);
3419
3420       result->rank = 2;
3421       result->shape = gfc_get_shape (result->rank);
3422       mpz_init_set_si (result->shape[0], result_rows);
3423       mpz_init_set_si (result->shape[1], result_columns);
3424     }
3425   else
3426     gcc_unreachable();
3427
3428   ma_ctor = matrix_a->value.constructor;
3429   mb_ctor = matrix_b->value.constructor;
3430
3431   for (col = 0; col < result_columns; ++col)
3432     {
3433       ma_ctor = matrix_a->value.constructor;
3434
3435       for (row = 0; row < result_rows; ++row)
3436         {
3437           gfc_expr *e;
3438           e = compute_dot_product (ma_ctor, stride_a,
3439                                    mb_ctor, 1);
3440
3441           gfc_append_constructor (result, e);
3442
3443           ADVANCE (ma_ctor, 1);
3444         }
3445
3446       ADVANCE (mb_ctor, stride_b);
3447     }
3448
3449   return result;
3450 }
3451
3452
3453 gfc_expr *
3454 gfc_simplify_merge (gfc_expr *tsource, gfc_expr *fsource, gfc_expr *mask)
3455 {
3456   if (tsource->expr_type != EXPR_CONSTANT
3457       || fsource->expr_type != EXPR_CONSTANT
3458       || mask->expr_type != EXPR_CONSTANT)
3459     return NULL;
3460
3461   return gfc_copy_expr (mask->value.logical ? tsource : fsource);
3462 }
3463
3464
3465 /* Selects bewteen current value and extremum for simplify_min_max
3466    and simplify_minval_maxval.  */
3467 static void
3468 min_max_choose (gfc_expr *arg, gfc_expr *extremum, int sign)
3469 {
3470   switch (arg->ts.type)
3471     {
3472       case BT_INTEGER:
3473         if (mpz_cmp (arg->value.integer,
3474                         extremum->value.integer) * sign > 0)
3475         mpz_set (extremum->value.integer, arg->value.integer);
3476         break;
3477
3478       case BT_REAL:
3479         /* We need to use mpfr_min and mpfr_max to treat NaN properly.  */
3480         if (sign > 0)
3481           mpfr_max (extremum->value.real, extremum->value.real,
3482                       arg->value.real, GFC_RND_MODE);
3483         else
3484           mpfr_min (extremum->value.real, extremum->value.real,
3485                       arg->value.real, GFC_RND_MODE);
3486         break;
3487
3488       case BT_CHARACTER:
3489 #define LENGTH(x) ((x)->value.character.length)
3490 #define STRING(x) ((x)->value.character.string)
3491         if (LENGTH(extremum) < LENGTH(arg))
3492           {
3493             gfc_char_t *tmp = STRING(extremum);
3494
3495             STRING(extremum) = gfc_get_wide_string (LENGTH(arg) + 1);
3496             memcpy (STRING(extremum), tmp,
3497                       LENGTH(extremum) * sizeof (gfc_char_t));
3498             gfc_wide_memset (&STRING(extremum)[LENGTH(extremum)], ' ',
3499                                LENGTH(arg) - LENGTH(extremum));
3500             STRING(extremum)[LENGTH(arg)] = '\0';  /* For debugger  */
3501             LENGTH(extremum) = LENGTH(arg);
3502             gfc_free (tmp);
3503           }
3504
3505         if (gfc_compare_string (arg, extremum) * sign > 0)
3506           {
3507             gfc_free (STRING(extremum));
3508             STRING(extremum) = gfc_get_wide_string (LENGTH(extremum) + 1);
3509             memcpy (STRING(extremum), STRING(arg),
3510                       LENGTH(arg) * sizeof (gfc_char_t));
3511             gfc_wide_memset (&STRING(extremum)[LENGTH(arg)], ' ',
3512                                LENGTH(extremum) - LENGTH(arg));
3513             STRING(extremum)[LENGTH(extremum)] = '\0';  /* For debugger  */
3514           }
3515 #undef LENGTH
3516 #undef STRING
3517         break;
3518               
3519       default:
3520         gfc_internal_error ("simplify_min_max(): Bad type in arglist");
3521     }
3522 }
3523
3524
3525 /* This function is special since MAX() can take any number of
3526    arguments.  The simplified expression is a rewritten version of the
3527    argument list containing at most one constant element.  Other
3528    constant elements are deleted.  Because the argument list has
3529    already been checked, this function always succeeds.  sign is 1 for
3530    MAX(), -1 for MIN().  */
3531
3532 static gfc_expr *
3533 simplify_min_max (gfc_expr *expr, int sign)
3534 {
3535   gfc_actual_arglist *arg, *last, *extremum;
3536   gfc_intrinsic_sym * specific;
3537
3538   last = NULL;
3539   extremum = NULL;
3540   specific = expr->value.function.isym;
3541
3542   arg = expr->value.function.actual;
3543
3544   for (; arg; last = arg, arg = arg->next)
3545     {
3546       if (arg->expr->expr_type != EXPR_CONSTANT)
3547         continue;
3548
3549       if (extremum == NULL)
3550         {
3551           extremum = arg;
3552           continue;
3553         }
3554
3555       min_max_choose (arg->expr, extremum->expr, sign);
3556
3557       /* Delete the extra constant argument.  */
3558       if (last == NULL)
3559         expr->value.function.actual = arg->next;
3560       else
3561         last->next = arg->next;
3562
3563       arg->next = NULL;
3564       gfc_free_actual_arglist (arg);
3565       arg = last;
3566     }
3567
3568   /* If there is one value left, replace the function call with the
3569      expression.  */
3570   if (expr->value.function.actual->next != NULL)
3571     return NULL;
3572
3573   /* Convert to the correct type and kind.  */
3574   if (expr->ts.type != BT_UNKNOWN) 
3575     return gfc_convert_constant (expr->value.function.actual->expr,
3576         expr->ts.type, expr->ts.kind);
3577
3578   if (specific->ts.type != BT_UNKNOWN) 
3579     return gfc_convert_constant (expr->value.function.actual->expr,
3580         specific->ts.type, specific->ts.kind); 
3581  
3582   return gfc_copy_expr (expr->value.function.actual->expr);
3583 }
3584
3585
3586 gfc_expr *
3587 gfc_simplify_min (gfc_expr *e)
3588 {
3589   return simplify_min_max (e, -1);
3590 }
3591
3592
3593 gfc_expr *
3594 gfc_simplify_max (gfc_expr *e)
3595 {
3596   return simplify_min_max (e, 1);
3597 }
3598
3599
3600 /* This is a simplified version of simplify_min_max to provide
3601    simplification of minval and maxval for a vector.  */
3602
3603 static gfc_expr *
3604 simplify_minval_maxval (gfc_expr *expr, int sign)
3605 {
3606   gfc_constructor *ctr, *extremum;
3607   gfc_intrinsic_sym * specific;
3608
3609   extremum = NULL;
3610   specific = expr->value.function.isym;
3611
3612   ctr = expr->value.constructor;
3613
3614   for (; ctr; ctr = ctr->next)
3615     {
3616       if (ctr->expr->expr_type != EXPR_CONSTANT)
3617         return NULL;
3618
3619       if (extremum == NULL)
3620         {
3621           extremum = ctr;
3622           continue;
3623         }
3624
3625       min_max_choose (ctr->expr, extremum->expr, sign);
3626      }
3627
3628   if (extremum == NULL)
3629     return NULL;
3630
3631   /* Convert to the correct type and kind.  */
3632   if (expr->ts.type != BT_UNKNOWN) 
3633     return gfc_convert_constant (extremum->expr,
3634         expr->ts.type, expr->ts.kind);
3635
3636   if (specific->ts.type != BT_UNKNOWN) 
3637     return gfc_convert_constant (extremum->expr,
3638         specific->ts.type, specific->ts.kind); 
3639  
3640   return gfc_copy_expr (extremum->expr);
3641 }
3642
3643
3644 gfc_expr *
3645 gfc_simplify_minval (gfc_expr *array, gfc_expr* dim, gfc_expr *mask)
3646 {
3647   if (array->expr_type != EXPR_ARRAY || array->rank != 1 || dim || mask)
3648     return NULL;
3649   
3650   return simplify_minval_maxval (array, -1);
3651 }
3652
3653
3654 gfc_expr *
3655 gfc_simplify_maxval (gfc_expr *array, gfc_expr* dim, gfc_expr *mask)
3656 {
3657   if (array->expr_type != EXPR_ARRAY || array->rank != 1 || dim || mask)
3658     return NULL;
3659   return simplify_minval_maxval (array, 1);
3660 }
3661
3662
3663 gfc_expr *
3664 gfc_simplify_maxexponent (gfc_expr *x)
3665 {
3666   gfc_expr *result;
3667   int i;
3668
3669   i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
3670
3671   result = gfc_int_expr (gfc_real_kinds[i].max_exponent);
3672   result->where = x->where;
3673
3674   return result;
3675 }
3676
3677
3678 gfc_expr *
3679 gfc_simplify_minexponent (gfc_expr *x)
3680 {
3681   gfc_expr *result;
3682   int i;
3683
3684   i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
3685
3686   result = gfc_int_expr (gfc_real_kinds[i].min_exponent);
3687   result->where = x->where;
3688
3689   return result;
3690 }
3691
3692
3693 gfc_expr *
3694 gfc_simplify_mod (gfc_expr *a, gfc_expr *p)
3695 {
3696   gfc_expr *result;
3697   mpfr_t tmp;
3698   int kind;
3699
3700   if (a->expr_type != EXPR_CONSTANT || p->expr_type != EXPR_CONSTANT)
3701     return NULL;
3702
3703   kind = a->ts.kind > p->ts.kind ? a->ts.kind : p->ts.kind;
3704   result = gfc_constant_result (a->ts.type, kind, &a->where);
3705
3706   switch (a->ts.type)
3707     {
3708     case BT_INTEGER:
3709       if (mpz_cmp_ui (p->value.integer, 0) == 0)
3710         {
3711           /* Result is processor-dependent.  */
3712           gfc_error ("Second argument MOD at %L is zero", &a->where);
3713           gfc_free_expr (result);
3714           return &gfc_bad_expr;
3715         }
3716       mpz_tdiv_r (result->value.integer, a->value.integer, p->value.integer);
3717       break;
3718
3719     case BT_REAL:
3720       if (mpfr_cmp_ui (p->value.real, 0) == 0)
3721         {
3722           /* Result is processor-dependent.  */
3723           gfc_error ("Second argument of MOD at %L is zero", &p->where);
3724           gfc_free_expr (result);
3725           return &gfc_bad_expr;
3726         }
3727
3728       gfc_set_model_kind (kind);
3729       mpfr_init (tmp);
3730       mpfr_div (tmp, a->value.real, p->value.real, GFC_RND_MODE);
3731       mpfr_trunc (tmp, tmp);
3732       mpfr_mul (tmp, tmp, p->value.real, GFC_RND_MODE);
3733       mpfr_sub (result->value.real, a->value.real, tmp, GFC_RND_MODE);
3734       mpfr_clear (tmp);
3735       break;
3736
3737     default:
3738       gfc_internal_error ("gfc_simplify_mod(): Bad arguments");
3739     }
3740
3741   return range_check (result, "MOD");
3742 }
3743
3744
3745 gfc_expr *
3746 gfc_simplify_modulo (gfc_expr *a, gfc_expr *p)
3747 {
3748   gfc_expr *result;
3749   mpfr_t tmp;
3750   int kind;
3751
3752   if (a->expr_type != EXPR_CONSTANT || p->expr_type != EXPR_CONSTANT)
3753     return NULL;
3754
3755   kind = a->ts.kind > p->ts.kind ? a->ts.kind : p->ts.kind;
3756   result = gfc_constant_result (a->ts.type, kind, &a->where);
3757
3758   switch (a->ts.type)
3759     {
3760     case BT_INTEGER:
3761       if (mpz_cmp_ui (p->value.integer, 0) == 0)
3762         {
3763           /* Result is processor-dependent. This processor just opts
3764              to not handle it at all.  */
3765           gfc_error ("Second argument of MODULO at %L is zero", &a->where);
3766           gfc_free_expr (result);
3767           return &gfc_bad_expr;
3768         }
3769       mpz_fdiv_r (result->value.integer, a->value.integer, p->value.integer);
3770
3771       break;
3772
3773     case BT_REAL:
3774       if (mpfr_cmp_ui (p->value.real, 0) == 0)
3775         {
3776           /* Result is processor-dependent.  */
3777           gfc_error ("Second argument of MODULO at %L is zero", &p->where);
3778           gfc_free_expr (result);
3779           return &gfc_bad_expr;
3780         }
3781
3782       gfc_set_model_kind (kind);
3783       mpfr_init (tmp);
3784       mpfr_div (tmp, a->value.real, p->value.real, GFC_RND_MODE);
3785       mpfr_floor (tmp, tmp);
3786       mpfr_mul (tmp, tmp, p->value.real, GFC_RND_MODE);
3787       mpfr_sub (result->value.real, a->value.real, tmp, GFC_RND_MODE);
3788       mpfr_clear (tmp);
3789       break;
3790
3791     default:
3792       gfc_internal_error ("gfc_simplify_modulo(): Bad arguments");
3793     }
3794
3795   return range_check (result, "MODULO");
3796 }
3797
3798
3799 /* Exists for the sole purpose of consistency with other intrinsics.  */
3800 gfc_expr *
3801 gfc_simplify_mvbits (gfc_expr *f  ATTRIBUTE_UNUSED,
3802                      gfc_expr *fp ATTRIBUTE_UNUSED,
3803                      gfc_expr *l  ATTRIBUTE_UNUSED,
3804                      gfc_expr *to ATTRIBUTE_UNUSED,
3805                      gfc_expr *tp ATTRIBUTE_UNUSED)
3806 {
3807   return NULL;
3808 }
3809
3810
3811 gfc_expr *
3812 gfc_simplify_nearest (gfc_expr *x, gfc_expr *s)
3813 {
3814   gfc_expr *result;
3815   mp_exp_t emin, emax;
3816   int kind;
3817
3818   if (x->expr_type != EXPR_CONSTANT || s->expr_type != EXPR_CONSTANT)
3819     return NULL;
3820
3821   if (mpfr_sgn (s->value.real) == 0)
3822     {
3823       gfc_error ("Second argument of NEAREST at %L shall not be zero",
3824                  &s->where);
3825       return &gfc_bad_expr;
3826     }
3827
3828   result = gfc_copy_expr (x);
3829
3830   /* Save current values of emin and emax.  */
3831   emin = mpfr_get_emin ();
3832   emax = mpfr_get_emax ();
3833
3834   /* Set emin and emax for the current model number.  */
3835   kind = gfc_validate_kind (BT_REAL, x->ts.kind, 0);
3836   mpfr_set_emin ((mp_exp_t) gfc_real_kinds[kind].min_exponent -
3837                 mpfr_get_prec(result->value.real) + 1);
3838   mpfr_set_emax ((mp_exp_t) gfc_real_kinds[kind].max_exponent - 1);
3839   mpfr_check_range (result->value.real, 0, GMP_RNDU);
3840
3841   if (mpfr_sgn (s->value.real) > 0)
3842     {
3843       mpfr_nextabove (result->value.real);
3844       mpfr_subnormalize (result->value.real, 0, GMP_RNDU);
3845     }
3846   else
3847     {
3848       mpfr_nextbelow (result->value.real);
3849       mpfr_subnormalize (result->value.real, 0, GMP_RNDD);
3850     }
3851
3852   mpfr_set_emin (emin);
3853   mpfr_set_emax (emax);
3854
3855   /* Only NaN can occur. Do not use range check as it gives an
3856      error for denormal numbers.  */
3857   if (mpfr_nan_p (result->value.real) && gfc_option.flag_range_check)
3858     {
3859       gfc_error ("Result of NEAREST is NaN at %L", &result->where);
3860       gfc_free_expr (result);
3861       return &gfc_bad_expr;
3862     }
3863
3864   return result;
3865 }
3866
3867
3868 static gfc_expr *
3869 simplify_nint (const char *name, gfc_expr *e, gfc_expr *k)
3870 {
3871   gfc_expr *itrunc, *result;
3872   int kind;
3873
3874   kind = get_kind (BT_INTEGER, k, name, gfc_default_integer_kind);
3875   if (kind == -1)
3876     return &gfc_bad_expr;
3877
3878   if (e->expr_type != EXPR_CONSTANT)
3879     return NULL;
3880
3881   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
3882
3883   itrunc = gfc_copy_expr (e);
3884
3885   mpfr_round (itrunc->value.real, e->value.real);
3886
3887   gfc_mpfr_to_mpz (result->value.integer, itrunc->value.real, &e->where);
3888
3889   gfc_free_expr (itrunc);
3890
3891   return range_check (result, name);
3892 }
3893
3894
3895 gfc_expr *
3896 gfc_simplify_new_line (gfc_expr *e)
3897 {
3898   gfc_expr *result;
3899
3900   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
3901   result->value.character.string = gfc_get_wide_string (2);
3902   result->value.character.length = 1;
3903   result->value.character.string[0] = '\n';
3904   result->value.character.string[1] = '\0';     /* For debugger */
3905   return result;
3906 }
3907
3908
3909 gfc_expr *
3910 gfc_simplify_nint (gfc_expr *e, gfc_expr *k)
3911 {
3912   return simplify_nint ("NINT", e, k);
3913 }
3914
3915
3916 gfc_expr *
3917 gfc_simplify_idnint (gfc_expr *e)
3918 {
3919   return simplify_nint ("IDNINT", e, NULL);
3920 }
3921
3922
3923 gfc_expr *
3924 gfc_simplify_not (gfc_expr *e)
3925 {
3926   gfc_expr *result;
3927
3928   if (e->expr_type != EXPR_CONSTANT)
3929     return NULL;
3930
3931   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
3932
3933   mpz_com (result->value.integer, e->value.integer);
3934
3935   return range_check (result, "NOT");
3936 }
3937
3938
3939 gfc_expr *
3940 gfc_simplify_null (gfc_expr *mold)
3941 {
3942   gfc_expr *result;
3943
3944   if (mold == NULL)
3945     {
3946       result = gfc_get_expr ();
3947       result->ts.type = BT_UNKNOWN;
3948     }
3949   else
3950     result = gfc_copy_expr (mold);
3951   result->expr_type = EXPR_NULL;
3952
3953   return result;
3954 }
3955
3956
3957 gfc_expr *
3958 gfc_simplify_or (gfc_expr *x, gfc_expr *y)
3959 {
3960   gfc_expr *result;
3961   int kind;
3962
3963   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
3964     return NULL;
3965
3966   kind = x->ts.kind > y->ts.kind ? x->ts.kind : y->ts.kind;
3967   if (x->ts.type == BT_INTEGER)
3968     {
3969       result = gfc_constant_result (BT_INTEGER, kind, &x->where);
3970       mpz_ior (result->value.integer, x->value.integer, y->value.integer);
3971       return range_check (result, "OR");
3972     }
3973   else /* BT_LOGICAL */
3974     {
3975       result = gfc_constant_result (BT_LOGICAL, kind, &x->where);
3976       result->value.logical = x->value.logical || y->value.logical;
3977       return result;
3978     }
3979 }
3980
3981
3982 gfc_expr *
3983 gfc_simplify_pack (gfc_expr *array, gfc_expr *mask, gfc_expr *vector)
3984 {
3985   gfc_expr *result;
3986   gfc_constructor *array_ctor, *mask_ctor, *vector_ctor;
3987
3988   if (!is_constant_array_expr(array)
3989       || !is_constant_array_expr(vector)
3990       || (!gfc_is_constant_expr (mask)
3991           && !is_constant_array_expr(mask)))
3992     return NULL;
3993
3994   result = gfc_start_constructor (array->ts.type, 
3995                                   array->ts.kind,
3996                                   &array->where);
3997
3998   array_ctor = array->value.constructor;
3999   vector_ctor = vector ? vector->value.constructor : NULL;
4000
4001   if (mask->expr_type == EXPR_CONSTANT
4002       && mask->value.logical)
4003     {
4004       /* Copy all elements of ARRAY to RESULT.  */
4005       while (array_ctor)
4006         {
4007           gfc_append_constructor (result, 
4008                                   gfc_copy_expr (array_ctor->expr));
4009
4010           ADVANCE (array_ctor, 1);
4011           ADVANCE (vector_ctor, 1);
4012         }
4013     }
4014   else if (mask->expr_type == EXPR_ARRAY)
4015     {
4016       /* Copy only those elements of ARRAY to RESULT whose 
4017          MASK equals .TRUE..  */
4018       mask_ctor = mask->value.constructor;
4019       while (mask_ctor)
4020         {
4021           if (mask_ctor->expr->value.logical)
4022             {
4023               gfc_append_constructor (result, 
4024                                       gfc_copy_expr (array_ctor->expr)); 
4025               ADVANCE (vector_ctor, 1);
4026             }
4027
4028           ADVANCE (array_ctor, 1);
4029           ADVANCE (mask_ctor, 1);
4030         }
4031     }
4032
4033   /* Append any left-over elements from VECTOR to RESULT.  */
4034   while (vector_ctor)
4035     {
4036       gfc_append_constructor (result, 
4037                               gfc_copy_expr (vector_ctor->expr));
4038       ADVANCE (vector_ctor, 1);
4039     }
4040
4041   result->shape = gfc_get_shape (1);
4042   gfc_array_size (result, &result->shape[0]);
4043
4044   if (array->ts.type == BT_CHARACTER)
4045     result->ts.cl = array->ts.cl;
4046
4047   return result;
4048 }
4049
4050
4051 gfc_expr *
4052 gfc_simplify_precision (gfc_expr *e)
4053 {
4054   gfc_expr *result;
4055   int i;
4056
4057   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
4058
4059   result = gfc_int_expr (gfc_real_kinds[i].precision);
4060   result->where = e->where;
4061
4062   return result;
4063 }
4064
4065
4066 gfc_expr *
4067 gfc_simplify_product (gfc_expr *array, gfc_expr *dim, gfc_expr *mask)
4068 {
4069   gfc_expr *result;
4070
4071   if (!is_constant_array_expr (array)
4072       || !gfc_is_constant_expr (dim))
4073     return NULL;
4074
4075   if (mask
4076       && !is_constant_array_expr (mask)
4077       && mask->expr_type != EXPR_CONSTANT)
4078     return NULL;
4079
4080   result = transformational_result (array, dim, array->ts.type,
4081                                     array->ts.kind, &array->where);
4082   init_result_expr (result, 1, NULL);
4083
4084   return !dim || array->rank == 1 ?
4085     simplify_transformation_to_scalar (result, array, mask, gfc_multiply) :
4086     simplify_transformation_to_array (result, array, dim, mask, gfc_multiply);
4087 }
4088
4089
4090 gfc_expr *
4091 gfc_simplify_radix (gfc_expr *e)
4092 {
4093   gfc_expr *result;
4094   int i;
4095
4096   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
4097   switch (e->ts.type)
4098     {
4099     case BT_INTEGER:
4100       i = gfc_integer_kinds[i].radix;
4101       break;
4102
4103     case BT_REAL:
4104       i = gfc_real_kinds[i].radix;
4105       break;
4106
4107     default:
4108       gcc_unreachable ();
4109     }
4110
4111   result = gfc_int_expr (i);
4112   result->where = e->where;
4113
4114   return result;
4115 }
4116
4117
4118 gfc_expr *
4119 gfc_simplify_range (gfc_expr *e)
4120 {
4121   gfc_expr *result;
4122   int i;
4123   long j;
4124
4125   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
4126
4127   switch (e->ts.type)
4128     {
4129     case BT_INTEGER:
4130       j = gfc_integer_kinds[i].range;
4131       break;
4132
4133     case BT_REAL:
4134     case BT_COMPLEX:
4135       j = gfc_real_kinds[i].range;
4136       break;
4137
4138     default:
4139       gcc_unreachable ();
4140     }
4141
4142   result = gfc_int_expr (j);
4143   result->where = e->where;
4144
4145   return result;
4146 }
4147
4148
4149 gfc_expr *
4150 gfc_simplify_real (gfc_expr *e, gfc_expr *k)
4151 {
4152   gfc_expr *result = NULL;
4153   int kind;
4154
4155   if (e->ts.type == BT_COMPLEX)
4156     kind = get_kind (BT_REAL, k, "REAL", e->ts.kind);
4157   else
4158     kind = get_kind (BT_REAL, k, "REAL", gfc_default_real_kind);
4159
4160   if (kind == -1)
4161     return &gfc_bad_expr;
4162
4163   if (e->expr_type != EXPR_CONSTANT)
4164     return NULL;
4165
4166   switch (e->ts.type)
4167     {
4168     case BT_INTEGER:
4169       if (!e->is_boz)
4170         result = gfc_int2real (e, kind);
4171       break;
4172
4173     case BT_REAL:
4174       result = gfc_real2real (e, kind);
4175       break;
4176
4177     case BT_COMPLEX:
4178       result = gfc_complex2real (e, kind);
4179       break;
4180
4181     default:
4182       gfc_internal_error ("bad type in REAL");
4183       /* Not reached */
4184     }
4185
4186   if (e->ts.type == BT_INTEGER && e->is_boz)
4187     {
4188       gfc_typespec ts;
4189       gfc_clear_ts (&ts);
4190       ts.type = BT_REAL;
4191       ts.kind = kind;
4192       result = gfc_copy_expr (e);
4193       if (!gfc_convert_boz (result, &ts))
4194         {
4195           gfc_free_expr (result);
4196           return &gfc_bad_expr;
4197         }
4198     }
4199
4200   return range_check (result, "REAL");
4201 }
4202
4203
4204 gfc_expr *
4205 gfc_simplify_realpart (gfc_expr *e)
4206 {
4207   gfc_expr *result;
4208
4209   if (e->expr_type != EXPR_CONSTANT)
4210     return NULL;
4211
4212   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
4213 #ifdef HAVE_mpc
4214   mpc_real (result->value.real, e->value.complex, GFC_RND_MODE);
4215 #else
4216   mpfr_set (result->value.real, e->value.complex.r, GFC_RND_MODE);
4217 #endif
4218
4219   return range_check (result, "REALPART");
4220 }
4221
4222 gfc_expr *
4223 gfc_simplify_repeat (gfc_expr *e, gfc_expr *n)
4224 {
4225   gfc_expr *result;
4226   int i, j, len, ncop, nlen;
4227   mpz_t ncopies;
4228   bool have_length = false;
4229
4230   /* If NCOPIES isn't a constant, there's nothing we can do.  */
4231   if (n->expr_type != EXPR_CONSTANT)
4232     return NULL;
4233
4234   /* If NCOPIES is negative, it's an error.  */
4235   if (mpz_sgn (n->value.integer) < 0)
4236     {
4237       gfc_error ("Argument NCOPIES of REPEAT intrinsic is negative at %L",
4238                  &n->where);
4239       return &gfc_bad_expr;
4240     }
4241
4242   /* If we don't know the character length, we can do no more.  */
4243   if (e->ts.cl && e->ts.cl->length
4244         && e->ts.cl->length->expr_type == EXPR_CONSTANT)
4245     {
4246       len = mpz_get_si (e->ts.cl->length->value.integer);
4247       have_length = true;
4248     }
4249   else if (e->expr_type == EXPR_CONSTANT
4250              && (e->ts.cl == NULL || e->ts.cl->length == NULL))
4251     {
4252       len = e->value.character.length;
4253     }
4254   else
4255     return NULL;
4256
4257   /* If the source length is 0, any value of NCOPIES is valid
4258      and everything behaves as if NCOPIES == 0.  */
4259   mpz_init (ncopies);
4260   if (len == 0)
4261     mpz_set_ui (ncopies, 0);
4262   else
4263     mpz_set (ncopies, n->value.integer);
4264
4265   /* Check that NCOPIES isn't too large.  */
4266   if (len)
4267     {
4268       mpz_t max, mlen;
4269       int i;
4270
4271       /* Compute the maximum value allowed for NCOPIES: huge(cl) / len.  */
4272       mpz_init (max);
4273       i = gfc_validate_kind (BT_INTEGER, gfc_charlen_int_kind, false);
4274
4275       if (have_length)
4276         {
4277           mpz_tdiv_q (max, gfc_integer_kinds[i].huge,
4278                       e->ts.cl->length->value.integer);
4279         }
4280       else
4281         {
4282           mpz_init_set_si (mlen, len);
4283           mpz_tdiv_q (max, gfc_integer_kinds[i].huge, mlen);
4284           mpz_clear (mlen);
4285         }
4286
4287       /* The check itself.  */
4288       if (mpz_cmp (ncopies, max) > 0)
4289         {
4290           mpz_clear (max);
4291           mpz_clear (ncopies);
4292           gfc_error ("Argument NCOPIES of REPEAT intrinsic is too large at %L",
4293                      &n->where);
4294           return &gfc_bad_expr;
4295         }
4296
4297       mpz_clear (max);
4298     }
4299   mpz_clear (ncopies);
4300
4301   /* For further simplification, we need the character string to be
4302      constant.  */
4303   if (e->expr_type != EXPR_CONSTANT)
4304     return NULL;
4305
4306   if (len || 
4307       (e->ts.cl->length && 
4308        mpz_sgn (e->ts.cl->length->value.integer)) != 0)
4309     {
4310       const char *res = gfc_extract_int (n, &ncop);
4311       gcc_assert (res == NULL);
4312     }
4313   else
4314     ncop = 0;
4315
4316   len = e->value.character.length;
4317   nlen = ncop * len;
4318
4319   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
4320
4321   if (ncop == 0)
4322     {
4323       result->value.character.string = gfc_get_wide_string (1);
4324       result->value.character.length = 0;
4325       result->value.character.string[0] = '\0';
4326       return result;
4327     }
4328
4329   result->value.character.length = nlen;
4330   result->value.character.string = gfc_get_wide_string (nlen + 1);
4331
4332   for (i = 0; i < ncop; i++)
4333     for (j = 0; j < len; j++)
4334       result->value.character.string[j+i*len]= e->value.character.string[j];
4335
4336   result->value.character.string[nlen] = '\0';  /* For debugger */
4337   return result;
4338 }
4339
4340
4341 /* This one is a bear, but mainly has to do with shuffling elements.  */
4342
4343 gfc_expr *
4344 gfc_simplify_reshape (gfc_expr *source, gfc_expr *shape_exp,
4345                       gfc_expr *pad, gfc_expr *order_exp)
4346 {
4347   int order[GFC_MAX_DIMENSIONS], shape[GFC_MAX_DIMENSIONS];
4348   int i, rank, npad, x[GFC_MAX_DIMENSIONS];
4349   gfc_constructor *head, *tail;
4350   mpz_t index, size;
4351   unsigned long j;
4352   size_t nsource;
4353   gfc_expr *e;
4354
4355   /* Check that argument expression types are OK.  */
4356   if (!is_constant_array_expr (source)
4357       || !is_constant_array_expr (shape_exp)
4358       || !is_constant_array_expr (pad)
4359       || !is_constant_array_expr (order_exp))
4360     return NULL;
4361
4362   /* Proceed with simplification, unpacking the array.  */
4363
4364   mpz_init (index);
4365   rank = 0;
4366   head = tail = NULL;
4367
4368   for (;;)
4369     {
4370       e = gfc_get_array_element (shape_exp, rank);
4371       if (e == NULL)
4372         break;
4373
4374       gfc_extract_int (e, &shape[rank]);
4375
4376       gcc_assert (rank >= 0 && rank < GFC_MAX_DIMENSIONS);
4377       gcc_assert (shape[rank] >= 0);
4378
4379       gfc_free_expr (e);
4380       rank++;
4381     }
4382
4383   gcc_assert (rank > 0);
4384
4385   /* Now unpack the order array if present.  */
4386   if (order_exp == NULL)
4387     {
4388       for (i = 0; i < rank; i++)
4389         order[i] = i;
4390     }
4391   else
4392     {
4393       for (i = 0; i < rank; i++)
4394         x[i] = 0;
4395
4396       for (i = 0; i < rank; i++)
4397         {
4398           e = gfc_get_array_element (order_exp, i);
4399           gcc_assert (e);
4400
4401           gfc_extract_int (e, &order[i]);
4402           gfc_free_expr (e);
4403
4404           gcc_assert (order[i] >= 1 && order[i] <= rank);
4405           order[i]--;
4406           gcc_assert (x[order[i]] == 0);
4407           x[order[i]] = 1;
4408         }
4409     }
4410
4411   /* Count the elements in the source and padding arrays.  */
4412
4413   npad = 0;
4414   if (pad != NULL)
4415     {
4416       gfc_array_size (pad, &size);
4417       npad = mpz_get_ui (size);
4418       mpz_clear (size);
4419     }
4420
4421   gfc_array_size (source, &size);
4422   nsource = mpz_get_ui (size);
4423   mpz_clear (size);
4424
4425   /* If it weren't for that pesky permutation we could just loop
4426      through the source and round out any shortage with pad elements.
4427      But no, someone just had to have the compiler do something the
4428      user should be doing.  */
4429
4430   for (i = 0; i < rank; i++)
4431     x[i] = 0;
4432
4433   while (nsource > 0 || npad > 0)
4434     {
4435       /* Figure out which element to extract.  */
4436       mpz_set_ui (index, 0);
4437
4438       for (i = rank - 1; i >= 0; i--)
4439         {
4440           mpz_add_ui (index, index, x[order[i]]);
4441           if (i != 0)
4442             mpz_mul_ui (index, index, shape[order[i - 1]]);
4443         }
4444
4445       if (mpz_cmp_ui (index, INT_MAX) > 0)
4446         gfc_internal_error ("Reshaped array too large at %C");
4447
4448       j = mpz_get_ui (index);
4449
4450       if (j < nsource)
4451         e = gfc_get_array_element (source, j);
4452       else
4453         {
4454           gcc_assert (npad > 0);
4455
4456           j = j - nsource;
4457           j = j % npad;
4458           e = gfc_get_array_element (pad, j);
4459         }
4460       gcc_assert (e);
4461
4462       if (head == NULL)
4463         head = tail = gfc_get_constructor ();
4464       else
4465         {
4466           tail->next = gfc_get_constructor ();
4467           tail = tail->next;
4468         }
4469
4470       tail->where = e->where;
4471       tail->expr = e;
4472
4473       /* Calculate the next element.  */
4474       i = 0;
4475
4476 inc:
4477       if (++x[i] < shape[i])
4478         continue;
4479       x[i++] = 0;
4480       if (i < rank)
4481         goto inc;
4482
4483       break;
4484     }
4485
4486   mpz_clear (index);
4487
4488   e = gfc_get_expr ();
4489   e->where = source->where;
4490   e->expr_type = EXPR_ARRAY;
4491   e->value.constructor = head;
4492   e->shape = gfc_get_shape (rank);
4493
4494   for (i = 0; i < rank; i++)
4495     mpz_init_set_ui (e->shape[i], shape[i]);
4496
4497   e->ts = source->ts;
4498   e->rank = rank;
4499
4500   return e;
4501 }
4502
4503
4504 gfc_expr *
4505 gfc_simplify_rrspacing (gfc_expr *x)
4506 {
4507   gfc_expr *result;
4508   int i;
4509   long int e, p;
4510
4511   if (x->expr_type != EXPR_CONSTANT)
4512     return NULL;
4513
4514   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
4515
4516   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
4517
4518   mpfr_abs (result->value.real, x->value.real, GFC_RND_MODE);
4519
4520   /* Special case x = -0 and 0.  */
4521   if (mpfr_sgn (result->value.real) == 0)
4522     {
4523       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
4524       return result;
4525     }
4526
4527   /* | x * 2**(-e) | * 2**p.  */
4528   e = - (long int) mpfr_get_exp (x->value.real);
4529   mpfr_mul_2si (result->value.real, result->value.real, e, GFC_RND_MODE);
4530
4531   p = (long int) gfc_real_kinds[i].digits;
4532   mpfr_mul_2si (result->value.real, result->value.real, p, GFC_RND_MODE);
4533
4534   return range_check (result, "RRSPACING");
4535 }
4536
4537
4538 gfc_expr *
4539 gfc_simplify_scale (gfc_expr *x, gfc_expr *i)
4540 {
4541   int k, neg_flag, power, exp_range;
4542   mpfr_t scale, radix;
4543   gfc_expr *result;
4544
4545   if (x->expr_type != EXPR_CONSTANT || i->expr_type != EXPR_CONSTANT)
4546     return NULL;
4547
4548   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
4549
4550   if (mpfr_sgn (x->value.real) == 0)
4551     {
4552       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
4553       return result;
4554     }
4555
4556   k = gfc_validate_kind (BT_REAL, x->ts.kind, false);
4557
4558   exp_range = gfc_real_kinds[k].max_exponent - gfc_real_kinds[k].min_exponent;
4559
4560   /* This check filters out values of i that would overflow an int.  */
4561   if (mpz_cmp_si (i->value.integer, exp_range + 2) > 0
4562       || mpz_cmp_si (i->value.integer, -exp_range - 2) < 0)
4563     {
4564       gfc_error ("Result of SCALE overflows its kind at %L", &result->where);
4565       gfc_free_expr (result);
4566       return &gfc_bad_expr;
4567     }
4568
4569   /* Compute scale = radix ** power.  */
4570   power = mpz_get_si (i->value.integer);
4571
4572   if (power >= 0)
4573     neg_flag = 0;
4574   else
4575     {
4576       neg_flag = 1;
4577       power = -power;
4578     }
4579
4580   gfc_set_model_kind (x->ts.kind);
4581   mpfr_init (scale);
4582   mpfr_init (radix);
4583   mpfr_set_ui (radix, gfc_real_kinds[k].radix, GFC_RND_MODE);
4584   mpfr_pow_ui (scale, radix, power, GFC_RND_MODE);
4585
4586   if (neg_flag)
4587     mpfr_div (result->value.real, x->value.real, scale, GFC_RND_MODE);
4588   else
4589     mpfr_mul (result->value.real, x->value.real, scale, GFC_RND_MODE);
4590
4591   mpfr_clears (scale, radix, NULL);
4592
4593   return range_check (result, "SCALE");
4594 }
4595
4596
4597 /* Variants of strspn and strcspn that operate on wide characters.  */
4598
4599 static size_t
4600 wide_strspn (const gfc_char_t *s1, const gfc_char_t *s2)
4601 {
4602   size_t i = 0;
4603   const gfc_char_t *c;
4604
4605   while (s1[i])
4606     {
4607       for (c = s2; *c; c++)
4608         {
4609           if (s1[i] == *c)
4610             break;
4611         }
4612       if (*c == '\0')
4613         break;
4614       i++;
4615     }
4616
4617   return i;
4618 }
4619
4620 static size_t
4621 wide_strcspn (const gfc_char_t *s1, const gfc_char_t *s2)
4622 {
4623   size_t i = 0;
4624   const gfc_char_t *c;
4625
4626   while (s1[i])
4627     {
4628       for (c = s2; *c; c++)
4629         {
4630           if (s1[i] == *c)
4631             break;
4632         }
4633       if (*c)
4634         break;
4635       i++;
4636     }
4637
4638   return i;
4639 }
4640
4641
4642 gfc_expr *
4643 gfc_simplify_scan (gfc_expr *e, gfc_expr *c, gfc_expr *b, gfc_expr *kind)
4644 {
4645   gfc_expr *result;
4646   int back;
4647   size_t i;
4648   size_t indx, len, lenc;
4649   int k = get_kind (BT_INTEGER, kind, "SCAN", gfc_default_integer_kind);
4650
4651   if (k == -1)
4652     return &gfc_bad_expr;
4653
4654   if (e->expr_type != EXPR_CONSTANT || c->expr_type != EXPR_CONSTANT)
4655     return NULL;
4656
4657   if (b != NULL && b->value.logical != 0)
4658     back = 1;
4659   else
4660     back = 0;
4661
4662   result = gfc_constant_result (BT_INTEGER, k, &e->where);
4663
4664   len = e->value.character.length;
4665   lenc = c->value.character.length;
4666
4667   if (len == 0 || lenc == 0)
4668     {
4669       indx = 0;
4670     }
4671   else
4672     {
4673       if (back == 0)
4674         {
4675           indx = wide_strcspn (e->value.character.string,
4676                                c->value.character.string) + 1;
4677           if (indx > len)
4678             indx = 0;
4679         }
4680       else
4681         {
4682           i = 0;
4683           for (indx = len; indx > 0; indx--)
4684             {
4685               for (i = 0; i < lenc; i++)
4686                 {
4687                   if (c->value.character.string[i]
4688                       == e->value.character.string[indx - 1])
4689                     break;
4690                 }
4691               if (i < lenc)
4692                 break;
4693             }
4694         }
4695     }
4696   mpz_set_ui (result->value.integer, indx);
4697   return range_check (result, "SCAN");
4698 }
4699
4700
4701 gfc_expr *
4702 gfc_simplify_selected_char_kind (gfc_expr *e)
4703 {
4704   int kind;
4705   gfc_expr *result;
4706
4707   if (e->expr_type != EXPR_CONSTANT)
4708     return NULL;
4709
4710   if (gfc_compare_with_Cstring (e, "ascii", false) == 0
4711       || gfc_compare_with_Cstring (e, "default", false) == 0)
4712     kind = 1;
4713   else if (gfc_compare_with_Cstring (e, "iso_10646", false) == 0)
4714     kind = 4;
4715   else
4716     kind = -1;
4717
4718   result = gfc_int_expr (kind);
4719   result->where = e->where;
4720
4721   return result;
4722 }
4723
4724
4725 gfc_expr *
4726 gfc_simplify_selected_int_kind (gfc_expr *e)
4727 {
4728   int i, kind, range;
4729   gfc_expr *result;
4730
4731   if (e->expr_type != EXPR_CONSTANT || gfc_extract_int (e, &range) != NULL)
4732     return NULL;
4733
4734   kind = INT_MAX;
4735
4736   for (i = 0; gfc_integer_kinds[i].kind != 0; i++)
4737     if (gfc_integer_kinds[i].range >= range
4738         && gfc_integer_kinds[i].kind < kind)
4739       kind = gfc_integer_kinds[i].kind;
4740
4741   if (kind == INT_MAX)
4742     kind = -1;
4743
4744   result = gfc_int_expr (kind);
4745   result->where = e->where;
4746
4747   return result;
4748 }
4749
4750
4751 gfc_expr *
4752 gfc_simplify_selected_real_kind (gfc_expr *p, gfc_expr *q)
4753 {
4754   int range, precision, i, kind, found_precision, found_range;
4755   gfc_expr *result;
4756
4757   if (p == NULL)
4758     precision = 0;
4759   else
4760     {
4761       if (p->expr_type != EXPR_CONSTANT
4762           || gfc_extract_int (p, &precision) != NULL)
4763         return NULL;
4764     }
4765
4766   if (q == NULL)
4767     range = 0;
4768   else
4769     {
4770       if (q->expr_type != EXPR_CONSTANT
4771           || gfc_extract_int (q, &range) != NULL)
4772         return NULL;
4773     }
4774
4775   kind = INT_MAX;
4776   found_precision = 0;
4777   found_range = 0;
4778
4779   for (i = 0; gfc_real_kinds[i].kind != 0; i++)
4780     {
4781       if (gfc_real_kinds[i].precision >= precision)
4782         found_precision = 1;
4783
4784       if (gfc_real_kinds[i].range >= range)
4785         found_range = 1;
4786
4787       if (gfc_real_kinds[i].precision >= precision
4788           && gfc_real_kinds[i].range >= range && gfc_real_kinds[i].kind < kind)
4789         kind = gfc_real_kinds[i].kind;
4790     }
4791
4792   if (kind == INT_MAX)
4793     {
4794       kind = 0;
4795
4796       if (!found_precision)
4797         kind = -1;
4798       if (!found_range)
4799         kind -= 2;
4800     }
4801
4802   result = gfc_int_expr (kind);
4803   result->where = (p != NULL) ? p->where : q->where;
4804
4805   return result;
4806 }
4807
4808
4809 gfc_expr *
4810 gfc_simplify_set_exponent (gfc_expr *x, gfc_expr *i)
4811 {
4812   gfc_expr *result;
4813   mpfr_t exp, absv, log2, pow2, frac;
4814   unsigned long exp2;
4815
4816   if (x->expr_type != EXPR_CONSTANT || i->expr_type != EXPR_CONSTANT)
4817     return NULL;
4818
4819   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
4820
4821   if (mpfr_sgn (x->value.real) == 0)
4822     {
4823       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
4824       return result;
4825     }
4826
4827   gfc_set_model_kind (x->ts.kind);
4828   mpfr_init (absv);
4829   mpfr_init (log2);
4830   mpfr_init (exp);
4831   mpfr_init (pow2);
4832   mpfr_init (frac);
4833
4834   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
4835   mpfr_log2 (log2, absv, GFC_RND_MODE);
4836
4837   mpfr_trunc (log2, log2);
4838   mpfr_add_ui (exp, log2, 1, GFC_RND_MODE);
4839
4840   /* Old exponent value, and fraction.  */
4841   mpfr_ui_pow (pow2, 2, exp, GFC_RND_MODE);
4842
4843   mpfr_div (frac, absv, pow2, GFC_RND_MODE);
4844
4845   /* New exponent.  */
4846   exp2 = (unsigned long) mpz_get_d (i->value.integer);
4847   mpfr_mul_2exp (result->value.real, frac, exp2, GFC_RND_MODE);
4848
4849   mpfr_clears (absv, log2, pow2, frac, NULL);
4850
4851   return range_check (result, "SET_EXPONENT");
4852 }
4853
4854
4855 gfc_expr *
4856 gfc_simplify_shape (gfc_expr *source)
4857 {
4858   mpz_t shape[GFC_MAX_DIMENSIONS];
4859   gfc_expr *result, *e, *f;
4860   gfc_array_ref *ar;
4861   int n;
4862   gfc_try t;
4863
4864   if (source->rank == 0)
4865     return gfc_start_constructor (BT_INTEGER, gfc_default_integer_kind,
4866                                   &source->where);
4867
4868   if (source->expr_type != EXPR_VARIABLE)
4869     return NULL;
4870
4871   result = gfc_start_constructor (BT_INTEGER, gfc_default_integer_kind,
4872                                   &source->where);
4873
4874   ar = gfc_find_array_ref (source);
4875
4876   t = gfc_array_ref_shape (ar, shape);
4877
4878   for (n = 0; n < source->rank; n++)
4879     {
4880       e = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
4881                                &source->where);
4882
4883       if (t == SUCCESS)
4884         {
4885           mpz_set (e->value.integer, shape[n]);
4886           mpz_clear (shape[n]);
4887         }
4888       else
4889         {
4890           mpz_set_ui (e->value.integer, n + 1);
4891
4892           f = gfc_simplify_size (source, e, NULL);
4893           gfc_free_expr (e);
4894           if (f == NULL)
4895             {
4896               gfc_free_expr (result);
4897               return NULL;
4898             }
4899           else
4900             {
4901               e = f;
4902             }
4903         }
4904
4905       gfc_append_constructor (result, e);
4906     }
4907
4908   return result;
4909 }
4910
4911
4912 gfc_expr *
4913 gfc_simplify_size (gfc_expr *array, gfc_expr *dim, gfc_expr *kind)
4914 {
4915   mpz_t size;
4916   gfc_expr *result;
4917   int d;
4918   int k = get_kind (BT_INTEGER, kind, "SIZE", gfc_default_integer_kind);
4919
4920   if (k == -1)
4921     return &gfc_bad_expr;
4922
4923   if (dim == NULL)
4924     {
4925       if (gfc_array_size (array, &size) == FAILURE)
4926         return NULL;
4927     }
4928   else
4929     {
4930       if (dim->expr_type != EXPR_CONSTANT)
4931         return NULL;
4932
4933       d = mpz_get_ui (dim->value.integer) - 1;
4934       if (gfc_array_dimen_size (array, d, &size) == FAILURE)
4935         return NULL;
4936     }
4937
4938   result = gfc_constant_result (BT_INTEGER, k, &array->where);
4939   mpz_set (result->value.integer, size);
4940   return result;
4941 }
4942
4943
4944 gfc_expr *
4945 gfc_simplify_sign (gfc_expr *x, gfc_expr *y)
4946 {
4947   gfc_expr *result;
4948
4949   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
4950     return NULL;
4951
4952   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
4953
4954   switch (x->ts.type)
4955     {
4956     case BT_INTEGER:
4957       mpz_abs (result->value.integer, x->value.integer);
4958       if (mpz_sgn (y->value.integer) < 0)
4959         mpz_neg (result->value.integer, result->value.integer);
4960       break;
4961
4962     case BT_REAL:
4963       if (gfc_option.flag_sign_zero)
4964         mpfr_copysign (result->value.real, x->value.real, y->value.real,
4965                        GFC_RND_MODE);
4966       else
4967         mpfr_setsign (result->value.real, x->value.real,
4968                       mpfr_sgn (y->value.real) < 0 ? 1 : 0, GFC_RND_MODE);
4969       break;
4970
4971     default:
4972       gfc_internal_error ("Bad type in gfc_simplify_sign");
4973     }
4974
4975   return result;
4976 }
4977
4978
4979 gfc_expr *
4980 gfc_simplify_sin (gfc_expr *x)
4981 {
4982   gfc_expr *result;
4983
4984   if (x->expr_type != EXPR_CONSTANT)
4985     return NULL;
4986
4987   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
4988
4989   switch (x->ts.type)
4990     {
4991     case BT_REAL:
4992       mpfr_sin (result->value.real, x->value.real, GFC_RND_MODE);
4993       break;
4994
4995     case BT_COMPLEX:
4996       gfc_set_model (x->value.real);
4997 #ifdef HAVE_mpc
4998       mpc_sin (result->value.complex, x->value.complex, GFC_MPC_RND_MODE);
4999 #else
5000     {
5001       mpfr_t xp, xq;
5002       mpfr_init (xp);
5003       mpfr_init (xq);
5004
5005       mpfr_sin  (xp, x->value.complex.r, GFC_RND_MODE);
5006       mpfr_cosh (xq, x->value.complex.i, GFC_RND_MODE);
5007       mpfr_mul (result->value.complex.r, xp, xq, GFC_RND_MODE);
5008
5009       mpfr_cos  (xp, x->value.complex.r, GFC_RND_MODE);
5010       mpfr_sinh (xq, x->value.complex.i, GFC_RND_MODE);
5011       mpfr_mul (result->value.complex.i, xp, xq, GFC_RND_MODE);
5012
5013       mpfr_clears (xp, xq, NULL);
5014     }
5015 #endif
5016       break;
5017
5018     default:
5019       gfc_internal_error ("in gfc_simplify_sin(): Bad type");
5020     }
5021
5022   return range_check (result, "SIN");
5023 }
5024
5025
5026 gfc_expr *
5027 gfc_simplify_sinh (gfc_expr *x)
5028 {
5029   gfc_expr *result;
5030
5031   if (x->expr_type != EXPR_CONSTANT)
5032     return NULL;
5033
5034   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
5035
5036   mpfr_sinh (result->value.real, x->value.real, GFC_RND_MODE);
5037
5038   return range_check (result, "SINH");
5039 }
5040
5041
5042 /* The argument is always a double precision real that is converted to
5043    single precision.  TODO: Rounding!  */
5044
5045 gfc_expr *
5046 gfc_simplify_sngl (gfc_expr *a)
5047 {
5048   gfc_expr *result;
5049
5050   if (a->expr_type != EXPR_CONSTANT)
5051     return NULL;
5052
5053   result = gfc_real2real (a, gfc_default_real_kind);
5054   return range_check (result, "SNGL");
5055 }
5056
5057
5058 gfc_expr *
5059 gfc_simplify_spacing (gfc_expr *x)
5060 {
5061   gfc_expr *result;
5062   int i;
5063   long int en, ep;
5064
5065   if (x->expr_type != EXPR_CONSTANT)
5066     return NULL;
5067
5068   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
5069
5070   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
5071
5072   /* Special case x = 0 and -0.  */
5073   mpfr_abs (result->value.real, x->value.real, GFC_RND_MODE);
5074   if (mpfr_sgn (result->value.real) == 0)
5075     {
5076       mpfr_set (result->value.real, gfc_real_kinds[i].tiny, GFC_RND_MODE);
5077       return result;
5078     }
5079
5080   /* In the Fortran 95 standard, the result is b**(e - p) where b, e, and p
5081      are the radix, exponent of x, and precision.  This excludes the 
5082      possibility of subnormal numbers.  Fortran 2003 states the result is
5083      b**max(e - p, emin - 1).  */
5084
5085   ep = (long int) mpfr_get_exp (x->value.real) - gfc_real_kinds[i].digits;
5086   en = (long int) gfc_real_kinds[i].min_exponent - 1;
5087   en = en > ep ? en : ep;
5088
5089   mpfr_set_ui (result->value.real, 1, GFC_RND_MODE);
5090   mpfr_mul_2si (result->value.real, result->value.real, en, GFC_RND_MODE);
5091
5092   return range_check (result, "SPACING");
5093 }
5094
5095
5096 gfc_expr *
5097 gfc_simplify_spread (gfc_expr *source, gfc_expr *dim_expr, gfc_expr *ncopies_expr)
5098 {
5099   gfc_expr *result = 0L;
5100   int i, j, dim, ncopies;
5101   mpz_t size;
5102
5103   if ((!gfc_is_constant_expr (source)
5104        && !is_constant_array_expr (source))
5105       || !gfc_is_constant_expr (dim_expr)
5106       || !gfc_is_constant_expr (ncopies_expr))
5107     return NULL;
5108
5109   gcc_assert (dim_expr->ts.type == BT_INTEGER);
5110   gfc_extract_int (dim_expr, &dim);
5111   dim -= 1;   /* zero-base DIM */
5112
5113   gcc_assert (ncopies_expr->ts.type == BT_INTEGER);
5114   gfc_extract_int (ncopies_expr, &ncopies);
5115   ncopies = MAX (ncopies, 0);
5116
5117   /* Do not allow the array size to exceed the limit for an array
5118      constructor.  */
5119   if (source->expr_type == EXPR_ARRAY)
5120     {
5121       if (gfc_array_size (source, &size) == FAILURE)
5122         gfc_internal_error ("Failure getting length of a constant array.");
5123     }
5124   else
5125     mpz_init_set_ui (size, 1);
5126
5127   if (mpz_get_si (size)*ncopies > gfc_option.flag_max_array_constructor)
5128     return NULL;
5129
5130   if (source->expr_type == EXPR_CONSTANT)
5131     {
5132       gcc_assert (dim == 0);
5133
5134       result = gfc_start_constructor (source->ts.type,
5135                                       source->ts.kind,
5136                                       &source->where);
5137       result->rank = 1;
5138       result->shape = gfc_get_shape (result->rank);
5139       mpz_init_set_si (result->shape[0], ncopies);
5140
5141       for (i = 0; i < ncopies; ++i)
5142         gfc_append_constructor (result, gfc_copy_expr (source));
5143     }
5144   else if (source->expr_type == EXPR_ARRAY)
5145     {
5146       int result_size, rstride[GFC_MAX_DIMENSIONS], extent[GFC_MAX_DIMENSIONS];
5147       gfc_constructor *ctor, *source_ctor, *result_ctor;
5148
5149       gcc_assert (source->rank < GFC_MAX_DIMENSIONS);
5150       gcc_assert (dim >= 0 && dim <= source->rank);
5151
5152       result = gfc_start_constructor (source->ts.type,
5153                                       source->ts.kind,
5154                                       &source->where);
5155       result->rank = source->rank + 1;
5156       result->shape = gfc_get_shape (result->rank);
5157
5158       result_size = 1;
5159       for (i = 0, j = 0; i < result->rank; ++i)
5160         {
5161           if (i != dim)
5162             mpz_init_set (result->shape[i], source->shape[j++]);
5163           else
5164             mpz_init_set_si (result->shape[i], ncopies);
5165
5166           extent[i] = mpz_get_si (result->shape[i]);
5167           rstride[i] = (i == 0) ? 1 : rstride[i-1] * extent[i-1];
5168           result_size *= extent[i];
5169         }
5170
5171       for (i = 0; i < result_size; ++i)
5172         gfc_append_constructor (result, NULL);
5173
5174       source_ctor = source->value.constructor;
5175       result_ctor = result->value.constructor;
5176       while (source_ctor)
5177         {
5178           ctor = result_ctor;
5179
5180           for (i = 0; i < ncopies; ++i)
5181           {
5182             ctor->expr = gfc_copy_expr (source_ctor->expr);
5183             ADVANCE (ctor, rstride[dim]);
5184           }
5185
5186           ADVANCE (result_ctor, (dim == 0 ? ncopies : 1));
5187           ADVANCE (source_ctor, 1);
5188         }
5189     }
5190   else
5191     /* FIXME: Returning here avoids a regression in array_simplify_1.f90.
5192        Replace NULL with gcc_unreachable() after implementing
5193        gfc_simplify_cshift(). */
5194     return NULL;
5195
5196   if (source->ts.type == BT_CHARACTER)
5197     result->ts.cl = source->ts.cl;
5198
5199   return result;
5200 }
5201
5202
5203 gfc_expr *
5204 gfc_simplify_sqrt (gfc_expr *e)
5205 {
5206   gfc_expr *result;
5207
5208   if (e->expr_type != EXPR_CONSTANT)
5209     return NULL;
5210
5211   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
5212
5213   switch (e->ts.type)
5214     {
5215     case BT_REAL:
5216       if (mpfr_cmp_si (e->value.real, 0) < 0)
5217         goto negative_arg;
5218       mpfr_sqrt (result->value.real, e->value.real, GFC_RND_MODE);
5219
5220       break;
5221
5222     case BT_COMPLEX:
5223       gfc_set_model (e->value.real);
5224 #ifdef HAVE_mpc
5225       mpc_sqrt (result->value.complex, e->value.complex, GFC_MPC_RND_MODE);
5226 #else
5227     {
5228       /* Formula taken from Numerical Recipes to avoid over- and
5229          underflow.  */
5230
5231       mpfr_t ac, ad, s, t, w;
5232       mpfr_init (ac);
5233       mpfr_init (ad);
5234       mpfr_init (s);
5235       mpfr_init (t);
5236       mpfr_init (w);
5237
5238       if (mpfr_cmp_ui (e->value.complex.r, 0) == 0
5239           && mpfr_cmp_ui (e->value.complex.i, 0) == 0)
5240         {
5241           mpfr_set_ui (result->value.complex.r, 0, GFC_RND_MODE);
5242           mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
5243           break;
5244         }
5245
5246       mpfr_abs (ac, e->value.complex.r, GFC_RND_MODE);
5247       mpfr_abs (ad, e->value.complex.i, GFC_RND_MODE);
5248
5249       if (mpfr_cmp (ac, ad) >= 0)
5250         {
5251           mpfr_div (t, e->value.complex.i, e->value.complex.r, GFC_RND_MODE);
5252           mpfr_mul (t, t, t, GFC_RND_MODE);
5253           mpfr_add_ui (t, t, 1, GFC_RND_MODE);
5254           mpfr_sqrt (t, t, GFC_RND_MODE);
5255           mpfr_add_ui (t, t, 1, GFC_RND_MODE);
5256           mpfr_div_ui (t, t, 2, GFC_RND_MODE);
5257           mpfr_sqrt (t, t, GFC_RND_MODE);
5258           mpfr_sqrt (s, ac, GFC_RND_MODE);
5259           mpfr_mul (w, s, t, GFC_RND_MODE);
5260         }
5261       else
5262         {
5263           mpfr_div (s, e->value.complex.r, e->value.complex.i, GFC_RND_MODE);
5264           mpfr_mul (t, s, s, GFC_RND_MODE);
5265           mpfr_add_ui (t, t, 1, GFC_RND_MODE);
5266           mpfr_sqrt (t, t, GFC_RND_MODE);
5267           mpfr_abs (s, s, GFC_RND_MODE);
5268           mpfr_add (t, t, s, GFC_RND_MODE);
5269           mpfr_div_ui (t, t, 2, GFC_RND_MODE);
5270           mpfr_sqrt (t, t, GFC_RND_MODE);
5271           mpfr_sqrt (s, ad, GFC_RND_MODE);
5272           mpfr_mul (w, s, t, GFC_RND_MODE);
5273         }
5274
5275       if (mpfr_cmp_ui (w, 0) != 0 && mpfr_cmp_ui (e->value.complex.r, 0) >= 0)
5276         {
5277           mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
5278           mpfr_div (result->value.complex.i, e->value.complex.i, t, GFC_RND_MODE);
5279           mpfr_set (result->value.complex.r, w, GFC_RND_MODE);
5280         }
5281       else if (mpfr_cmp_ui (w, 0) != 0
5282                && mpfr_cmp_ui (e->value.complex.r, 0) < 0
5283                && mpfr_cmp_ui (e->value.complex.i, 0) >= 0)
5284         {
5285           mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
5286           mpfr_div (result->value.complex.r, e->value.complex.i, t, GFC_RND_MODE);
5287           mpfr_set (result->value.complex.i, w, GFC_RND_MODE);
5288         }
5289       else if (mpfr_cmp_ui (w, 0) != 0
5290                && mpfr_cmp_ui (e->value.complex.r, 0) < 0
5291                && mpfr_cmp_ui (e->value.complex.i, 0) < 0)
5292         {
5293           mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
5294           mpfr_div (result->value.complex.r, ad, t, GFC_RND_MODE);
5295           mpfr_neg (w, w, GFC_RND_MODE);
5296           mpfr_set (result->value.complex.i, w, GFC_RND_MODE);
5297         }
5298       else
5299         gfc_internal_error ("invalid complex argument of SQRT at %L",
5300                             &e->where);
5301
5302       mpfr_clears (s, t, ac, ad, w, NULL);
5303     }
5304 #endif
5305       break;
5306
5307     default:
5308       gfc_internal_error ("invalid argument of SQRT at %L", &e->where);
5309     }
5310
5311   return range_check (result, "SQRT");
5312
5313 negative_arg:
5314   gfc_free_expr (result);
5315   gfc_error ("Argument of SQRT at %L has a negative value", &e->where);
5316   return &gfc_bad_expr;
5317 }
5318
5319
5320 gfc_expr *
5321 gfc_simplify_sum (gfc_expr *array, gfc_expr *dim, gfc_expr *mask)
5322 {
5323   gfc_expr *result;
5324
5325   if (!is_constant_array_expr (array)
5326       || !gfc_is_constant_expr (dim))
5327     return NULL;
5328
5329   if (mask
5330       && !is_constant_array_expr (mask)
5331       && mask->expr_type != EXPR_CONSTANT)
5332     return NULL;
5333
5334   result = transformational_result (array, dim, array->ts.type,
5335                                     array->ts.kind, &array->where);
5336   init_result_expr (result, 0, NULL);
5337
5338   return !dim || array->rank == 1 ?
5339     simplify_transformation_to_scalar (result, array, mask, gfc_add) :
5340     simplify_transformation_to_array (result, array, dim, mask, gfc_add);
5341 }
5342
5343
5344 gfc_expr *
5345 gfc_simplify_tan (gfc_expr *x)
5346 {
5347   int i;
5348   gfc_expr *result;
5349
5350   if (x->expr_type != EXPR_CONSTANT)
5351     return NULL;
5352
5353   i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
5354
5355   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
5356
5357   mpfr_tan (result->value.real, x->value.real, GFC_RND_MODE);
5358
5359   return range_check (result, "TAN");
5360 }
5361
5362
5363 gfc_expr *
5364 gfc_simplify_tanh (gfc_expr *x)
5365 {
5366   gfc_expr *result;
5367
5368   if (x->expr_type != EXPR_CONSTANT)
5369     return NULL;
5370
5371   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
5372
5373   mpfr_tanh (result->value.real, x->value.real, GFC_RND_MODE);
5374
5375   return range_check (result, "TANH");
5376
5377 }
5378
5379
5380 gfc_expr *
5381 gfc_simplify_tiny (gfc_expr *e)
5382 {
5383   gfc_expr *result;
5384   int i;
5385
5386   i = gfc_validate_kind (BT_REAL, e->ts.kind, false);
5387
5388   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
5389   mpfr_set (result->value.real, gfc_real_kinds[i].tiny, GFC_RND_MODE);
5390
5391   return result;
5392 }
5393
5394
5395 gfc_expr *
5396 gfc_simplify_trailz (gfc_expr *e)
5397 {
5398   gfc_expr *result;
5399   unsigned long tz, bs;
5400   int i;
5401
5402   if (e->expr_type != EXPR_CONSTANT)
5403     return NULL;
5404
5405   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
5406   bs = gfc_integer_kinds[i].bit_size;
5407   tz = mpz_scan1 (e->value.integer, 0);
5408
5409   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind, &e->where);
5410   mpz_set_ui (result->value.integer, MIN (tz, bs));
5411
5412   return result;
5413 }
5414
5415
5416 gfc_expr *
5417 gfc_simplify_transfer (gfc_expr *source, gfc_expr *mold, gfc_expr *size)
5418 {
5419   gfc_expr *result;
5420   gfc_expr *mold_element;
5421   size_t source_size;
5422   size_t result_size;
5423   size_t result_elt_size;
5424   size_t buffer_size;
5425   mpz_t tmp;
5426   unsigned char *buffer;
5427
5428   if (!gfc_is_constant_expr (source)
5429         || (gfc_init_expr && !gfc_is_constant_expr (mold))
5430         || !gfc_is_constant_expr (size))
5431     return NULL;
5432
5433   if (source->expr_type == EXPR_FUNCTION)
5434     return NULL;
5435
5436   /* Calculate the size of the source.  */
5437   if (source->expr_type == EXPR_ARRAY
5438       && gfc_array_size (source, &tmp) == FAILURE)
5439     gfc_internal_error ("Failure getting length of a constant array.");
5440
5441   source_size = gfc_target_expr_size (source);
5442
5443   /* Create an empty new expression with the appropriate characteristics.  */
5444   result = gfc_constant_result (mold->ts.type, mold->ts.kind,
5445                                 &source->where);
5446   result->ts = mold->ts;
5447
5448   mold_element = mold->expr_type == EXPR_ARRAY
5449                  ? mold->value.constructor->expr
5450                  : mold;
5451
5452   /* Set result character length, if needed.  Note that this needs to be
5453      set even for array expressions, in order to pass this information into 
5454      gfc_target_interpret_expr.  */
5455   if (result->ts.type == BT_CHARACTER && gfc_is_constant_expr (mold_element))
5456     result->value.character.length = mold_element->value.character.length;
5457   
5458   /* Set the number of elements in the result, and determine its size.  */
5459   result_elt_size = gfc_target_expr_size (mold_element);
5460   if (result_elt_size == 0)
5461     {
5462       gfc_free_expr (result);
5463       return NULL;
5464     }
5465
5466   if (mold->expr_type == EXPR_ARRAY || mold->rank || size)
5467     {
5468       int result_length;
5469
5470       result->expr_type = EXPR_ARRAY;
5471       result->rank = 1;
5472
5473       if (size)
5474         result_length = (size_t)mpz_get_ui (size->value.integer);
5475       else
5476         {
5477           result_length = source_size / result_elt_size;
5478           if (result_length * result_elt_size < source_size)
5479             result_length += 1;
5480         }
5481
5482       result->shape = gfc_get_shape (1);
5483       mpz_init_set_ui (result->shape[0], result_length);
5484
5485       result_size = result_length * result_elt_size;
5486     }
5487   else
5488     {
5489       result->rank = 0;
5490       result_size = result_elt_size;
5491     }
5492
5493   if (gfc_option.warn_surprising && source_size < result_size)
5494     gfc_warning("Intrinsic TRANSFER at %L has partly undefined result: "
5495                 "source size %ld < result size %ld", &source->where,
5496                 (long) source_size, (long) result_size);
5497
5498   /* Allocate the buffer to store the binary version of the source.  */
5499   buffer_size = MAX (source_size, result_size);
5500   buffer = (unsigned char*)alloca (buffer_size);
5501   memset (buffer, 0, buffer_size);
5502
5503   /* Now write source to the buffer.  */
5504   gfc_target_encode_expr (source, buffer, buffer_size);
5505
5506   /* And read the buffer back into the new expression.  */
5507   gfc_target_interpret_expr (buffer, buffer_size, result);
5508
5509   return result;
5510 }
5511
5512
5513 gfc_expr *
5514 gfc_simplify_transpose (gfc_expr *matrix)
5515 {
5516   int i, matrix_rows;
5517   gfc_expr *result;
5518   gfc_constructor *matrix_ctor;
5519
5520   if (!is_constant_array_expr (matrix))
5521     return NULL;
5522
5523   gcc_assert (matrix->rank == 2);
5524
5525   result = gfc_start_constructor (matrix->ts.type, matrix->ts.kind, &matrix->where);
5526   result->rank = 2;
5527   result->shape = gfc_get_shape (result->rank);
5528   mpz_set (result->shape[0], matrix->shape[1]);
5529   mpz_set (result->shape[1], matrix->shape[0]);
5530
5531   if (matrix->ts.type == BT_CHARACTER)
5532     result->ts.cl = matrix->ts.cl;
5533
5534   matrix_rows = mpz_get_si (matrix->shape[0]);
5535   matrix_ctor = matrix->value.constructor;
5536   for (i = 0; i < matrix_rows; ++i)
5537     {
5538       gfc_constructor *column_ctor = matrix_ctor;
5539       while (column_ctor)
5540         {
5541           gfc_append_constructor (result, 
5542                                   gfc_copy_expr (column_ctor->expr));
5543
5544           ADVANCE (column_ctor, matrix_rows);
5545         }
5546
5547       ADVANCE (matrix_ctor, 1);
5548     }
5549
5550   return result;
5551 }
5552
5553
5554 gfc_expr *
5555 gfc_simplify_trim (gfc_expr *e)
5556 {
5557   gfc_expr *result;
5558   int count, i, len, lentrim;
5559
5560   if (e->expr_type != EXPR_CONSTANT)
5561     return NULL;
5562
5563   len = e->value.character.length;
5564
5565   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
5566
5567   for (count = 0, i = 1; i <= len; ++i)
5568     {
5569       if (e->value.character.string[len - i] == ' ')
5570         count++;
5571       else
5572         break;
5573     }
5574
5575   lentrim = len - count;
5576
5577   result->value.character.length = lentrim;
5578   result->value.character.string = gfc_get_wide_string (lentrim + 1);
5579
5580   for (i = 0; i < lentrim; i++)
5581     result->value.character.string[i] = e->value.character.string[i];
5582
5583   result->value.character.string[lentrim] = '\0';       /* For debugger */
5584
5585   return result;
5586 }
5587
5588
5589 gfc_expr *
5590 gfc_simplify_ubound (gfc_expr *array, gfc_expr *dim, gfc_expr *kind)
5591 {
5592   return simplify_bound (array, dim, kind, 1);
5593 }
5594
5595
5596 gfc_expr *
5597 gfc_simplify_unpack (gfc_expr *vector, gfc_expr *mask, gfc_expr *field)
5598 {
5599   gfc_expr *result, *e;
5600   gfc_constructor *vector_ctor, *mask_ctor, *field_ctor;
5601
5602   if (!is_constant_array_expr (vector)
5603       || !is_constant_array_expr (mask)
5604       || (!gfc_is_constant_expr (field)
5605           && !is_constant_array_expr(field)))
5606     return NULL;
5607
5608   result = gfc_start_constructor (vector->ts.type,
5609                                   vector->ts.kind,
5610                                   &vector->where);
5611   result->rank = mask->rank;
5612   result->shape = gfc_copy_shape (mask->shape, mask->rank);
5613
5614   if (vector->ts.type == BT_CHARACTER)
5615     result->ts.cl = vector->ts.cl;
5616
5617   vector_ctor = vector->value.constructor;
5618   mask_ctor = mask->value.constructor;
5619   field_ctor = field->expr_type == EXPR_ARRAY ? field->value.constructor : NULL;
5620
5621   while (mask_ctor)
5622     {
5623       if (mask_ctor->expr->value.logical)
5624         {
5625           gcc_assert (vector_ctor);
5626           e = gfc_copy_expr (vector_ctor->expr);
5627           ADVANCE (vector_ctor, 1);
5628         }
5629       else if (field->expr_type == EXPR_ARRAY)
5630         e = gfc_copy_expr (field_ctor->expr);
5631       else
5632         e = gfc_copy_expr (field);
5633
5634       gfc_append_constructor (result, e);
5635
5636       ADVANCE (mask_ctor, 1);
5637       ADVANCE (field_ctor, 1);
5638     }
5639
5640   return result;
5641 }
5642
5643
5644 gfc_expr *
5645 gfc_simplify_verify (gfc_expr *s, gfc_expr *set, gfc_expr *b, gfc_expr *kind)
5646 {
5647   gfc_expr *result;
5648   int back;
5649   size_t index, len, lenset;
5650   size_t i;
5651   int k = get_kind (BT_INTEGER, kind, "VERIFY", gfc_default_integer_kind);
5652
5653   if (k == -1)
5654     return &gfc_bad_expr;
5655
5656   if (s->expr_type != EXPR_CONSTANT || set->expr_type != EXPR_CONSTANT)
5657     return NULL;
5658
5659   if (b != NULL && b->value.logical != 0)
5660     back = 1;
5661   else
5662     back = 0;
5663
5664   result = gfc_constant_result (BT_INTEGER, k, &s->where);
5665
5666   len = s->value.character.length;
5667   lenset = set->value.character.length;
5668
5669   if (len == 0)
5670     {
5671       mpz_set_ui (result->value.integer, 0);
5672       return result;
5673     }
5674
5675   if (back == 0)
5676     {
5677       if (lenset == 0)
5678         {
5679           mpz_set_ui (result->value.integer, 1);
5680           return result;
5681         }
5682
5683       index = wide_strspn (s->value.character.string,
5684                            set->value.character.string) + 1;
5685       if (index > len)
5686         index = 0;
5687
5688     }
5689   else
5690     {
5691       if (lenset == 0)
5692         {
5693           mpz_set_ui (result->value.integer, len);
5694           return result;
5695         }
5696       for (index = len; index > 0; index --)
5697         {
5698           for (i = 0; i < lenset; i++)
5699             {
5700               if (s->value.character.string[index - 1]
5701                   == set->value.character.string[i])
5702                 break;
5703             }
5704           if (i == lenset)
5705             break;
5706         }
5707     }
5708
5709   mpz_set_ui (result->value.integer, index);
5710   return result;
5711 }
5712
5713
5714 gfc_expr *
5715 gfc_simplify_xor (gfc_expr *x, gfc_expr *y)
5716 {
5717   gfc_expr *result;
5718   int kind;
5719
5720   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
5721     return NULL;
5722
5723   kind = x->ts.kind > y->ts.kind ? x->ts.kind : y->ts.kind;
5724   if (x->ts.type == BT_INTEGER)
5725     {
5726       result = gfc_constant_result (BT_INTEGER, kind, &x->where);
5727       mpz_xor (result->value.integer, x->value.integer, y->value.integer);
5728       return range_check (result, "XOR");
5729     }
5730   else /* BT_LOGICAL */
5731     {
5732       result = gfc_constant_result (BT_LOGICAL, kind, &x->where);
5733       result->value.logical = (x->value.logical && !y->value.logical)
5734                               || (!x->value.logical && y->value.logical);
5735       return result;
5736     }
5737
5738 }
5739
5740
5741 /****************** Constant simplification *****************/
5742
5743 /* Master function to convert one constant to another.  While this is
5744    used as a simplification function, it requires the destination type
5745    and kind information which is supplied by a special case in
5746    do_simplify().  */
5747
5748 gfc_expr *
5749 gfc_convert_constant (gfc_expr *e, bt type, int kind)
5750 {
5751   gfc_expr *g, *result, *(*f) (gfc_expr *, int);
5752   gfc_constructor *head, *c, *tail = NULL;
5753
5754   switch (e->ts.type)
5755     {
5756     case BT_INTEGER:
5757       switch (type)
5758         {
5759         case BT_INTEGER:
5760           f = gfc_int2int;
5761           break;
5762         case BT_REAL:
5763           f = gfc_int2real;
5764           break;
5765         case BT_COMPLEX:
5766           f = gfc_int2complex;
5767           break;
5768         case BT_LOGICAL:
5769           f = gfc_int2log;
5770           break;
5771         default:
5772           goto oops;
5773         }
5774       break;
5775
5776     case BT_REAL:
5777       switch (type)
5778         {
5779         case BT_INTEGER:
5780           f = gfc_real2int;
5781           break;
5782         case BT_REAL:
5783           f = gfc_real2real;
5784           break;
5785         case BT_COMPLEX:
5786           f = gfc_real2complex;
5787           break;
5788         default:
5789           goto oops;
5790         }
5791       break;
5792
5793     case BT_COMPLEX:
5794       switch (type)
5795         {
5796         case BT_INTEGER:
5797           f = gfc_complex2int;
5798           break;
5799         case BT_REAL:
5800           f = gfc_complex2real;
5801           break;
5802         case BT_COMPLEX:
5803           f = gfc_complex2complex;
5804           break;
5805
5806         default:
5807           goto oops;
5808         }
5809       break;
5810
5811     case BT_LOGICAL:
5812       switch (type)
5813         {
5814         case BT_INTEGER:
5815           f = gfc_log2int;
5816           break;
5817         case BT_LOGICAL:
5818           f = gfc_log2log;
5819           break;
5820         default:
5821           goto oops;
5822         }
5823       break;
5824
5825     case BT_HOLLERITH:
5826       switch (type)
5827         {
5828         case BT_INTEGER:
5829           f = gfc_hollerith2int;
5830           break;
5831
5832         case BT_REAL:
5833           f = gfc_hollerith2real;
5834           break;
5835
5836         case BT_COMPLEX:
5837           f = gfc_hollerith2complex;
5838           break;
5839
5840         case BT_CHARACTER:
5841           f = gfc_hollerith2character;
5842           break;
5843
5844         case BT_LOGICAL:
5845           f = gfc_hollerith2logical;
5846           break;
5847
5848         default:
5849           goto oops;
5850         }
5851       break;
5852
5853     default:
5854     oops:
5855       gfc_internal_error ("gfc_convert_constant(): Unexpected type");
5856     }
5857
5858   result = NULL;
5859
5860   switch (e->expr_type)
5861     {
5862     case EXPR_CONSTANT:
5863       result = f (e, kind);
5864       if (result == NULL)
5865         return &gfc_bad_expr;
5866       break;
5867
5868     case EXPR_ARRAY:
5869       if (!gfc_is_constant_expr (e))
5870         break;
5871
5872       head = NULL;
5873
5874       for (c = e->value.constructor; c; c = c->next)
5875         {
5876           if (head == NULL)
5877             head = tail = gfc_get_constructor ();
5878           else
5879             {
5880               tail->next = gfc_get_constructor ();
5881               tail = tail->next;
5882             }
5883
5884           tail->where = c->where;
5885
5886           if (c->iterator == NULL)
5887             tail->expr = f (c->expr, kind);
5888           else
5889             {
5890               g = gfc_convert_constant (c->expr, type, kind);
5891               if (g == &gfc_bad_expr)
5892                 return g;
5893               tail->expr = g;
5894             }
5895
5896           if (tail->expr == NULL)
5897             {
5898               gfc_free_constructor (head);
5899               return NULL;
5900             }
5901         }
5902
5903       result = gfc_get_expr ();
5904       result->ts.type = type;
5905       result->ts.kind = kind;
5906       result->expr_type = EXPR_ARRAY;
5907       result->value.constructor = head;
5908       result->shape = gfc_copy_shape (e->shape, e->rank);
5909       result->where = e->where;
5910       result->rank = e->rank;
5911       break;
5912
5913     default:
5914       break;
5915     }
5916
5917   return result;
5918 }
5919
5920
5921 /* Function for converting character constants.  */
5922 gfc_expr *
5923 gfc_convert_char_constant (gfc_expr *e, bt type ATTRIBUTE_UNUSED, int kind)
5924 {
5925   gfc_expr *result;
5926   int i;
5927
5928   if (!gfc_is_constant_expr (e))
5929     return NULL;
5930
5931   if (e->expr_type == EXPR_CONSTANT)
5932     {
5933       /* Simple case of a scalar.  */
5934       result = gfc_constant_result (BT_CHARACTER, kind, &e->where);
5935       if (result == NULL)
5936         return &gfc_bad_expr;
5937
5938       result->value.character.length = e->value.character.length;
5939       result->value.character.string
5940         = gfc_get_wide_string (e->value.character.length + 1);
5941       memcpy (result->value.character.string, e->value.character.string,
5942               (e->value.character.length + 1) * sizeof (gfc_char_t));
5943
5944       /* Check we only have values representable in the destination kind.  */
5945       for (i = 0; i < result->value.character.length; i++)
5946         if (!gfc_check_character_range (result->value.character.string[i],
5947                                         kind))
5948           {
5949             gfc_error ("Character '%s' in string at %L cannot be converted "
5950                        "into character kind %d",
5951                        gfc_print_wide_char (result->value.character.string[i]),
5952                        &e->where, kind);
5953             return &gfc_bad_expr;
5954           }
5955
5956       return result;
5957     }
5958   else if (e->expr_type == EXPR_ARRAY)
5959     {
5960       /* For an array constructor, we convert each constructor element.  */
5961       gfc_constructor *head = NULL, *tail = NULL, *c;
5962
5963       for (c = e->value.constructor; c; c = c->next)
5964         {
5965           if (head == NULL)
5966             head = tail = gfc_get_constructor ();
5967           else
5968             {
5969               tail->next = gfc_get_constructor ();
5970               tail = tail->next;
5971             }
5972
5973           tail->where = c->where;
5974           tail->expr = gfc_convert_char_constant (c->expr, type, kind);
5975           if (tail->expr == &gfc_bad_expr)
5976             {
5977               tail->expr = NULL;
5978               return &gfc_bad_expr;
5979             }
5980
5981           if (tail->expr == NULL)
5982             {
5983               gfc_free_constructor (head);
5984               return NULL;
5985             }
5986         }
5987
5988       result = gfc_get_expr ();
5989       result->ts.type = type;
5990       result->ts.kind = kind;
5991       result->expr_type = EXPR_ARRAY;
5992       result->value.constructor = head;
5993       result->shape = gfc_copy_shape (e->shape, e->rank);
5994       result->where = e->where;
5995       result->rank = e->rank;
5996       result->ts.cl = e->ts.cl;
5997
5998       return result;
5999     }
6000   else
6001     return NULL;
6002 }