OSDN Git Service

* simplify.c (gfc_simplify_anint): Use mpfr_round()
[pf3gnuchains/gcc-fork.git] / gcc / fortran / simplify.c
1 /* Simplify intrinsic functions at compile-time.
2    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005 Free Software Foundation,
3    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 2, 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 COPYING.  If not, write to the Free
20 Software Foundation, 59 Temple Place - Suite 330, Boston, MA
21 02111-1307, USA.  */
22
23 #include "config.h"
24 #include "system.h"
25 #include "flags.h"
26 #include "gfortran.h"
27 #include "arith.h"
28 #include "intrinsic.h"
29
30 gfc_expr gfc_bad_expr;
31
32
33 /* Note that 'simplification' is not just transforming expressions.
34    For functions that are not simplified at compile time, range
35    checking is done if possible.
36
37    The return convention is that each simplification function returns:
38
39      A new expression node corresponding to the simplified arguments.
40      The original arguments are destroyed by the caller, and must not
41      be a part of the new expression.
42
43      NULL pointer indicating that no simplification was possible and
44      the original expression should remain intact.  If the
45      simplification function sets the type and/or the function name
46      via the pointer gfc_simple_expression, then this type is
47      retained.
48
49      An expression pointer to gfc_bad_expr (a static placeholder)
50      indicating that some error has prevented simplification.  For
51      example, sqrt(-1.0).  The error is generated within the function
52      and should be propagated upwards
53
54    By the time a simplification function gets control, it has been
55    decided that the function call is really supposed to be the
56    intrinsic.  No type checking is strictly necessary, since only
57    valid types will be passed on.  On the other hand, a simplification
58    subroutine may have to look at the type of an argument as part of
59    its processing.
60
61    Array arguments are never passed to these subroutines.
62
63    The functions in this file don't have much comment with them, but
64    everything is reasonably straight-forward.  The Standard, chapter 13
65    is the best comment you'll find for this file anyway.  */
66
67 /* Static table for converting non-ascii character sets to ascii.
68    The xascii_table[] is the inverse table.  */
69
70 static int ascii_table[256] = {
71   '\0', '\0', '\0', '\0', '\0', '\0', '\0', '\0',
72   '\b', '\t', '\n', '\v', '\0', '\r', '\0', '\0',
73   '\0', '\0', '\0', '\0', '\0', '\0', '\0', '\0',
74   '\0', '\0', '\0', '\0', '\0', '\0', '\0', '\0',
75   ' ', '!', '\'', '#', '$', '%', '&', '\'',
76   '(', ')', '*', '+', ',', '-', '.', '/',
77   '0', '1', '2', '3', '4', '5', '6', '7',
78   '8', '9', ':', ';', '<', '=', '>', '?',
79   '@', 'A', 'B', 'C', 'D', 'E', 'F', 'G',
80   'H', 'I', 'J', 'K', 'L', 'M', 'N', 'O',
81   'P', 'Q', 'R', 'S', 'T', 'U', 'V', 'W',
82   'X', 'Y', 'Z', '[', '\\', ']', '^', '_',
83   '`', 'a', 'b', 'c', 'd', 'e', 'f', 'g',
84   'h', 'i', 'j', 'k', 'l', 'm', 'n', 'o',
85   'p', 'q', 'r', 's', 't', 'u', 'v', 'w',
86   'x', 'y', 'z', '{', '|', '}', '~', '\?'
87 };
88
89 static int xascii_table[256];
90
91
92 /* Range checks an expression node.  If all goes well, returns the
93    node, otherwise returns &gfc_bad_expr and frees the node.  */
94
95 static gfc_expr *
96 range_check (gfc_expr * result, const char *name)
97 {
98   if (gfc_range_check (result) == ARITH_OK)
99     return result;
100
101   gfc_error ("Result of %s overflows its kind at %L", name, &result->where);
102   gfc_free_expr (result);
103   return &gfc_bad_expr;
104 }
105
106
107 /* A helper function that gets an optional and possibly missing
108    kind parameter.  Returns the kind, -1 if something went wrong.  */
109
110 static int
111 get_kind (bt type, gfc_expr * k, const char *name, int default_kind)
112 {
113   int kind;
114
115   if (k == NULL)
116     return default_kind;
117
118   if (k->expr_type != EXPR_CONSTANT)
119     {
120       gfc_error ("KIND parameter of %s at %L must be an initialization "
121                  "expression", name, &k->where);
122
123       return -1;
124     }
125
126   if (gfc_extract_int (k, &kind) != NULL
127       || gfc_validate_kind (type, kind, true) < 0)
128     {
129
130       gfc_error ("Invalid KIND parameter of %s at %L", name, &k->where);
131       return -1;
132     }
133
134   return kind;
135 }
136
137
138 /* Checks if X, which is assumed to represent a two's complement
139    integer of binary width BITSIZE, has the signbit set.  If so, makes 
140    X the corresponding negative number.  */
141
142 static void
143 twos_complement (mpz_t x, int bitsize)
144 {
145   mpz_t mask;
146
147   if (mpz_tstbit (x, bitsize - 1) == 1)
148     {
149       mpz_init_set_ui(mask, 1);
150       mpz_mul_2exp(mask, mask, bitsize);
151       mpz_sub_ui(mask, mask, 1);
152
153       /* We negate the number by hand, zeroing the high bits, that is
154         make it the corresponding positive number, and then have it
155         negated by GMP, giving the correct representation of the
156         negative number.  */
157       mpz_com (x, x);
158       mpz_add_ui (x, x, 1);
159       mpz_and (x, x, mask);
160
161       mpz_neg (x, x);
162
163       mpz_clear (mask);
164     }
165 }
166
167
168 /********************** Simplification functions *****************************/
169
170 gfc_expr *
171 gfc_simplify_abs (gfc_expr * e)
172 {
173   gfc_expr *result;
174
175   if (e->expr_type != EXPR_CONSTANT)
176     return NULL;
177
178   switch (e->ts.type)
179     {
180     case BT_INTEGER:
181       result = gfc_constant_result (BT_INTEGER, e->ts.kind, &e->where);
182
183       mpz_abs (result->value.integer, e->value.integer);
184
185       result = range_check (result, "IABS");
186       break;
187
188     case BT_REAL:
189       result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
190
191       mpfr_abs (result->value.real, e->value.real, GFC_RND_MODE);
192
193       result = range_check (result, "ABS");
194       break;
195
196     case BT_COMPLEX:
197       result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
198
199       gfc_set_model_kind (e->ts.kind);
200
201       mpfr_hypot (result->value.real, e->value.complex.r, 
202                   e->value.complex.i, GFC_RND_MODE);
203       result = range_check (result, "CABS");
204       break;
205
206     default:
207       gfc_internal_error ("gfc_simplify_abs(): Bad type");
208     }
209
210   return result;
211 }
212
213
214 gfc_expr *
215 gfc_simplify_achar (gfc_expr * e)
216 {
217   gfc_expr *result;
218   int index;
219
220   if (e->expr_type != EXPR_CONSTANT)
221     return NULL;
222
223   /* We cannot assume that the native character set is ASCII in this
224      function.  */
225   if (gfc_extract_int (e, &index) != NULL || index < 0 || index > 127)
226     {
227       gfc_error ("Extended ASCII not implemented: argument of ACHAR at %L "
228                  "must be between 0 and 127", &e->where);
229       return &gfc_bad_expr;
230     }
231
232   result = gfc_constant_result (BT_CHARACTER, gfc_default_character_kind,
233                                 &e->where);
234
235   result->value.character.string = gfc_getmem (2);
236
237   result->value.character.length = 1;
238   result->value.character.string[0] = ascii_table[index];
239   result->value.character.string[1] = '\0';     /* For debugger */
240   return result;
241 }
242
243
244 gfc_expr *
245 gfc_simplify_acos (gfc_expr * x)
246 {
247   gfc_expr *result;
248
249   if (x->expr_type != EXPR_CONSTANT)
250     return NULL;
251
252   if (mpfr_cmp_si (x->value.real, 1) > 0 || mpfr_cmp_si (x->value.real, -1) < 0)
253     {
254       gfc_error ("Argument of ACOS at %L must be between -1 and 1",
255                  &x->where);
256       return &gfc_bad_expr;
257     }
258
259   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
260
261   mpfr_acos (result->value.real, x->value.real, GFC_RND_MODE);
262
263   return range_check (result, "ACOS");
264 }
265
266
267 gfc_expr *
268 gfc_simplify_adjustl (gfc_expr * e)
269 {
270   gfc_expr *result;
271   int count, i, len;
272   char ch;
273
274   if (e->expr_type != EXPR_CONSTANT)
275     return NULL;
276
277   len = e->value.character.length;
278
279   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
280
281   result->value.character.length = len;
282   result->value.character.string = gfc_getmem (len + 1);
283
284   for (count = 0, i = 0; i < len; ++i)
285     {
286       ch = e->value.character.string[i];
287       if (ch != ' ')
288         break;
289       ++count;
290     }
291
292   for (i = 0; i < len - count; ++i)
293     {
294       result->value.character.string[i] =
295         e->value.character.string[count + i];
296     }
297
298   for (i = len - count; i < len; ++i)
299     {
300       result->value.character.string[i] = ' ';
301     }
302
303   result->value.character.string[len] = '\0';   /* For debugger */
304
305   return result;
306 }
307
308
309 gfc_expr *
310 gfc_simplify_adjustr (gfc_expr * e)
311 {
312   gfc_expr *result;
313   int count, i, len;
314   char ch;
315
316   if (e->expr_type != EXPR_CONSTANT)
317     return NULL;
318
319   len = e->value.character.length;
320
321   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
322
323   result->value.character.length = len;
324   result->value.character.string = gfc_getmem (len + 1);
325
326   for (count = 0, i = len - 1; i >= 0; --i)
327     {
328       ch = e->value.character.string[i];
329       if (ch != ' ')
330         break;
331       ++count;
332     }
333
334   for (i = 0; i < count; ++i)
335     {
336       result->value.character.string[i] = ' ';
337     }
338
339   for (i = count; i < len; ++i)
340     {
341       result->value.character.string[i] =
342         e->value.character.string[i - count];
343     }
344
345   result->value.character.string[len] = '\0';   /* For debugger */
346
347   return result;
348 }
349
350
351 gfc_expr *
352 gfc_simplify_aimag (gfc_expr * e)
353 {
354   gfc_expr *result;
355
356   if (e->expr_type != EXPR_CONSTANT)
357     return NULL;
358
359   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
360   mpfr_set (result->value.real, e->value.complex.i, GFC_RND_MODE);
361
362   return range_check (result, "AIMAG");
363 }
364
365
366 gfc_expr *
367 gfc_simplify_aint (gfc_expr * e, gfc_expr * k)
368 {
369   gfc_expr *rtrunc, *result;
370   int kind;
371
372   kind = get_kind (BT_REAL, k, "AINT", e->ts.kind);
373   if (kind == -1)
374     return &gfc_bad_expr;
375
376   if (e->expr_type != EXPR_CONSTANT)
377     return NULL;
378
379   rtrunc = gfc_copy_expr (e);
380
381   mpfr_trunc (rtrunc->value.real, e->value.real);
382
383   result = gfc_real2real (rtrunc, kind);
384   gfc_free_expr (rtrunc);
385
386   return range_check (result, "AINT");
387 }
388
389
390 gfc_expr *
391 gfc_simplify_dint (gfc_expr * e)
392 {
393   gfc_expr *rtrunc, *result;
394
395   if (e->expr_type != EXPR_CONSTANT)
396     return NULL;
397
398   rtrunc = gfc_copy_expr (e);
399
400   mpfr_trunc (rtrunc->value.real, e->value.real);
401
402   result = gfc_real2real (rtrunc, gfc_default_double_kind);
403   gfc_free_expr (rtrunc);
404
405   return range_check (result, "DINT");
406 }
407
408
409 gfc_expr *
410 gfc_simplify_anint (gfc_expr * e, gfc_expr * k)
411 {
412   gfc_expr *result;
413   int kind;
414
415   kind = get_kind (BT_REAL, k, "ANINT", e->ts.kind);
416   if (kind == -1)
417     return &gfc_bad_expr;
418
419   if (e->expr_type != EXPR_CONSTANT)
420     return NULL;
421
422   result = gfc_constant_result (e->ts.type, kind, &e->where);
423
424   mpfr_round (result->value.real, e->value.real);
425
426   return range_check (result, "ANINT");
427 }
428
429
430 gfc_expr *
431 gfc_simplify_dnint (gfc_expr * e)
432 {
433   gfc_expr *result;
434
435   if (e->expr_type != EXPR_CONSTANT)
436     return NULL;
437
438   result = gfc_constant_result (BT_REAL, gfc_default_double_kind, &e->where);
439
440   mpfr_round (result->value.real, e->value.real);
441
442   return range_check (result, "DNINT");
443 }
444
445
446 gfc_expr *
447 gfc_simplify_asin (gfc_expr * x)
448 {
449   gfc_expr *result;
450
451   if (x->expr_type != EXPR_CONSTANT)
452     return NULL;
453
454   if (mpfr_cmp_si (x->value.real, 1) > 0 || mpfr_cmp_si (x->value.real, -1) < 0)
455     {
456       gfc_error ("Argument of ASIN at %L must be between -1 and 1",
457                  &x->where);
458       return &gfc_bad_expr;
459     }
460
461   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
462
463   mpfr_asin(result->value.real, x->value.real, GFC_RND_MODE);
464
465   return range_check (result, "ASIN");
466 }
467
468
469 gfc_expr *
470 gfc_simplify_atan (gfc_expr * x)
471 {
472   gfc_expr *result;
473
474   if (x->expr_type != EXPR_CONSTANT)
475     return NULL;
476
477   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
478
479   mpfr_atan(result->value.real, x->value.real, GFC_RND_MODE);
480
481   return range_check (result, "ATAN");
482
483 }
484
485
486 gfc_expr *
487 gfc_simplify_atan2 (gfc_expr * y, gfc_expr * x)
488 {
489   gfc_expr *result;
490
491   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
492     return NULL;
493
494   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
495
496   if (mpfr_sgn (y->value.real) == 0 && mpfr_sgn (x->value.real) == 0)
497     {
498       gfc_error
499         ("If first argument of ATAN2 %L is zero, then the second argument "
500           "must not be zero", &x->where);
501       gfc_free_expr (result);
502       return &gfc_bad_expr;
503     }
504
505   arctangent2 (y->value.real, x->value.real, result->value.real);
506
507   return range_check (result, "ATAN2");
508
509 }
510
511
512 gfc_expr *
513 gfc_simplify_bit_size (gfc_expr * e)
514 {
515   gfc_expr *result;
516   int i;
517
518   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
519   result = gfc_constant_result (BT_INTEGER, e->ts.kind, &e->where);
520   mpz_set_ui (result->value.integer, gfc_integer_kinds[i].bit_size);
521
522   return result;
523 }
524
525
526 gfc_expr *
527 gfc_simplify_btest (gfc_expr * e, gfc_expr * bit)
528 {
529   int b;
530
531   if (e->expr_type != EXPR_CONSTANT || bit->expr_type != EXPR_CONSTANT)
532     return NULL;
533
534   if (gfc_extract_int (bit, &b) != NULL || b < 0)
535     return gfc_logical_expr (0, &e->where);
536
537   return gfc_logical_expr (mpz_tstbit (e->value.integer, b), &e->where);
538 }
539
540
541 gfc_expr *
542 gfc_simplify_ceiling (gfc_expr * e, gfc_expr * k)
543 {
544   gfc_expr *ceil, *result;
545   int kind;
546
547   kind = get_kind (BT_INTEGER, k, "CEILING", gfc_default_integer_kind);
548   if (kind == -1)
549     return &gfc_bad_expr;
550
551   if (e->expr_type != EXPR_CONSTANT)
552     return NULL;
553
554   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
555
556   ceil = gfc_copy_expr (e);
557
558   mpfr_ceil (ceil->value.real, e->value.real);
559   gfc_mpfr_to_mpz(result->value.integer, ceil->value.real);
560
561   gfc_free_expr (ceil);
562
563   return range_check (result, "CEILING");
564 }
565
566
567 gfc_expr *
568 gfc_simplify_char (gfc_expr * e, gfc_expr * k)
569 {
570   gfc_expr *result;
571   int c, kind;
572
573   kind = get_kind (BT_CHARACTER, k, "CHAR", gfc_default_character_kind);
574   if (kind == -1)
575     return &gfc_bad_expr;
576
577   if (e->expr_type != EXPR_CONSTANT)
578     return NULL;
579
580   if (gfc_extract_int (e, &c) != NULL || c < 0 || c > 255)
581     {
582       gfc_error ("Bad character in CHAR function at %L", &e->where);
583       return &gfc_bad_expr;
584     }
585
586   result = gfc_constant_result (BT_CHARACTER, kind, &e->where);
587
588   result->value.character.length = 1;
589   result->value.character.string = gfc_getmem (2);
590
591   result->value.character.string[0] = c;
592   result->value.character.string[1] = '\0';     /* For debugger */
593
594   return result;
595 }
596
597
598 /* Common subroutine for simplifying CMPLX and DCMPLX.  */
599
600 static gfc_expr *
601 simplify_cmplx (const char *name, gfc_expr * x, gfc_expr * y, int kind)
602 {
603   gfc_expr *result;
604
605   result = gfc_constant_result (BT_COMPLEX, kind, &x->where);
606
607   mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
608
609   switch (x->ts.type)
610     {
611     case BT_INTEGER:
612       mpfr_set_z (result->value.complex.r, x->value.integer, GFC_RND_MODE);
613       break;
614
615     case BT_REAL:
616       mpfr_set (result->value.complex.r, x->value.real, GFC_RND_MODE);
617       break;
618
619     case BT_COMPLEX:
620       mpfr_set (result->value.complex.r, x->value.complex.r, GFC_RND_MODE);
621       mpfr_set (result->value.complex.i, x->value.complex.i, GFC_RND_MODE);
622       break;
623
624     default:
625       gfc_internal_error ("gfc_simplify_dcmplx(): Bad type (x)");
626     }
627
628   if (y != NULL)
629     {
630       switch (y->ts.type)
631         {
632         case BT_INTEGER:
633           mpfr_set_z (result->value.complex.i, y->value.integer, GFC_RND_MODE);
634           break;
635
636         case BT_REAL:
637           mpfr_set (result->value.complex.i, y->value.real, GFC_RND_MODE);
638           break;
639
640         default:
641           gfc_internal_error ("gfc_simplify_dcmplx(): Bad type (y)");
642         }
643     }
644
645   return range_check (result, name);
646 }
647
648
649 gfc_expr *
650 gfc_simplify_cmplx (gfc_expr * x, gfc_expr * y, gfc_expr * k)
651 {
652   int kind;
653
654   if (x->expr_type != EXPR_CONSTANT
655       || (y != NULL && y->expr_type != EXPR_CONSTANT))
656     return NULL;
657
658   kind = get_kind (BT_REAL, k, "CMPLX", gfc_default_real_kind);
659   if (kind == -1)
660     return &gfc_bad_expr;
661
662   return simplify_cmplx ("CMPLX", x, y, kind);
663 }
664
665
666 gfc_expr *
667 gfc_simplify_conjg (gfc_expr * e)
668 {
669   gfc_expr *result;
670
671   if (e->expr_type != EXPR_CONSTANT)
672     return NULL;
673
674   result = gfc_copy_expr (e);
675   mpfr_neg (result->value.complex.i, result->value.complex.i, GFC_RND_MODE);
676
677   return range_check (result, "CONJG");
678 }
679
680
681 gfc_expr *
682 gfc_simplify_cos (gfc_expr * x)
683 {
684   gfc_expr *result;
685   mpfr_t xp, xq;
686
687   if (x->expr_type != EXPR_CONSTANT)
688     return NULL;
689
690   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
691
692   switch (x->ts.type)
693     {
694     case BT_REAL:
695       mpfr_cos (result->value.real, x->value.real, GFC_RND_MODE);
696       break;
697     case BT_COMPLEX:
698       gfc_set_model_kind (x->ts.kind);
699       mpfr_init (xp);
700       mpfr_init (xq);
701
702       mpfr_cos  (xp, x->value.complex.r, GFC_RND_MODE);
703       mpfr_cosh (xq, x->value.complex.i, GFC_RND_MODE);
704       mpfr_mul(result->value.complex.r, xp, xq, GFC_RND_MODE);
705
706       mpfr_sin  (xp, x->value.complex.r, GFC_RND_MODE);
707       mpfr_sinh (xq, x->value.complex.i, GFC_RND_MODE);
708       mpfr_mul (xp, xp, xq, GFC_RND_MODE);
709       mpfr_neg (result->value.complex.i, xp, GFC_RND_MODE );
710
711       mpfr_clear (xp);
712       mpfr_clear (xq);
713       break;
714     default:
715       gfc_internal_error ("in gfc_simplify_cos(): Bad type");
716     }
717
718   return range_check (result, "COS");
719
720 }
721
722
723 gfc_expr *
724 gfc_simplify_cosh (gfc_expr * x)
725 {
726   gfc_expr *result;
727
728   if (x->expr_type != EXPR_CONSTANT)
729     return NULL;
730
731   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
732
733   mpfr_cosh (result->value.real, x->value.real, GFC_RND_MODE);
734
735   return range_check (result, "COSH");
736 }
737
738
739 gfc_expr *
740 gfc_simplify_dcmplx (gfc_expr * x, gfc_expr * y)
741 {
742
743   if (x->expr_type != EXPR_CONSTANT
744       || (y != NULL && y->expr_type != EXPR_CONSTANT))
745     return NULL;
746
747   return simplify_cmplx ("DCMPLX", x, y, gfc_default_double_kind);
748 }
749
750
751 gfc_expr *
752 gfc_simplify_dble (gfc_expr * e)
753 {
754   gfc_expr *result;
755
756   if (e->expr_type != EXPR_CONSTANT)
757     return NULL;
758
759   switch (e->ts.type)
760     {
761     case BT_INTEGER:
762       result = gfc_int2real (e, gfc_default_double_kind);
763       break;
764
765     case BT_REAL:
766       result = gfc_real2real (e, gfc_default_double_kind);
767       break;
768
769     case BT_COMPLEX:
770       result = gfc_complex2real (e, gfc_default_double_kind);
771       break;
772
773     default:
774       gfc_internal_error ("gfc_simplify_dble(): bad type at %L", &e->where);
775     }
776
777   return range_check (result, "DBLE");
778 }
779
780
781 gfc_expr *
782 gfc_simplify_digits (gfc_expr * x)
783 {
784   int i, digits;
785
786   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
787   switch (x->ts.type)
788     {
789     case BT_INTEGER:
790       digits = gfc_integer_kinds[i].digits;
791       break;
792
793     case BT_REAL:
794     case BT_COMPLEX:
795       digits = gfc_real_kinds[i].digits;
796       break;
797
798     default:
799       gcc_unreachable ();
800     }
801
802   return gfc_int_expr (digits);
803 }
804
805
806 gfc_expr *
807 gfc_simplify_dim (gfc_expr * x, gfc_expr * y)
808 {
809   gfc_expr *result;
810
811   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
812     return NULL;
813
814   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
815
816   switch (x->ts.type)
817     {
818     case BT_INTEGER:
819       if (mpz_cmp (x->value.integer, y->value.integer) > 0)
820         mpz_sub (result->value.integer, x->value.integer, y->value.integer);
821       else
822         mpz_set_ui (result->value.integer, 0);
823
824       break;
825
826     case BT_REAL:
827       if (mpfr_cmp (x->value.real, y->value.real) > 0)
828         mpfr_sub (result->value.real, x->value.real, y->value.real, GFC_RND_MODE);
829       else
830         mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
831
832       break;
833
834     default:
835       gfc_internal_error ("gfc_simplify_dim(): Bad type");
836     }
837
838   return range_check (result, "DIM");
839 }
840
841
842 gfc_expr *
843 gfc_simplify_dprod (gfc_expr * x, gfc_expr * y)
844 {
845   gfc_expr *a1, *a2, *result;
846
847   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
848     return NULL;
849
850   result =
851     gfc_constant_result (BT_REAL, gfc_default_double_kind, &x->where);
852
853   a1 = gfc_real2real (x, gfc_default_double_kind);
854   a2 = gfc_real2real (y, gfc_default_double_kind);
855
856   mpfr_mul (result->value.real, a1->value.real, a2->value.real, GFC_RND_MODE);
857
858   gfc_free_expr (a1);
859   gfc_free_expr (a2);
860
861   return range_check (result, "DPROD");
862 }
863
864
865 gfc_expr *
866 gfc_simplify_epsilon (gfc_expr * e)
867 {
868   gfc_expr *result;
869   int i;
870
871   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
872
873   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
874
875   mpfr_set (result->value.real, gfc_real_kinds[i].epsilon, GFC_RND_MODE);
876
877   return range_check (result, "EPSILON");
878 }
879
880
881 gfc_expr *
882 gfc_simplify_exp (gfc_expr * x)
883 {
884   gfc_expr *result;
885   mpfr_t xp, xq;
886
887   if (x->expr_type != EXPR_CONSTANT)
888     return NULL;
889
890   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
891
892   switch (x->ts.type)
893     {
894     case BT_REAL:
895       mpfr_exp(result->value.real, x->value.real, GFC_RND_MODE);
896       break;
897
898     case BT_COMPLEX:
899       gfc_set_model_kind (x->ts.kind);
900       mpfr_init (xp);
901       mpfr_init (xq);
902       mpfr_exp (xq, x->value.complex.r, GFC_RND_MODE);
903       mpfr_cos (xp, x->value.complex.i, GFC_RND_MODE);
904       mpfr_mul (result->value.complex.r, xq, xp, GFC_RND_MODE);
905       mpfr_sin (xp, x->value.complex.i, GFC_RND_MODE);
906       mpfr_mul (result->value.complex.i, xq, xp, GFC_RND_MODE);
907       mpfr_clear (xp);
908       mpfr_clear (xq);
909       break;
910
911     default:
912       gfc_internal_error ("in gfc_simplify_exp(): Bad type");
913     }
914
915   return range_check (result, "EXP");
916 }
917
918 /* FIXME:  MPFR should be able to do this better */
919 gfc_expr *
920 gfc_simplify_exponent (gfc_expr * x)
921 {
922   int i;
923   mpfr_t tmp;
924   gfc_expr *result;
925
926   if (x->expr_type != EXPR_CONSTANT)
927     return NULL;
928
929   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
930                                 &x->where);
931
932   gfc_set_model (x->value.real);
933
934   if (mpfr_sgn (x->value.real) == 0)
935     {
936       mpz_set_ui (result->value.integer, 0);
937       return result;
938     }
939
940   mpfr_init (tmp);
941
942   mpfr_abs (tmp, x->value.real, GFC_RND_MODE);
943   mpfr_log2 (tmp, tmp, GFC_RND_MODE);
944
945   gfc_mpfr_to_mpz (result->value.integer, tmp);
946
947   /* The model number for tiny(x) is b**(emin - 1) where b is the base and emin
948      is the smallest exponent value.  So, we need to add 1 if x is tiny(x).  */
949   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
950   if (mpfr_cmp (x->value.real, gfc_real_kinds[i].tiny) == 0)
951     mpz_add_ui (result->value.integer,result->value.integer, 1);
952
953   mpfr_clear (tmp);
954
955   return range_check (result, "EXPONENT");
956 }
957
958
959 gfc_expr *
960 gfc_simplify_float (gfc_expr * a)
961 {
962   gfc_expr *result;
963
964   if (a->expr_type != EXPR_CONSTANT)
965     return NULL;
966
967   result = gfc_int2real (a, gfc_default_real_kind);
968   return range_check (result, "FLOAT");
969 }
970
971
972 gfc_expr *
973 gfc_simplify_floor (gfc_expr * e, gfc_expr * k)
974 {
975   gfc_expr *result;
976   mpfr_t floor;
977   int kind;
978
979   kind = get_kind (BT_INTEGER, k, "FLOOR", gfc_default_integer_kind);
980   if (kind == -1)
981     gfc_internal_error ("gfc_simplify_floor(): Bad kind");
982
983   if (e->expr_type != EXPR_CONSTANT)
984     return NULL;
985
986   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
987
988   gfc_set_model_kind (kind);
989   mpfr_init (floor);
990   mpfr_floor (floor, e->value.real);
991
992   gfc_mpfr_to_mpz (result->value.integer, floor);
993
994   mpfr_clear (floor);
995
996   return range_check (result, "FLOOR");
997 }
998
999
1000 gfc_expr *
1001 gfc_simplify_fraction (gfc_expr * x)
1002 {
1003   gfc_expr *result;
1004   mpfr_t absv, exp, pow2;
1005
1006   if (x->expr_type != EXPR_CONSTANT)
1007     return NULL;
1008
1009   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
1010
1011   gfc_set_model_kind (x->ts.kind);
1012
1013   if (mpfr_sgn (x->value.real) == 0)
1014     {
1015       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
1016       return result;
1017     }
1018
1019   mpfr_init (exp);
1020   mpfr_init (absv);
1021   mpfr_init (pow2);
1022
1023   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
1024   mpfr_log2 (exp, absv, GFC_RND_MODE);
1025
1026   mpfr_trunc (exp, exp);
1027   mpfr_add_ui (exp, exp, 1, GFC_RND_MODE);
1028
1029   mpfr_ui_pow (pow2, 2, exp, GFC_RND_MODE);
1030
1031   mpfr_div (result->value.real, absv, pow2, GFC_RND_MODE);
1032
1033   mpfr_clear (exp);
1034   mpfr_clear (absv);
1035   mpfr_clear (pow2);
1036
1037   return range_check (result, "FRACTION");
1038 }
1039
1040
1041 gfc_expr *
1042 gfc_simplify_huge (gfc_expr * e)
1043 {
1044   gfc_expr *result;
1045   int i;
1046
1047   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
1048
1049   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
1050
1051   switch (e->ts.type)
1052     {
1053     case BT_INTEGER:
1054       mpz_set (result->value.integer, gfc_integer_kinds[i].huge);
1055       break;
1056
1057     case BT_REAL:
1058       mpfr_set (result->value.real, gfc_real_kinds[i].huge, GFC_RND_MODE);
1059       break;
1060
1061     default:
1062       gcc_unreachable ();
1063     }
1064
1065   return result;
1066 }
1067
1068
1069 gfc_expr *
1070 gfc_simplify_iachar (gfc_expr * e)
1071 {
1072   gfc_expr *result;
1073   int index;
1074
1075   if (e->expr_type != EXPR_CONSTANT)
1076     return NULL;
1077
1078   if (e->value.character.length != 1)
1079     {
1080       gfc_error ("Argument of IACHAR at %L must be of length one", &e->where);
1081       return &gfc_bad_expr;
1082     }
1083
1084   index = xascii_table[(int) e->value.character.string[0] & 0xFF];
1085
1086   result = gfc_int_expr (index);
1087   result->where = e->where;
1088
1089   return range_check (result, "IACHAR");
1090 }
1091
1092
1093 gfc_expr *
1094 gfc_simplify_iand (gfc_expr * x, gfc_expr * y)
1095 {
1096   gfc_expr *result;
1097
1098   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1099     return NULL;
1100
1101   result = gfc_constant_result (BT_INTEGER, x->ts.kind, &x->where);
1102
1103   mpz_and (result->value.integer, x->value.integer, y->value.integer);
1104
1105   return range_check (result, "IAND");
1106 }
1107
1108
1109 gfc_expr *
1110 gfc_simplify_ibclr (gfc_expr * x, gfc_expr * y)
1111 {
1112   gfc_expr *result;
1113   int k, pos;
1114
1115   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1116     return NULL;
1117
1118   if (gfc_extract_int (y, &pos) != NULL || pos < 0)
1119     {
1120       gfc_error ("Invalid second argument of IBCLR at %L", &y->where);
1121       return &gfc_bad_expr;
1122     }
1123
1124   k = gfc_validate_kind (x->ts.type, x->ts.kind, false);
1125
1126   if (pos > gfc_integer_kinds[k].bit_size)
1127     {
1128       gfc_error ("Second argument of IBCLR exceeds bit size at %L",
1129                  &y->where);
1130       return &gfc_bad_expr;
1131     }
1132
1133   result = gfc_copy_expr (x);
1134
1135   mpz_clrbit (result->value.integer, pos);
1136   return range_check (result, "IBCLR");
1137 }
1138
1139
1140 gfc_expr *
1141 gfc_simplify_ibits (gfc_expr * x, gfc_expr * y, gfc_expr * z)
1142 {
1143   gfc_expr *result;
1144   int pos, len;
1145   int i, k, bitsize;
1146   int *bits;
1147
1148   if (x->expr_type != EXPR_CONSTANT
1149       || y->expr_type != EXPR_CONSTANT
1150       || z->expr_type != EXPR_CONSTANT)
1151     return NULL;
1152
1153   if (gfc_extract_int (y, &pos) != NULL || pos < 0)
1154     {
1155       gfc_error ("Invalid second argument of IBITS at %L", &y->where);
1156       return &gfc_bad_expr;
1157     }
1158
1159   if (gfc_extract_int (z, &len) != NULL || len < 0)
1160     {
1161       gfc_error ("Invalid third argument of IBITS at %L", &z->where);
1162       return &gfc_bad_expr;
1163     }
1164
1165   k = gfc_validate_kind (BT_INTEGER, x->ts.kind, false);
1166
1167   bitsize = gfc_integer_kinds[k].bit_size;
1168
1169   if (pos + len > bitsize)
1170     {
1171       gfc_error
1172         ("Sum of second and third arguments of IBITS exceeds bit size "
1173          "at %L", &y->where);
1174       return &gfc_bad_expr;
1175     }
1176
1177   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1178
1179   bits = gfc_getmem (bitsize * sizeof (int));
1180
1181   for (i = 0; i < bitsize; i++)
1182     bits[i] = 0;
1183
1184   for (i = 0; i < len; i++)
1185     bits[i] = mpz_tstbit (x->value.integer, i + pos);
1186
1187   for (i = 0; i < bitsize; i++)
1188     {
1189       if (bits[i] == 0)
1190         {
1191           mpz_clrbit (result->value.integer, i);
1192         }
1193       else if (bits[i] == 1)
1194         {
1195           mpz_setbit (result->value.integer, i);
1196         }
1197       else
1198         {
1199           gfc_internal_error ("IBITS: Bad bit");
1200         }
1201     }
1202
1203   gfc_free (bits);
1204
1205   return range_check (result, "IBITS");
1206 }
1207
1208
1209 gfc_expr *
1210 gfc_simplify_ibset (gfc_expr * x, gfc_expr * y)
1211 {
1212   gfc_expr *result;
1213   int k, pos;
1214
1215   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1216     return NULL;
1217
1218   if (gfc_extract_int (y, &pos) != NULL || pos < 0)
1219     {
1220       gfc_error ("Invalid second argument of IBSET at %L", &y->where);
1221       return &gfc_bad_expr;
1222     }
1223
1224   k = gfc_validate_kind (x->ts.type, x->ts.kind, false);
1225
1226   if (pos > gfc_integer_kinds[k].bit_size)
1227     {
1228       gfc_error ("Second argument of IBSET exceeds bit size at %L",
1229                  &y->where);
1230       return &gfc_bad_expr;
1231     }
1232
1233   result = gfc_copy_expr (x);
1234
1235   mpz_setbit (result->value.integer, pos);
1236   return range_check (result, "IBSET");
1237 }
1238
1239
1240 gfc_expr *
1241 gfc_simplify_ichar (gfc_expr * e)
1242 {
1243   gfc_expr *result;
1244   int index;
1245
1246   if (e->expr_type != EXPR_CONSTANT)
1247     return NULL;
1248
1249   if (e->value.character.length != 1)
1250     {
1251       gfc_error ("Argument of ICHAR at %L must be of length one", &e->where);
1252       return &gfc_bad_expr;
1253     }
1254
1255   index = (int) e->value.character.string[0];
1256
1257   if (index < CHAR_MIN || index > CHAR_MAX)
1258     {
1259       gfc_error ("Argument of ICHAR at %L out of range of this processor",
1260                  &e->where);
1261       return &gfc_bad_expr;
1262     }
1263
1264   result = gfc_int_expr (index);
1265   result->where = e->where;
1266   return range_check (result, "ICHAR");
1267 }
1268
1269
1270 gfc_expr *
1271 gfc_simplify_ieor (gfc_expr * x, gfc_expr * y)
1272 {
1273   gfc_expr *result;
1274
1275   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1276     return NULL;
1277
1278   result = gfc_constant_result (BT_INTEGER, x->ts.kind, &x->where);
1279
1280   mpz_xor (result->value.integer, x->value.integer, y->value.integer);
1281
1282   return range_check (result, "IEOR");
1283 }
1284
1285
1286 gfc_expr *
1287 gfc_simplify_index (gfc_expr * x, gfc_expr * y, gfc_expr * b)
1288 {
1289   gfc_expr *result;
1290   int back, len, lensub;
1291   int i, j, k, count, index = 0, start;
1292
1293   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1294     return NULL;
1295
1296   if (b != NULL && b->value.logical != 0)
1297     back = 1;
1298   else
1299     back = 0;
1300
1301   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1302                                 &x->where);
1303
1304   len = x->value.character.length;
1305   lensub = y->value.character.length;
1306
1307   if (len < lensub)
1308     {
1309       mpz_set_si (result->value.integer, 0);
1310       return result;
1311     }
1312
1313   if (back == 0)
1314     {
1315
1316       if (lensub == 0)
1317         {
1318           mpz_set_si (result->value.integer, 1);
1319           return result;
1320         }
1321       else if (lensub == 1)
1322         {
1323           for (i = 0; i < len; i++)
1324             {
1325               for (j = 0; j < lensub; j++)
1326                 {
1327                   if (y->value.character.string[j] ==
1328                       x->value.character.string[i])
1329                     {
1330                       index = i + 1;
1331                       goto done;
1332                     }
1333                 }
1334             }
1335         }
1336       else
1337         {
1338           for (i = 0; i < len; i++)
1339             {
1340               for (j = 0; j < lensub; j++)
1341                 {
1342                   if (y->value.character.string[j] ==
1343                       x->value.character.string[i])
1344                     {
1345                       start = i;
1346                       count = 0;
1347
1348                       for (k = 0; k < lensub; k++)
1349                         {
1350                           if (y->value.character.string[k] ==
1351                               x->value.character.string[k + start])
1352                             count++;
1353                         }
1354
1355                       if (count == lensub)
1356                         {
1357                           index = start + 1;
1358                           goto done;
1359                         }
1360                     }
1361                 }
1362             }
1363         }
1364
1365     }
1366   else
1367     {
1368
1369       if (lensub == 0)
1370         {
1371           mpz_set_si (result->value.integer, len + 1);
1372           return result;
1373         }
1374       else if (lensub == 1)
1375         {
1376           for (i = 0; i < len; i++)
1377             {
1378               for (j = 0; j < lensub; j++)
1379                 {
1380                   if (y->value.character.string[j] ==
1381                       x->value.character.string[len - i])
1382                     {
1383                       index = len - i + 1;
1384                       goto done;
1385                     }
1386                 }
1387             }
1388         }
1389       else
1390         {
1391           for (i = 0; i < len; i++)
1392             {
1393               for (j = 0; j < lensub; j++)
1394                 {
1395                   if (y->value.character.string[j] ==
1396                       x->value.character.string[len - i])
1397                     {
1398                       start = len - i;
1399                       if (start <= len - lensub)
1400                         {
1401                           count = 0;
1402                           for (k = 0; k < lensub; k++)
1403                             if (y->value.character.string[k] ==
1404                                 x->value.character.string[k + start])
1405                               count++;
1406
1407                           if (count == lensub)
1408                             {
1409                               index = start + 1;
1410                               goto done;
1411                             }
1412                         }
1413                       else
1414                         {
1415                           continue;
1416                         }
1417                     }
1418                 }
1419             }
1420         }
1421     }
1422
1423 done:
1424   mpz_set_si (result->value.integer, index);
1425   return range_check (result, "INDEX");
1426 }
1427
1428
1429 gfc_expr *
1430 gfc_simplify_int (gfc_expr * e, gfc_expr * k)
1431 {
1432   gfc_expr *rpart, *rtrunc, *result;
1433   int kind;
1434
1435   kind = get_kind (BT_INTEGER, k, "INT", gfc_default_integer_kind);
1436   if (kind == -1)
1437     return &gfc_bad_expr;
1438
1439   if (e->expr_type != EXPR_CONSTANT)
1440     return NULL;
1441
1442   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
1443
1444   switch (e->ts.type)
1445     {
1446     case BT_INTEGER:
1447       mpz_set (result->value.integer, e->value.integer);
1448       break;
1449
1450     case BT_REAL:
1451       rtrunc = gfc_copy_expr (e);
1452       mpfr_trunc (rtrunc->value.real, e->value.real);
1453       gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real);
1454       gfc_free_expr (rtrunc);
1455       break;
1456
1457     case BT_COMPLEX:
1458       rpart = gfc_complex2real (e, kind);
1459       rtrunc = gfc_copy_expr (rpart);
1460       mpfr_trunc (rtrunc->value.real, rpart->value.real);
1461       gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real);
1462       gfc_free_expr (rpart);
1463       gfc_free_expr (rtrunc);
1464       break;
1465
1466     default:
1467       gfc_error ("Argument of INT at %L is not a valid type", &e->where);
1468       gfc_free_expr (result);
1469       return &gfc_bad_expr;
1470     }
1471
1472   return range_check (result, "INT");
1473 }
1474
1475
1476 gfc_expr *
1477 gfc_simplify_ifix (gfc_expr * e)
1478 {
1479   gfc_expr *rtrunc, *result;
1480
1481   if (e->expr_type != EXPR_CONSTANT)
1482     return NULL;
1483
1484   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1485                                 &e->where);
1486
1487   rtrunc = gfc_copy_expr (e);
1488
1489   mpfr_trunc (rtrunc->value.real, e->value.real);
1490   gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real);
1491
1492   gfc_free_expr (rtrunc);
1493   return range_check (result, "IFIX");
1494 }
1495
1496
1497 gfc_expr *
1498 gfc_simplify_idint (gfc_expr * e)
1499 {
1500   gfc_expr *rtrunc, *result;
1501
1502   if (e->expr_type != EXPR_CONSTANT)
1503     return NULL;
1504
1505   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1506                                 &e->where);
1507
1508   rtrunc = gfc_copy_expr (e);
1509
1510   mpfr_trunc (rtrunc->value.real, e->value.real);
1511   gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real);
1512
1513   gfc_free_expr (rtrunc);
1514   return range_check (result, "IDINT");
1515 }
1516
1517
1518 gfc_expr *
1519 gfc_simplify_ior (gfc_expr * x, gfc_expr * y)
1520 {
1521   gfc_expr *result;
1522
1523   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1524     return NULL;
1525
1526   result = gfc_constant_result (BT_INTEGER, x->ts.kind, &x->where);
1527
1528   mpz_ior (result->value.integer, x->value.integer, y->value.integer);
1529   return range_check (result, "IOR");
1530 }
1531
1532
1533 gfc_expr *
1534 gfc_simplify_ishft (gfc_expr * e, gfc_expr * s)
1535 {
1536   gfc_expr *result;
1537   int shift, ashift, isize, k, *bits, i;
1538
1539   if (e->expr_type != EXPR_CONSTANT || s->expr_type != EXPR_CONSTANT)
1540     return NULL;
1541
1542   if (gfc_extract_int (s, &shift) != NULL)
1543     {
1544       gfc_error ("Invalid second argument of ISHFT at %L", &s->where);
1545       return &gfc_bad_expr;
1546     }
1547
1548   k = gfc_validate_kind (BT_INTEGER, e->ts.kind, false);
1549
1550   isize = gfc_integer_kinds[k].bit_size;
1551
1552   if (shift >= 0)
1553     ashift = shift;
1554   else
1555     ashift = -shift;
1556
1557   if (ashift > isize)
1558     {
1559       gfc_error
1560         ("Magnitude of second argument of ISHFT exceeds bit size at %L",
1561          &s->where);
1562       return &gfc_bad_expr;
1563     }
1564
1565   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
1566
1567   if (shift == 0)
1568     {
1569       mpz_set (result->value.integer, e->value.integer);
1570       return range_check (result, "ISHFT");
1571     }
1572   
1573   bits = gfc_getmem (isize * sizeof (int));
1574
1575   for (i = 0; i < isize; i++)
1576     bits[i] = mpz_tstbit (e->value.integer, i);
1577
1578   if (shift > 0)
1579     {
1580       for (i = 0; i < shift; i++)
1581         mpz_clrbit (result->value.integer, i);
1582
1583       for (i = 0; i < isize - shift; i++)
1584         {
1585           if (bits[i] == 0)
1586             mpz_clrbit (result->value.integer, i + shift);
1587           else
1588             mpz_setbit (result->value.integer, i + shift);
1589         }
1590     }
1591   else
1592     {
1593       for (i = isize - 1; i >= isize - ashift; i--)
1594         mpz_clrbit (result->value.integer, i);
1595
1596       for (i = isize - 1; i >= ashift; i--)
1597         {
1598           if (bits[i] == 0)
1599             mpz_clrbit (result->value.integer, i - ashift);
1600           else
1601             mpz_setbit (result->value.integer, i - ashift);
1602         }
1603     }
1604
1605   twos_complement (result->value.integer, isize);
1606
1607   gfc_free (bits);
1608   return result;
1609 }
1610
1611
1612 gfc_expr *
1613 gfc_simplify_ishftc (gfc_expr * e, gfc_expr * s, gfc_expr * sz)
1614 {
1615   gfc_expr *result;
1616   int shift, ashift, isize, delta, k;
1617   int i, *bits;
1618
1619   if (e->expr_type != EXPR_CONSTANT || s->expr_type != EXPR_CONSTANT)
1620     return NULL;
1621
1622   if (gfc_extract_int (s, &shift) != NULL)
1623     {
1624       gfc_error ("Invalid second argument of ISHFTC at %L", &s->where);
1625       return &gfc_bad_expr;
1626     }
1627
1628   k = gfc_validate_kind (e->ts.type, e->ts.kind, false);
1629
1630   if (sz != NULL)
1631     {
1632       if (gfc_extract_int (sz, &isize) != NULL || isize < 0)
1633         {
1634           gfc_error ("Invalid third argument of ISHFTC at %L", &sz->where);
1635           return &gfc_bad_expr;
1636         }
1637     }
1638   else
1639     isize = gfc_integer_kinds[k].bit_size;
1640
1641   if (shift >= 0)
1642     ashift = shift;
1643   else
1644     ashift = -shift;
1645
1646   if (ashift > isize)
1647     {
1648       gfc_error
1649         ("Magnitude of second argument of ISHFTC exceeds third argument "
1650          "at %L", &s->where);
1651       return &gfc_bad_expr;
1652     }
1653
1654   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
1655
1656   if (shift == 0)
1657     {
1658       mpz_set (result->value.integer, e->value.integer);
1659       return result;
1660     }
1661
1662   bits = gfc_getmem (isize * sizeof (int));
1663
1664   for (i = 0; i < isize; i++)
1665     bits[i] = mpz_tstbit (e->value.integer, i);
1666
1667   delta = isize - ashift;
1668
1669   if (shift > 0)
1670     {
1671       for (i = 0; i < delta; i++)
1672         {
1673           if (bits[i] == 0)
1674             mpz_clrbit (result->value.integer, i + shift);
1675           else
1676             mpz_setbit (result->value.integer, i + shift);
1677         }
1678
1679       for (i = delta; i < isize; i++)
1680         {
1681           if (bits[i] == 0)
1682             mpz_clrbit (result->value.integer, i - delta);
1683           else
1684             mpz_setbit (result->value.integer, i - delta);
1685         }
1686     }
1687   else
1688     {
1689       for (i = 0; i < ashift; i++)
1690         {
1691           if (bits[i] == 0)
1692             mpz_clrbit (result->value.integer, i + delta);
1693           else
1694             mpz_setbit (result->value.integer, i + delta);
1695         }
1696
1697       for (i = ashift; i < isize; i++)
1698         {
1699           if (bits[i] == 0)
1700             mpz_clrbit (result->value.integer, i + shift);
1701           else
1702             mpz_setbit (result->value.integer, i + shift);
1703         }
1704     }
1705
1706   twos_complement (result->value.integer, isize);
1707
1708   gfc_free (bits);
1709   return result;
1710 }
1711
1712
1713 gfc_expr *
1714 gfc_simplify_kind (gfc_expr * e)
1715 {
1716
1717   if (e->ts.type == BT_DERIVED)
1718     {
1719       gfc_error ("Argument of KIND at %L is a DERIVED type", &e->where);
1720       return &gfc_bad_expr;
1721     }
1722
1723   return gfc_int_expr (e->ts.kind);
1724 }
1725
1726
1727 static gfc_expr *
1728 simplify_bound (gfc_expr * array, gfc_expr * dim, int upper)
1729 {
1730   gfc_ref *ref;
1731   gfc_array_spec *as;
1732   gfc_expr *e;
1733   int d;
1734
1735   if (array->expr_type != EXPR_VARIABLE)
1736     return NULL;
1737
1738   if (dim == NULL)
1739     /* TODO: Simplify constant multi-dimensional bounds.  */
1740     return NULL;
1741
1742   if (dim->expr_type != EXPR_CONSTANT)
1743     return NULL;
1744
1745   /* Follow any component references.  */
1746   as = array->symtree->n.sym->as;
1747   for (ref = array->ref; ref; ref = ref->next)
1748     {
1749       switch (ref->type)
1750         {
1751         case REF_ARRAY:
1752           switch (ref->u.ar.type)
1753             {
1754             case AR_ELEMENT:
1755               as = NULL;
1756               continue;
1757
1758             case AR_FULL:
1759               /* We're done because 'as' has already been set in the
1760                  previous iteration.  */
1761               goto done;
1762
1763             case AR_SECTION:
1764             case AR_UNKNOWN:
1765               return NULL;
1766             }
1767
1768           gcc_unreachable ();
1769
1770         case REF_COMPONENT:
1771           as = ref->u.c.component->as;
1772           continue;
1773
1774         case REF_SUBSTRING:
1775           continue;
1776         }
1777     }
1778
1779   gcc_unreachable ();
1780
1781  done:
1782   if (as->type == AS_DEFERRED || as->type == AS_ASSUMED_SHAPE)
1783     return NULL;
1784
1785   d = mpz_get_si (dim->value.integer);
1786
1787   if (d < 1 || d > as->rank
1788       || (d == as->rank && as->type == AS_ASSUMED_SIZE && upper))
1789     {
1790       gfc_error ("DIM argument at %L is out of bounds", &dim->where);
1791       return &gfc_bad_expr;
1792     }
1793
1794   e = upper ? as->upper[d-1] : as->lower[d-1];
1795
1796   if (e->expr_type != EXPR_CONSTANT)
1797     return NULL;
1798
1799   return gfc_copy_expr (e);
1800 }
1801
1802
1803 gfc_expr *
1804 gfc_simplify_lbound (gfc_expr * array, gfc_expr * dim)
1805 {
1806   return simplify_bound (array, dim, 0);
1807 }
1808
1809
1810 gfc_expr *
1811 gfc_simplify_len (gfc_expr * e)
1812 {
1813   gfc_expr *result;
1814
1815   if (e->expr_type != EXPR_CONSTANT)
1816     return NULL;
1817
1818   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1819                                 &e->where);
1820
1821   mpz_set_si (result->value.integer, e->value.character.length);
1822   return range_check (result, "LEN");
1823 }
1824
1825
1826 gfc_expr *
1827 gfc_simplify_len_trim (gfc_expr * e)
1828 {
1829   gfc_expr *result;
1830   int count, len, lentrim, i;
1831
1832   if (e->expr_type != EXPR_CONSTANT)
1833     return NULL;
1834
1835   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1836                                 &e->where);
1837
1838   len = e->value.character.length;
1839
1840   for (count = 0, i = 1; i <= len; i++)
1841     if (e->value.character.string[len - i] == ' ')
1842       count++;
1843     else
1844       break;
1845
1846   lentrim = len - count;
1847
1848   mpz_set_si (result->value.integer, lentrim);
1849   return range_check (result, "LEN_TRIM");
1850 }
1851
1852
1853 gfc_expr *
1854 gfc_simplify_lge (gfc_expr * a, gfc_expr * b)
1855 {
1856
1857   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
1858     return NULL;
1859
1860   return gfc_logical_expr (gfc_compare_string (a, b, xascii_table) >= 0,
1861                            &a->where);
1862 }
1863
1864
1865 gfc_expr *
1866 gfc_simplify_lgt (gfc_expr * a, gfc_expr * b)
1867 {
1868
1869   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
1870     return NULL;
1871
1872   return gfc_logical_expr (gfc_compare_string (a, b, xascii_table) > 0,
1873                            &a->where);
1874 }
1875
1876
1877 gfc_expr *
1878 gfc_simplify_lle (gfc_expr * a, gfc_expr * b)
1879 {
1880
1881   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
1882     return NULL;
1883
1884   return gfc_logical_expr (gfc_compare_string (a, b, xascii_table) <= 0,
1885                            &a->where);
1886 }
1887
1888
1889 gfc_expr *
1890 gfc_simplify_llt (gfc_expr * a, gfc_expr * b)
1891 {
1892
1893   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
1894     return NULL;
1895
1896   return gfc_logical_expr (gfc_compare_string (a, b, xascii_table) < 0,
1897                            &a->where);
1898 }
1899
1900
1901 gfc_expr *
1902 gfc_simplify_log (gfc_expr * x)
1903 {
1904   gfc_expr *result;
1905   mpfr_t xr, xi;
1906
1907   if (x->expr_type != EXPR_CONSTANT)
1908     return NULL;
1909
1910   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1911
1912   gfc_set_model_kind (x->ts.kind);
1913
1914   switch (x->ts.type)
1915     {
1916     case BT_REAL:
1917       if (mpfr_sgn (x->value.real) <= 0)
1918         {
1919           gfc_error
1920             ("Argument of LOG at %L cannot be less than or equal to zero",
1921              &x->where);
1922           gfc_free_expr (result);
1923           return &gfc_bad_expr;
1924         }
1925
1926       mpfr_log(result->value.real, x->value.real, GFC_RND_MODE);
1927       break;
1928
1929     case BT_COMPLEX:
1930       if ((mpfr_sgn (x->value.complex.r) == 0)
1931           && (mpfr_sgn (x->value.complex.i) == 0))
1932         {
1933           gfc_error ("Complex argument of LOG at %L cannot be zero",
1934                      &x->where);
1935           gfc_free_expr (result);
1936           return &gfc_bad_expr;
1937         }
1938
1939       mpfr_init (xr);
1940       mpfr_init (xi);
1941
1942       arctangent2 (x->value.complex.i, x->value.complex.r,
1943                    result->value.complex.i);
1944
1945       mpfr_mul (xr, x->value.complex.r, x->value.complex.r, GFC_RND_MODE);
1946       mpfr_mul (xi, x->value.complex.i, x->value.complex.i, GFC_RND_MODE);
1947       mpfr_add (xr, xr, xi, GFC_RND_MODE);
1948       mpfr_sqrt (xr, xr, GFC_RND_MODE);
1949       mpfr_log (result->value.complex.r, xr, GFC_RND_MODE);
1950
1951       mpfr_clear (xr);
1952       mpfr_clear (xi);
1953
1954       break;
1955
1956     default:
1957       gfc_internal_error ("gfc_simplify_log: bad type");
1958     }
1959
1960   return range_check (result, "LOG");
1961 }
1962
1963
1964 gfc_expr *
1965 gfc_simplify_log10 (gfc_expr * x)
1966 {
1967   gfc_expr *result;
1968
1969   if (x->expr_type != EXPR_CONSTANT)
1970     return NULL;
1971
1972   gfc_set_model_kind (x->ts.kind);
1973
1974   if (mpfr_sgn (x->value.real) <= 0)
1975     {
1976       gfc_error
1977         ("Argument of LOG10 at %L cannot be less than or equal to zero",
1978          &x->where);
1979       return &gfc_bad_expr;
1980     }
1981
1982   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1983
1984   mpfr_log10 (result->value.real, x->value.real, GFC_RND_MODE);
1985
1986   return range_check (result, "LOG10");
1987 }
1988
1989
1990 gfc_expr *
1991 gfc_simplify_logical (gfc_expr * e, gfc_expr * k)
1992 {
1993   gfc_expr *result;
1994   int kind;
1995
1996   kind = get_kind (BT_LOGICAL, k, "LOGICAL", gfc_default_logical_kind);
1997   if (kind < 0)
1998     return &gfc_bad_expr;
1999
2000   if (e->expr_type != EXPR_CONSTANT)
2001     return NULL;
2002
2003   result = gfc_constant_result (BT_LOGICAL, kind, &e->where);
2004
2005   result->value.logical = e->value.logical;
2006
2007   return result;
2008 }
2009
2010
2011 /* This function is special since MAX() can take any number of
2012    arguments.  The simplified expression is a rewritten version of the
2013    argument list containing at most one constant element.  Other
2014    constant elements are deleted.  Because the argument list has
2015    already been checked, this function always succeeds.  sign is 1 for
2016    MAX(), -1 for MIN().  */
2017
2018 static gfc_expr *
2019 simplify_min_max (gfc_expr * expr, int sign)
2020 {
2021   gfc_actual_arglist *arg, *last, *extremum;
2022   gfc_intrinsic_sym * specific;
2023
2024   last = NULL;
2025   extremum = NULL;
2026   specific = expr->value.function.isym;
2027
2028   arg = expr->value.function.actual;
2029
2030   for (; arg; last = arg, arg = arg->next)
2031     {
2032       if (arg->expr->expr_type != EXPR_CONSTANT)
2033         continue;
2034
2035       if (extremum == NULL)
2036         {
2037           extremum = arg;
2038           continue;
2039         }
2040
2041       switch (arg->expr->ts.type)
2042         {
2043         case BT_INTEGER:
2044           if (mpz_cmp (arg->expr->value.integer,
2045                        extremum->expr->value.integer) * sign > 0)
2046             mpz_set (extremum->expr->value.integer, arg->expr->value.integer);
2047
2048           break;
2049
2050         case BT_REAL:
2051           if (mpfr_cmp (arg->expr->value.real, extremum->expr->value.real) *
2052               sign > 0)
2053             mpfr_set (extremum->expr->value.real, arg->expr->value.real,
2054                       GFC_RND_MODE);
2055
2056           break;
2057
2058         default:
2059           gfc_internal_error ("gfc_simplify_max(): Bad type in arglist");
2060         }
2061
2062       /* Delete the extra constant argument.  */
2063       if (last == NULL)
2064         expr->value.function.actual = arg->next;
2065       else
2066         last->next = arg->next;
2067
2068       arg->next = NULL;
2069       gfc_free_actual_arglist (arg);
2070       arg = last;
2071     }
2072
2073   /* If there is one value left, replace the function call with the
2074      expression.  */
2075   if (expr->value.function.actual->next != NULL)
2076     return NULL;
2077
2078   /* Convert to the correct type and kind.  */
2079   if (expr->ts.type != BT_UNKNOWN) 
2080     return gfc_convert_constant (expr->value.function.actual->expr,
2081         expr->ts.type, expr->ts.kind);
2082
2083   if (specific->ts.type != BT_UNKNOWN) 
2084     return gfc_convert_constant (expr->value.function.actual->expr,
2085         specific->ts.type, specific->ts.kind); 
2086  
2087   return gfc_copy_expr (expr->value.function.actual->expr);
2088 }
2089
2090
2091 gfc_expr *
2092 gfc_simplify_min (gfc_expr * e)
2093 {
2094   return simplify_min_max (e, -1);
2095 }
2096
2097
2098 gfc_expr *
2099 gfc_simplify_max (gfc_expr * e)
2100 {
2101   return simplify_min_max (e, 1);
2102 }
2103
2104
2105 gfc_expr *
2106 gfc_simplify_maxexponent (gfc_expr * x)
2107 {
2108   gfc_expr *result;
2109   int i;
2110
2111   i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
2112
2113   result = gfc_int_expr (gfc_real_kinds[i].max_exponent);
2114   result->where = x->where;
2115
2116   return result;
2117 }
2118
2119
2120 gfc_expr *
2121 gfc_simplify_minexponent (gfc_expr * x)
2122 {
2123   gfc_expr *result;
2124   int i;
2125
2126   i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
2127
2128   result = gfc_int_expr (gfc_real_kinds[i].min_exponent);
2129   result->where = x->where;
2130
2131   return result;
2132 }
2133
2134
2135 gfc_expr *
2136 gfc_simplify_mod (gfc_expr * a, gfc_expr * p)
2137 {
2138   gfc_expr *result;
2139   mpfr_t quot, iquot, term;
2140
2141   if (a->expr_type != EXPR_CONSTANT || p->expr_type != EXPR_CONSTANT)
2142     return NULL;
2143
2144   result = gfc_constant_result (a->ts.type, a->ts.kind, &a->where);
2145
2146   switch (a->ts.type)
2147     {
2148     case BT_INTEGER:
2149       if (mpz_cmp_ui (p->value.integer, 0) == 0)
2150         {
2151           /* Result is processor-dependent.  */
2152           gfc_error ("Second argument MOD at %L is zero", &a->where);
2153           gfc_free_expr (result);
2154           return &gfc_bad_expr;
2155         }
2156       mpz_tdiv_r (result->value.integer, a->value.integer, p->value.integer);
2157       break;
2158
2159     case BT_REAL:
2160       if (mpfr_cmp_ui (p->value.real, 0) == 0)
2161         {
2162           /* Result is processor-dependent.  */
2163           gfc_error ("Second argument of MOD at %L is zero", &p->where);
2164           gfc_free_expr (result);
2165           return &gfc_bad_expr;
2166         }
2167
2168       gfc_set_model_kind (a->ts.kind);
2169       mpfr_init (quot);
2170       mpfr_init (iquot);
2171       mpfr_init (term);
2172
2173       mpfr_div (quot, a->value.real, p->value.real, GFC_RND_MODE);
2174       mpfr_trunc (iquot, quot);
2175       mpfr_mul (term, iquot, p->value.real, GFC_RND_MODE);
2176       mpfr_sub (result->value.real, a->value.real, term, GFC_RND_MODE);
2177
2178       mpfr_clear (quot);
2179       mpfr_clear (iquot);
2180       mpfr_clear (term);
2181       break;
2182
2183     default:
2184       gfc_internal_error ("gfc_simplify_mod(): Bad arguments");
2185     }
2186
2187   return range_check (result, "MOD");
2188 }
2189
2190
2191 gfc_expr *
2192 gfc_simplify_modulo (gfc_expr * a, gfc_expr * p)
2193 {
2194   gfc_expr *result;
2195   mpfr_t quot, iquot, term;
2196
2197   if (a->expr_type != EXPR_CONSTANT || p->expr_type != EXPR_CONSTANT)
2198     return NULL;
2199
2200   result = gfc_constant_result (a->ts.type, a->ts.kind, &a->where);
2201
2202   switch (a->ts.type)
2203     {
2204     case BT_INTEGER:
2205       if (mpz_cmp_ui (p->value.integer, 0) == 0)
2206         {
2207           /* Result is processor-dependent. This processor just opts
2208              to not handle it at all.  */
2209           gfc_error ("Second argument of MODULO at %L is zero", &a->where);
2210           gfc_free_expr (result);
2211           return &gfc_bad_expr;
2212         }
2213       mpz_fdiv_r (result->value.integer, a->value.integer, p->value.integer);
2214
2215       break;
2216
2217     case BT_REAL:
2218       if (mpfr_cmp_ui (p->value.real, 0) == 0)
2219         {
2220           /* Result is processor-dependent.  */
2221           gfc_error ("Second argument of MODULO at %L is zero", &p->where);
2222           gfc_free_expr (result);
2223           return &gfc_bad_expr;
2224         }
2225
2226       gfc_set_model_kind (a->ts.kind);
2227       mpfr_init (quot);
2228       mpfr_init (iquot);
2229       mpfr_init (term);
2230
2231       mpfr_div (quot, a->value.real, p->value.real, GFC_RND_MODE);
2232       mpfr_floor (iquot, quot);
2233       mpfr_mul (term, iquot, p->value.real, GFC_RND_MODE);
2234
2235       mpfr_clear (quot);
2236       mpfr_clear (iquot);
2237       mpfr_clear (term);
2238
2239       mpfr_sub (result->value.real, a->value.real, term, GFC_RND_MODE);
2240       break;
2241
2242     default:
2243       gfc_internal_error ("gfc_simplify_modulo(): Bad arguments");
2244     }
2245
2246   return range_check (result, "MODULO");
2247 }
2248
2249
2250 /* Exists for the sole purpose of consistency with other intrinsics.  */
2251 gfc_expr *
2252 gfc_simplify_mvbits (gfc_expr * f  ATTRIBUTE_UNUSED,
2253                      gfc_expr * fp ATTRIBUTE_UNUSED,
2254                      gfc_expr * l  ATTRIBUTE_UNUSED,
2255                      gfc_expr * to ATTRIBUTE_UNUSED,
2256                      gfc_expr * tp ATTRIBUTE_UNUSED)
2257 {
2258   return NULL;
2259 }
2260
2261
2262 gfc_expr *
2263 gfc_simplify_nearest (gfc_expr * x, gfc_expr * s)
2264 {
2265   gfc_expr *result;
2266   float rval;
2267   double val, eps;
2268   int p, i, k, match_float;
2269
2270   /* FIXME: This implementation is dopey and probably not quite right,
2271      but it's a start.  */
2272
2273   if (x->expr_type != EXPR_CONSTANT)
2274     return NULL;
2275
2276   k = gfc_validate_kind (x->ts.type, x->ts.kind, false);
2277
2278   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
2279
2280   val = mpfr_get_d (x->value.real, GFC_RND_MODE);
2281   p = gfc_real_kinds[k].digits;
2282
2283   eps = 1.;
2284   for (i = 1; i < p; ++i)
2285     {
2286       eps = eps / 2.;
2287     }
2288
2289   /* TODO we should make sure that 'float' matches kind 4 */
2290   match_float = gfc_real_kinds[k].kind == 4;
2291   if (mpfr_cmp_ui (s->value.real, 0) > 0)
2292     {
2293       if (match_float)
2294         {
2295           rval = (float) val;
2296           rval = rval + eps;
2297           mpfr_set_d (result->value.real, rval, GFC_RND_MODE);
2298         }
2299       else
2300         {
2301           val = val + eps;
2302           mpfr_set_d (result->value.real, val, GFC_RND_MODE);
2303         }
2304     }
2305   else if (mpfr_cmp_ui (s->value.real, 0) < 0)
2306     {
2307       if (match_float)
2308         {
2309           rval = (float) val;
2310           rval = rval - eps;
2311           mpfr_set_d (result->value.real, rval, GFC_RND_MODE);
2312         }
2313       else
2314         {
2315           val = val - eps;
2316           mpfr_set_d (result->value.real, val, GFC_RND_MODE);
2317         }
2318     }
2319   else
2320     {
2321       gfc_error ("Invalid second argument of NEAREST at %L", &s->where);
2322       gfc_free (result);
2323       return &gfc_bad_expr;
2324     }
2325
2326   return range_check (result, "NEAREST");
2327 }
2328
2329
2330 static gfc_expr *
2331 simplify_nint (const char *name, gfc_expr * e, gfc_expr * k)
2332 {
2333   gfc_expr *itrunc, *result;
2334   int kind;
2335
2336   kind = get_kind (BT_INTEGER, k, name, gfc_default_integer_kind);
2337   if (kind == -1)
2338     return &gfc_bad_expr;
2339
2340   if (e->expr_type != EXPR_CONSTANT)
2341     return NULL;
2342
2343   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
2344
2345   itrunc = gfc_copy_expr (e);
2346
2347   mpfr_round(itrunc->value.real, e->value.real);
2348
2349   gfc_mpfr_to_mpz (result->value.integer, itrunc->value.real);
2350
2351   gfc_free_expr (itrunc);
2352
2353   return range_check (result, name);
2354 }
2355
2356
2357 gfc_expr *
2358 gfc_simplify_nint (gfc_expr * e, gfc_expr * k)
2359 {
2360   return simplify_nint ("NINT", e, k);
2361 }
2362
2363
2364 gfc_expr *
2365 gfc_simplify_idnint (gfc_expr * e)
2366 {
2367   return simplify_nint ("IDNINT", e, NULL);
2368 }
2369
2370
2371 gfc_expr *
2372 gfc_simplify_not (gfc_expr * e)
2373 {
2374   gfc_expr *result;
2375   int i;
2376
2377   if (e->expr_type != EXPR_CONSTANT)
2378     return NULL;
2379
2380   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
2381
2382   mpz_com (result->value.integer, e->value.integer);
2383
2384   /* Because of how GMP handles numbers, the result must be ANDed with
2385      the max_int mask.  For radices <> 2, this will require change.  */
2386
2387   i = gfc_validate_kind (BT_INTEGER, e->ts.kind, false);
2388
2389   mpz_and (result->value.integer, result->value.integer,
2390            gfc_integer_kinds[i].max_int);
2391
2392   return range_check (result, "NOT");
2393 }
2394
2395
2396 gfc_expr *
2397 gfc_simplify_null (gfc_expr * mold)
2398 {
2399   gfc_expr *result;
2400
2401   result = gfc_get_expr ();
2402   result->expr_type = EXPR_NULL;
2403
2404   if (mold == NULL)
2405     result->ts.type = BT_UNKNOWN;
2406   else
2407     {
2408       result->ts = mold->ts;
2409       result->where = mold->where;
2410     }
2411
2412   return result;
2413 }
2414
2415
2416 gfc_expr *
2417 gfc_simplify_precision (gfc_expr * e)
2418 {
2419   gfc_expr *result;
2420   int i;
2421
2422   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
2423
2424   result = gfc_int_expr (gfc_real_kinds[i].precision);
2425   result->where = e->where;
2426
2427   return result;
2428 }
2429
2430
2431 gfc_expr *
2432 gfc_simplify_radix (gfc_expr * e)
2433 {
2434   gfc_expr *result;
2435   int i;
2436
2437   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
2438   switch (e->ts.type)
2439     {
2440     case BT_INTEGER:
2441       i = gfc_integer_kinds[i].radix;
2442       break;
2443
2444     case BT_REAL:
2445       i = gfc_real_kinds[i].radix;
2446       break;
2447
2448     default:
2449       gcc_unreachable ();
2450     }
2451
2452   result = gfc_int_expr (i);
2453   result->where = e->where;
2454
2455   return result;
2456 }
2457
2458
2459 gfc_expr *
2460 gfc_simplify_range (gfc_expr * e)
2461 {
2462   gfc_expr *result;
2463   int i;
2464   long j;
2465
2466   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
2467
2468   switch (e->ts.type)
2469     {
2470     case BT_INTEGER:
2471       j = gfc_integer_kinds[i].range;
2472       break;
2473
2474     case BT_REAL:
2475     case BT_COMPLEX:
2476       j = gfc_real_kinds[i].range;
2477       break;
2478
2479     default:
2480       gcc_unreachable ();
2481     }
2482
2483   result = gfc_int_expr (j);
2484   result->where = e->where;
2485
2486   return result;
2487 }
2488
2489
2490 gfc_expr *
2491 gfc_simplify_real (gfc_expr * e, gfc_expr * k)
2492 {
2493   gfc_expr *result;
2494   int kind;
2495
2496   if (e->ts.type == BT_COMPLEX)
2497     kind = get_kind (BT_REAL, k, "REAL", e->ts.kind);
2498   else
2499     kind = get_kind (BT_REAL, k, "REAL", gfc_default_real_kind);
2500
2501   if (kind == -1)
2502     return &gfc_bad_expr;
2503
2504   if (e->expr_type != EXPR_CONSTANT)
2505     return NULL;
2506
2507   switch (e->ts.type)
2508     {
2509     case BT_INTEGER:
2510       result = gfc_int2real (e, kind);
2511       break;
2512
2513     case BT_REAL:
2514       result = gfc_real2real (e, kind);
2515       break;
2516
2517     case BT_COMPLEX:
2518       result = gfc_complex2real (e, kind);
2519       break;
2520
2521     default:
2522       gfc_internal_error ("bad type in REAL");
2523       /* Not reached */
2524     }
2525
2526   return range_check (result, "REAL");
2527 }
2528
2529 gfc_expr *
2530 gfc_simplify_repeat (gfc_expr * e, gfc_expr * n)
2531 {
2532   gfc_expr *result;
2533   int i, j, len, ncopies, nlen;
2534
2535   if (e->expr_type != EXPR_CONSTANT || n->expr_type != EXPR_CONSTANT)
2536     return NULL;
2537
2538   if (n != NULL && (gfc_extract_int (n, &ncopies) != NULL || ncopies < 0))
2539     {
2540       gfc_error ("Invalid second argument of REPEAT at %L", &n->where);
2541       return &gfc_bad_expr;
2542     }
2543
2544   len = e->value.character.length;
2545   nlen = ncopies * len;
2546
2547   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
2548
2549   if (ncopies == 0)
2550     {
2551       result->value.character.string = gfc_getmem (1);
2552       result->value.character.length = 0;
2553       result->value.character.string[0] = '\0';
2554       return result;
2555     }
2556
2557   result->value.character.length = nlen;
2558   result->value.character.string = gfc_getmem (nlen + 1);
2559
2560   for (i = 0; i < ncopies; i++)
2561     for (j = 0; j < len; j++)
2562       result->value.character.string[j + i * len] =
2563         e->value.character.string[j];
2564
2565   result->value.character.string[nlen] = '\0';  /* For debugger */
2566   return result;
2567 }
2568
2569
2570 /* This one is a bear, but mainly has to do with shuffling elements.  */
2571
2572 gfc_expr *
2573 gfc_simplify_reshape (gfc_expr * source, gfc_expr * shape_exp,
2574                       gfc_expr * pad, gfc_expr * order_exp)
2575 {
2576
2577   int order[GFC_MAX_DIMENSIONS], shape[GFC_MAX_DIMENSIONS];
2578   int i, rank, npad, x[GFC_MAX_DIMENSIONS];
2579   gfc_constructor *head, *tail;
2580   mpz_t index, size;
2581   unsigned long j;
2582   size_t nsource;
2583   gfc_expr *e;
2584
2585   /* Unpack the shape array.  */
2586   if (source->expr_type != EXPR_ARRAY || !gfc_is_constant_expr (source))
2587     return NULL;
2588
2589   if (shape_exp->expr_type != EXPR_ARRAY || !gfc_is_constant_expr (shape_exp))
2590     return NULL;
2591
2592   if (pad != NULL
2593       && (pad->expr_type != EXPR_ARRAY
2594           || !gfc_is_constant_expr (pad)))
2595     return NULL;
2596
2597   if (order_exp != NULL
2598       && (order_exp->expr_type != EXPR_ARRAY
2599           || !gfc_is_constant_expr (order_exp)))
2600     return NULL;
2601
2602   mpz_init (index);
2603   rank = 0;
2604   head = tail = NULL;
2605
2606   for (;;)
2607     {
2608       e = gfc_get_array_element (shape_exp, rank);
2609       if (e == NULL)
2610         break;
2611
2612       if (gfc_extract_int (e, &shape[rank]) != NULL)
2613         {
2614           gfc_error ("Integer too large in shape specification at %L",
2615                      &e->where);
2616           gfc_free_expr (e);
2617           goto bad_reshape;
2618         }
2619
2620       gfc_free_expr (e);
2621
2622       if (rank >= GFC_MAX_DIMENSIONS)
2623         {
2624           gfc_error ("Too many dimensions in shape specification for RESHAPE "
2625                      "at %L", &e->where);
2626
2627           goto bad_reshape;
2628         }
2629
2630       if (shape[rank] < 0)
2631         {
2632           gfc_error ("Shape specification at %L cannot be negative",
2633                      &e->where);
2634           goto bad_reshape;
2635         }
2636
2637       rank++;
2638     }
2639
2640   if (rank == 0)
2641     {
2642       gfc_error ("Shape specification at %L cannot be the null array",
2643                  &shape_exp->where);
2644       goto bad_reshape;
2645     }
2646
2647   /* Now unpack the order array if present.  */
2648   if (order_exp == NULL)
2649     {
2650       for (i = 0; i < rank; i++)
2651         order[i] = i;
2652
2653     }
2654   else
2655     {
2656
2657       for (i = 0; i < rank; i++)
2658         x[i] = 0;
2659
2660       for (i = 0; i < rank; i++)
2661         {
2662           e = gfc_get_array_element (order_exp, i);
2663           if (e == NULL)
2664             {
2665               gfc_error
2666                 ("ORDER parameter of RESHAPE at %L is not the same size "
2667                  "as SHAPE parameter", &order_exp->where);
2668               goto bad_reshape;
2669             }
2670
2671           if (gfc_extract_int (e, &order[i]) != NULL)
2672             {
2673               gfc_error ("Error in ORDER parameter of RESHAPE at %L",
2674                          &e->where);
2675               gfc_free_expr (e);
2676               goto bad_reshape;
2677             }
2678
2679           gfc_free_expr (e);
2680
2681           if (order[i] < 1 || order[i] > rank)
2682             {
2683               gfc_error ("ORDER parameter of RESHAPE at %L is out of range",
2684                          &e->where);
2685               goto bad_reshape;
2686             }
2687
2688           order[i]--;
2689
2690           if (x[order[i]])
2691             {
2692               gfc_error ("Invalid permutation in ORDER parameter at %L",
2693                          &e->where);
2694               goto bad_reshape;
2695             }
2696
2697           x[order[i]] = 1;
2698         }
2699     }
2700
2701   /* Count the elements in the source and padding arrays.  */
2702
2703   npad = 0;
2704   if (pad != NULL)
2705     {
2706       gfc_array_size (pad, &size);
2707       npad = mpz_get_ui (size);
2708       mpz_clear (size);
2709     }
2710
2711   gfc_array_size (source, &size);
2712   nsource = mpz_get_ui (size);
2713   mpz_clear (size);
2714
2715   /* If it weren't for that pesky permutation we could just loop
2716      through the source and round out any shortage with pad elements.
2717      But no, someone just had to have the compiler do something the
2718      user should be doing.  */
2719
2720   for (i = 0; i < rank; i++)
2721     x[i] = 0;
2722
2723   for (;;)
2724     {
2725       /* Figure out which element to extract.  */
2726       mpz_set_ui (index, 0);
2727
2728       for (i = rank - 1; i >= 0; i--)
2729         {
2730           mpz_add_ui (index, index, x[order[i]]);
2731           if (i != 0)
2732             mpz_mul_ui (index, index, shape[order[i - 1]]);
2733         }
2734
2735       if (mpz_cmp_ui (index, INT_MAX) > 0)
2736         gfc_internal_error ("Reshaped array too large at %L", &e->where);
2737
2738       j = mpz_get_ui (index);
2739
2740       if (j < nsource)
2741         e = gfc_get_array_element (source, j);
2742       else
2743         {
2744           j = j - nsource;
2745
2746           if (npad == 0)
2747             {
2748               gfc_error
2749                 ("PAD parameter required for short SOURCE parameter at %L",
2750                  &source->where);
2751               goto bad_reshape;
2752             }
2753
2754           j = j % npad;
2755           e = gfc_get_array_element (pad, j);
2756         }
2757
2758       if (head == NULL)
2759         head = tail = gfc_get_constructor ();
2760       else
2761         {
2762           tail->next = gfc_get_constructor ();
2763           tail = tail->next;
2764         }
2765
2766       if (e == NULL)
2767         goto bad_reshape;
2768
2769       tail->where = e->where;
2770       tail->expr = e;
2771
2772       /* Calculate the next element.  */
2773       i = 0;
2774
2775 inc:
2776       if (++x[i] < shape[i])
2777         continue;
2778       x[i++] = 0;
2779       if (i < rank)
2780         goto inc;
2781
2782       break;
2783     }
2784
2785   mpz_clear (index);
2786
2787   e = gfc_get_expr ();
2788   e->where = source->where;
2789   e->expr_type = EXPR_ARRAY;
2790   e->value.constructor = head;
2791   e->shape = gfc_get_shape (rank);
2792
2793   for (i = 0; i < rank; i++)
2794     mpz_init_set_ui (e->shape[i], shape[i]);
2795
2796   e->ts = head->expr->ts;
2797   e->rank = rank;
2798
2799   return e;
2800
2801 bad_reshape:
2802   gfc_free_constructor (head);
2803   mpz_clear (index);
2804   return &gfc_bad_expr;
2805 }
2806
2807
2808 gfc_expr *
2809 gfc_simplify_rrspacing (gfc_expr * x)
2810 {
2811   gfc_expr *result;
2812   mpfr_t absv, log2, exp, frac, pow2;
2813   int i, p;
2814
2815   if (x->expr_type != EXPR_CONSTANT)
2816     return NULL;
2817
2818   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
2819
2820   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
2821
2822   p = gfc_real_kinds[i].digits;
2823
2824   gfc_set_model_kind (x->ts.kind);
2825
2826   if (mpfr_sgn (x->value.real) == 0)
2827     {
2828       mpfr_ui_div (result->value.real, 1, gfc_real_kinds[i].tiny, GFC_RND_MODE);
2829       return result;
2830     }
2831
2832   mpfr_init (log2);
2833   mpfr_init (absv);
2834   mpfr_init (frac);
2835   mpfr_init (pow2);
2836
2837   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
2838   mpfr_log2 (log2, absv, GFC_RND_MODE);
2839
2840   mpfr_trunc (log2, log2);
2841   mpfr_add_ui (exp, log2, 1, GFC_RND_MODE);
2842
2843   mpfr_ui_pow (pow2, 2, exp, GFC_RND_MODE);
2844   mpfr_div (frac, absv, pow2, GFC_RND_MODE);
2845
2846   mpfr_mul_2exp (result->value.real, frac, (unsigned long)p, GFC_RND_MODE);
2847
2848   mpfr_clear (log2);
2849   mpfr_clear (absv);
2850   mpfr_clear (frac);
2851   mpfr_clear (pow2);
2852
2853   return range_check (result, "RRSPACING");
2854 }
2855
2856
2857 gfc_expr *
2858 gfc_simplify_scale (gfc_expr * x, gfc_expr * i)
2859 {
2860   int k, neg_flag, power, exp_range;
2861   mpfr_t scale, radix;
2862   gfc_expr *result;
2863
2864   if (x->expr_type != EXPR_CONSTANT || i->expr_type != EXPR_CONSTANT)
2865     return NULL;
2866
2867   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
2868
2869   if (mpfr_sgn (x->value.real) == 0)
2870     {
2871       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
2872       return result;
2873     }
2874
2875   k = gfc_validate_kind (BT_REAL, x->ts.kind, false);
2876
2877   exp_range = gfc_real_kinds[k].max_exponent - gfc_real_kinds[k].min_exponent;
2878
2879   /* This check filters out values of i that would overflow an int.  */
2880   if (mpz_cmp_si (i->value.integer, exp_range + 2) > 0
2881       || mpz_cmp_si (i->value.integer, -exp_range - 2) < 0)
2882     {
2883       gfc_error ("Result of SCALE overflows its kind at %L", &result->where);
2884       return &gfc_bad_expr;
2885     }
2886
2887   /* Compute scale = radix ** power.  */
2888   power = mpz_get_si (i->value.integer);
2889
2890   if (power >= 0)
2891     neg_flag = 0;
2892   else
2893     {
2894       neg_flag = 1;
2895       power = -power;
2896     }
2897
2898   gfc_set_model_kind (x->ts.kind);
2899   mpfr_init (scale);
2900   mpfr_init (radix);
2901   mpfr_set_ui (radix, gfc_real_kinds[k].radix, GFC_RND_MODE);
2902   mpfr_pow_ui (scale, radix, power, GFC_RND_MODE);
2903
2904   if (neg_flag)
2905     mpfr_div (result->value.real, x->value.real, scale, GFC_RND_MODE);
2906   else
2907     mpfr_mul (result->value.real, x->value.real, scale, GFC_RND_MODE);
2908
2909   mpfr_clear (scale);
2910   mpfr_clear (radix);
2911
2912   return range_check (result, "SCALE");
2913 }
2914
2915
2916 gfc_expr *
2917 gfc_simplify_scan (gfc_expr * e, gfc_expr * c, gfc_expr * b)
2918 {
2919   gfc_expr *result;
2920   int back;
2921   size_t i;
2922   size_t indx, len, lenc;
2923
2924   if (e->expr_type != EXPR_CONSTANT || c->expr_type != EXPR_CONSTANT)
2925     return NULL;
2926
2927   if (b != NULL && b->value.logical != 0)
2928     back = 1;
2929   else
2930     back = 0;
2931
2932   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
2933                                 &e->where);
2934
2935   len = e->value.character.length;
2936   lenc = c->value.character.length;
2937
2938   if (len == 0 || lenc == 0)
2939     {
2940       indx = 0;
2941     }
2942   else
2943     {
2944       if (back == 0)
2945         {
2946           indx =
2947             strcspn (e->value.character.string, c->value.character.string) + 1;
2948           if (indx > len)
2949             indx = 0;
2950         }
2951       else
2952         {
2953           i = 0;
2954           for (indx = len; indx > 0; indx--)
2955             {
2956               for (i = 0; i < lenc; i++)
2957                 {
2958                   if (c->value.character.string[i]
2959                         == e->value.character.string[indx - 1])
2960                     break;
2961                 }
2962               if (i < lenc)
2963                 break;
2964             }
2965         }
2966     }
2967   mpz_set_ui (result->value.integer, indx);
2968   return range_check (result, "SCAN");
2969 }
2970
2971
2972 gfc_expr *
2973 gfc_simplify_selected_int_kind (gfc_expr * e)
2974 {
2975   int i, kind, range;
2976   gfc_expr *result;
2977
2978   if (e->expr_type != EXPR_CONSTANT || gfc_extract_int (e, &range) != NULL)
2979     return NULL;
2980
2981   kind = INT_MAX;
2982
2983   for (i = 0; gfc_integer_kinds[i].kind != 0; i++)
2984     if (gfc_integer_kinds[i].range >= range
2985         && gfc_integer_kinds[i].kind < kind)
2986       kind = gfc_integer_kinds[i].kind;
2987
2988   if (kind == INT_MAX)
2989     kind = -1;
2990
2991   result = gfc_int_expr (kind);
2992   result->where = e->where;
2993
2994   return result;
2995 }
2996
2997
2998 gfc_expr *
2999 gfc_simplify_selected_real_kind (gfc_expr * p, gfc_expr * q)
3000 {
3001   int range, precision, i, kind, found_precision, found_range;
3002   gfc_expr *result;
3003
3004   if (p == NULL)
3005     precision = 0;
3006   else
3007     {
3008       if (p->expr_type != EXPR_CONSTANT
3009           || gfc_extract_int (p, &precision) != NULL)
3010         return NULL;
3011     }
3012
3013   if (q == NULL)
3014     range = 0;
3015   else
3016     {
3017       if (q->expr_type != EXPR_CONSTANT
3018           || gfc_extract_int (q, &range) != NULL)
3019         return NULL;
3020     }
3021
3022   kind = INT_MAX;
3023   found_precision = 0;
3024   found_range = 0;
3025
3026   for (i = 0; gfc_real_kinds[i].kind != 0; i++)
3027     {
3028       if (gfc_real_kinds[i].precision >= precision)
3029         found_precision = 1;
3030
3031       if (gfc_real_kinds[i].range >= range)
3032         found_range = 1;
3033
3034       if (gfc_real_kinds[i].precision >= precision
3035           && gfc_real_kinds[i].range >= range && gfc_real_kinds[i].kind < kind)
3036         kind = gfc_real_kinds[i].kind;
3037     }
3038
3039   if (kind == INT_MAX)
3040     {
3041       kind = 0;
3042
3043       if (!found_precision)
3044         kind = -1;
3045       if (!found_range)
3046         kind -= 2;
3047     }
3048
3049   result = gfc_int_expr (kind);
3050   result->where = (p != NULL) ? p->where : q->where;
3051
3052   return result;
3053 }
3054
3055
3056 gfc_expr *
3057 gfc_simplify_set_exponent (gfc_expr * x, gfc_expr * i)
3058 {
3059   gfc_expr *result;
3060   mpfr_t exp, absv, log2, pow2, frac;
3061   unsigned long exp2;
3062
3063   if (x->expr_type != EXPR_CONSTANT || i->expr_type != EXPR_CONSTANT)
3064     return NULL;
3065
3066   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
3067
3068   gfc_set_model_kind (x->ts.kind);
3069
3070   if (mpfr_sgn (x->value.real) == 0)
3071     {
3072       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
3073       return result;
3074     }
3075
3076   mpfr_init (absv);
3077   mpfr_init (log2);
3078   mpfr_init (exp);
3079   mpfr_init (pow2);
3080   mpfr_init (frac);
3081
3082   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
3083   mpfr_log2 (log2, absv, GFC_RND_MODE);
3084
3085   mpfr_trunc (log2, log2);
3086   mpfr_add_ui (exp, log2, 1, GFC_RND_MODE);
3087
3088   /* Old exponent value, and fraction.  */
3089   mpfr_ui_pow (pow2, 2, exp, GFC_RND_MODE);
3090
3091   mpfr_div (frac, absv, pow2, GFC_RND_MODE);
3092
3093   /* New exponent.  */
3094   exp2 = (unsigned long) mpz_get_d (i->value.integer);
3095   mpfr_mul_2exp (result->value.real, frac, exp2, GFC_RND_MODE);
3096
3097   mpfr_clear (absv);
3098   mpfr_clear (log2);
3099   mpfr_clear (pow2);
3100   mpfr_clear (frac);
3101
3102   return range_check (result, "SET_EXPONENT");
3103 }
3104
3105
3106 gfc_expr *
3107 gfc_simplify_shape (gfc_expr * source)
3108 {
3109   mpz_t shape[GFC_MAX_DIMENSIONS];
3110   gfc_expr *result, *e, *f;
3111   gfc_array_ref *ar;
3112   int n;
3113   try t;
3114
3115   if (source->rank == 0 || source->expr_type != EXPR_VARIABLE)
3116     return NULL;
3117
3118   result = gfc_start_constructor (BT_INTEGER, gfc_default_integer_kind,
3119                                   &source->where);
3120
3121   ar = gfc_find_array_ref (source);
3122
3123   t = gfc_array_ref_shape (ar, shape);
3124
3125   for (n = 0; n < source->rank; n++)
3126     {
3127       e = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
3128                                &source->where);
3129
3130       if (t == SUCCESS)
3131         {
3132           mpz_set (e->value.integer, shape[n]);
3133           mpz_clear (shape[n]);
3134         }
3135       else
3136         {
3137           mpz_set_ui (e->value.integer, n + 1);
3138
3139           f = gfc_simplify_size (source, e);
3140           gfc_free_expr (e);
3141           if (f == NULL)
3142             {
3143               gfc_free_expr (result);
3144               return NULL;
3145             }
3146           else
3147             {
3148               e = f;
3149             }
3150         }
3151
3152       gfc_append_constructor (result, e);
3153     }
3154
3155   return result;
3156 }
3157
3158
3159 gfc_expr *
3160 gfc_simplify_size (gfc_expr * array, gfc_expr * dim)
3161 {
3162   mpz_t size;
3163   gfc_expr *result;
3164   int d;
3165
3166   if (dim == NULL)
3167     {
3168       if (gfc_array_size (array, &size) == FAILURE)
3169         return NULL;
3170     }
3171   else
3172     {
3173       if (dim->expr_type != EXPR_CONSTANT)
3174         return NULL;
3175
3176       d = mpz_get_ui (dim->value.integer) - 1;
3177       if (gfc_array_dimen_size (array, d, &size) == FAILURE)
3178         return NULL;
3179     }
3180
3181   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
3182                                 &array->where);
3183
3184   mpz_set (result->value.integer, size);
3185
3186   return result;
3187 }
3188
3189
3190 gfc_expr *
3191 gfc_simplify_sign (gfc_expr * x, gfc_expr * y)
3192 {
3193   gfc_expr *result;
3194
3195   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
3196     return NULL;
3197
3198   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3199
3200   switch (x->ts.type)
3201     {
3202     case BT_INTEGER:
3203       mpz_abs (result->value.integer, x->value.integer);
3204       if (mpz_sgn (y->value.integer) < 0)
3205         mpz_neg (result->value.integer, result->value.integer);
3206
3207       break;
3208
3209     case BT_REAL:
3210       /* TODO: Handle -0.0 and +0.0 correctly on machines that support
3211          it.  */
3212       mpfr_abs (result->value.real, x->value.real, GFC_RND_MODE);
3213       if (mpfr_sgn (y->value.real) < 0)
3214         mpfr_neg (result->value.real, result->value.real, GFC_RND_MODE);
3215
3216       break;
3217
3218     default:
3219       gfc_internal_error ("Bad type in gfc_simplify_sign");
3220     }
3221
3222   return result;
3223 }
3224
3225
3226 gfc_expr *
3227 gfc_simplify_sin (gfc_expr * x)
3228 {
3229   gfc_expr *result;
3230   mpfr_t xp, xq;
3231
3232   if (x->expr_type != EXPR_CONSTANT)
3233     return NULL;
3234
3235   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3236
3237   switch (x->ts.type)
3238     {
3239     case BT_REAL:
3240       mpfr_sin (result->value.real, x->value.real, GFC_RND_MODE);
3241       break;
3242
3243     case BT_COMPLEX:
3244       gfc_set_model (x->value.real);
3245       mpfr_init (xp);
3246       mpfr_init (xq);
3247
3248       mpfr_sin  (xp, x->value.complex.r, GFC_RND_MODE);
3249       mpfr_cosh (xq, x->value.complex.i, GFC_RND_MODE);
3250       mpfr_mul (result->value.complex.r, xp, xq, GFC_RND_MODE);
3251
3252       mpfr_cos  (xp, x->value.complex.r, GFC_RND_MODE);
3253       mpfr_sinh (xq, x->value.complex.i, GFC_RND_MODE);
3254       mpfr_mul (result->value.complex.i, xp, xq, GFC_RND_MODE);
3255
3256       mpfr_clear (xp);
3257       mpfr_clear (xq);
3258       break;
3259
3260     default:
3261       gfc_internal_error ("in gfc_simplify_sin(): Bad type");
3262     }
3263
3264   return range_check (result, "SIN");
3265 }
3266
3267
3268 gfc_expr *
3269 gfc_simplify_sinh (gfc_expr * x)
3270 {
3271   gfc_expr *result;
3272
3273   if (x->expr_type != EXPR_CONSTANT)
3274     return NULL;
3275
3276   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3277
3278   mpfr_sinh(result->value.real, x->value.real, GFC_RND_MODE);
3279
3280   return range_check (result, "SINH");
3281 }
3282
3283
3284 /* The argument is always a double precision real that is converted to
3285    single precision.  TODO: Rounding!  */
3286
3287 gfc_expr *
3288 gfc_simplify_sngl (gfc_expr * a)
3289 {
3290   gfc_expr *result;
3291
3292   if (a->expr_type != EXPR_CONSTANT)
3293     return NULL;
3294
3295   result = gfc_real2real (a, gfc_default_real_kind);
3296   return range_check (result, "SNGL");
3297 }
3298
3299
3300 gfc_expr *
3301 gfc_simplify_spacing (gfc_expr * x)
3302 {
3303   gfc_expr *result;
3304   mpfr_t absv, log2;
3305   long diff;
3306   int i, p;
3307
3308   if (x->expr_type != EXPR_CONSTANT)
3309     return NULL;
3310
3311   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
3312
3313   p = gfc_real_kinds[i].digits;
3314
3315   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
3316
3317   gfc_set_model_kind (x->ts.kind);
3318
3319   if (mpfr_sgn (x->value.real) == 0)
3320     {
3321       mpfr_set (result->value.real, gfc_real_kinds[i].tiny, GFC_RND_MODE);
3322       return result;
3323     }
3324
3325   mpfr_init (log2);
3326   mpfr_init (absv);
3327
3328   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
3329   mpfr_log2 (log2, absv, GFC_RND_MODE);
3330   mpfr_trunc (log2, log2);
3331
3332   mpfr_add_ui (log2, log2, 1, GFC_RND_MODE);
3333
3334   /* FIXME: We should be using mpfr_get_si here, but this function is
3335      not available with the version of mpfr distributed with gmp (as of
3336      2004-09-17). Replace once mpfr has been imported into the gcc cvs
3337      tree.  */
3338   diff = (long)mpfr_get_d (log2, GFC_RND_MODE) - (long)p;
3339   mpfr_set_ui (result->value.real, 1, GFC_RND_MODE);
3340   mpfr_mul_2si (result->value.real, result->value.real, diff, GFC_RND_MODE);
3341
3342   mpfr_clear (log2);
3343   mpfr_clear (absv);
3344
3345   if (mpfr_cmp (result->value.real, gfc_real_kinds[i].tiny) < 0)
3346     mpfr_set (result->value.real, gfc_real_kinds[i].tiny, GFC_RND_MODE);
3347
3348   return range_check (result, "SPACING");
3349 }
3350
3351
3352 gfc_expr *
3353 gfc_simplify_sqrt (gfc_expr * e)
3354 {
3355   gfc_expr *result;
3356   mpfr_t ac, ad, s, t, w;
3357
3358   if (e->expr_type != EXPR_CONSTANT)
3359     return NULL;
3360
3361   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
3362
3363   switch (e->ts.type)
3364     {
3365     case BT_REAL:
3366       if (mpfr_cmp_si (e->value.real, 0) < 0)
3367         goto negative_arg;
3368       mpfr_sqrt (result->value.real, e->value.real, GFC_RND_MODE);
3369
3370       break;
3371
3372     case BT_COMPLEX:
3373       /* Formula taken from Numerical Recipes to avoid over- and
3374          underflow.  */
3375
3376       gfc_set_model (e->value.real);
3377       mpfr_init (ac);
3378       mpfr_init (ad);
3379       mpfr_init (s);
3380       mpfr_init (t);
3381       mpfr_init (w);
3382
3383       if (mpfr_cmp_ui (e->value.complex.r, 0) == 0
3384           && mpfr_cmp_ui (e->value.complex.i, 0) == 0)
3385         {
3386
3387           mpfr_set_ui (result->value.complex.r, 0, GFC_RND_MODE);
3388           mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
3389           break;
3390         }
3391
3392       mpfr_abs (ac, e->value.complex.r, GFC_RND_MODE);
3393       mpfr_abs (ad, e->value.complex.i, GFC_RND_MODE);
3394
3395       if (mpfr_cmp (ac, ad) >= 0)
3396         {
3397           mpfr_div (t, e->value.complex.i, e->value.complex.r, GFC_RND_MODE);
3398           mpfr_mul (t, t, t, GFC_RND_MODE);
3399           mpfr_add_ui (t, t, 1, GFC_RND_MODE);
3400           mpfr_sqrt (t, t, GFC_RND_MODE);
3401           mpfr_add_ui (t, t, 1, GFC_RND_MODE);
3402           mpfr_div_ui (t, t, 2, GFC_RND_MODE);
3403           mpfr_sqrt (t, t, GFC_RND_MODE);
3404           mpfr_sqrt (s, ac, GFC_RND_MODE);
3405           mpfr_mul (w, s, t, GFC_RND_MODE);
3406         }
3407       else
3408         {
3409           mpfr_div (s, e->value.complex.r, e->value.complex.i, GFC_RND_MODE);
3410           mpfr_mul (t, s, s, GFC_RND_MODE);
3411           mpfr_add_ui (t, t, 1, GFC_RND_MODE);
3412           mpfr_sqrt (t, t, GFC_RND_MODE);
3413           mpfr_abs (s, s, GFC_RND_MODE);
3414           mpfr_add (t, t, s, GFC_RND_MODE);
3415           mpfr_div_ui (t, t, 2, GFC_RND_MODE);
3416           mpfr_sqrt (t, t, GFC_RND_MODE);
3417           mpfr_sqrt (s, ad, GFC_RND_MODE);
3418           mpfr_mul (w, s, t, GFC_RND_MODE);
3419         }
3420
3421       if (mpfr_cmp_ui (w, 0) != 0 && mpfr_cmp_ui (e->value.complex.r, 0) >= 0)
3422         {
3423           mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
3424           mpfr_div (result->value.complex.i, e->value.complex.i, t, GFC_RND_MODE);
3425           mpfr_set (result->value.complex.r, w, GFC_RND_MODE);
3426         }
3427       else if (mpfr_cmp_ui (w, 0) != 0
3428                && mpfr_cmp_ui (e->value.complex.r, 0) < 0
3429                && mpfr_cmp_ui (e->value.complex.i, 0) >= 0)
3430         {
3431           mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
3432           mpfr_div (result->value.complex.r, e->value.complex.i, t, GFC_RND_MODE);
3433           mpfr_set (result->value.complex.i, w, GFC_RND_MODE);
3434         }
3435       else if (mpfr_cmp_ui (w, 0) != 0
3436                && mpfr_cmp_ui (e->value.complex.r, 0) < 0
3437                && mpfr_cmp_ui (e->value.complex.i, 0) < 0)
3438         {
3439           mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
3440           mpfr_div (result->value.complex.r, ad, t, GFC_RND_MODE);
3441           mpfr_neg (w, w, GFC_RND_MODE);
3442           mpfr_set (result->value.complex.i, w, GFC_RND_MODE);
3443         }
3444       else
3445         gfc_internal_error ("invalid complex argument of SQRT at %L",
3446                             &e->where);
3447
3448       mpfr_clear (s);
3449       mpfr_clear (t);
3450       mpfr_clear (ac);
3451       mpfr_clear (ad);
3452       mpfr_clear (w);
3453
3454       break;
3455
3456     default:
3457       gfc_internal_error ("invalid argument of SQRT at %L", &e->where);
3458     }
3459
3460   return range_check (result, "SQRT");
3461
3462 negative_arg:
3463   gfc_free_expr (result);
3464   gfc_error ("Argument of SQRT at %L has a negative value", &e->where);
3465   return &gfc_bad_expr;
3466 }
3467
3468
3469 gfc_expr *
3470 gfc_simplify_tan (gfc_expr * x)
3471 {
3472   int i;
3473   gfc_expr *result;
3474
3475   if (x->expr_type != EXPR_CONSTANT)
3476     return NULL;
3477
3478   i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
3479
3480   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3481
3482   mpfr_tan (result->value.real, x->value.real, GFC_RND_MODE);
3483
3484   return range_check (result, "TAN");
3485 }
3486
3487
3488 gfc_expr *
3489 gfc_simplify_tanh (gfc_expr * x)
3490 {
3491   gfc_expr *result;
3492
3493   if (x->expr_type != EXPR_CONSTANT)
3494     return NULL;
3495
3496   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3497
3498   mpfr_tanh (result->value.real, x->value.real, GFC_RND_MODE);
3499
3500   return range_check (result, "TANH");
3501
3502 }
3503
3504
3505 gfc_expr *
3506 gfc_simplify_tiny (gfc_expr * e)
3507 {
3508   gfc_expr *result;
3509   int i;
3510
3511   i = gfc_validate_kind (BT_REAL, e->ts.kind, false);
3512
3513   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
3514   mpfr_set (result->value.real, gfc_real_kinds[i].tiny, GFC_RND_MODE);
3515
3516   return result;
3517 }
3518
3519
3520 gfc_expr *
3521 gfc_simplify_trim (gfc_expr * e)
3522 {
3523   gfc_expr *result;
3524   int count, i, len, lentrim;
3525
3526   if (e->expr_type != EXPR_CONSTANT)
3527     return NULL;
3528
3529   len = e->value.character.length;
3530
3531   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
3532
3533   for (count = 0, i = 1; i <= len; ++i)
3534     {
3535       if (e->value.character.string[len - i] == ' ')
3536         count++;
3537       else
3538         break;
3539     }
3540
3541   lentrim = len - count;
3542
3543   result->value.character.length = lentrim;
3544   result->value.character.string = gfc_getmem (lentrim + 1);
3545
3546   for (i = 0; i < lentrim; i++)
3547     result->value.character.string[i] = e->value.character.string[i];
3548
3549   result->value.character.string[lentrim] = '\0';       /* For debugger */
3550
3551   return result;
3552 }
3553
3554
3555 gfc_expr *
3556 gfc_simplify_ubound (gfc_expr * array, gfc_expr * dim)
3557 {
3558   return simplify_bound (array, dim, 1);
3559 }
3560
3561
3562 gfc_expr *
3563 gfc_simplify_verify (gfc_expr * s, gfc_expr * set, gfc_expr * b)
3564 {
3565   gfc_expr *result;
3566   int back;
3567   size_t index, len, lenset;
3568   size_t i;
3569
3570   if (s->expr_type != EXPR_CONSTANT || set->expr_type != EXPR_CONSTANT)
3571     return NULL;
3572
3573   if (b != NULL && b->value.logical != 0)
3574     back = 1;
3575   else
3576     back = 0;
3577
3578   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
3579                                 &s->where);
3580
3581   len = s->value.character.length;
3582   lenset = set->value.character.length;
3583
3584   if (len == 0)
3585     {
3586       mpz_set_ui (result->value.integer, 0);
3587       return result;
3588     }
3589
3590   if (back == 0)
3591     {
3592       if (lenset == 0)
3593         {
3594           mpz_set_ui (result->value.integer, len);
3595           return result;
3596         }
3597
3598       index =
3599         strspn (s->value.character.string, set->value.character.string) + 1;
3600       if (index > len)
3601         index = 0;
3602
3603     }
3604   else
3605     {
3606       if (lenset == 0)
3607         {
3608           mpz_set_ui (result->value.integer, 1);
3609           return result;
3610         }
3611       for (index = len; index > 0; index --)
3612         {
3613           for (i = 0; i < lenset; i++)
3614             {
3615               if (s->value.character.string[index - 1]
3616                     == set->value.character.string[i])
3617                 break;
3618             }
3619           if (i == lenset)
3620             break;
3621         }
3622     }
3623
3624   mpz_set_ui (result->value.integer, index);
3625   return result;
3626 }
3627
3628 /****************** Constant simplification *****************/
3629
3630 /* Master function to convert one constant to another.  While this is
3631    used as a simplification function, it requires the destination type
3632    and kind information which is supplied by a special case in
3633    do_simplify().  */
3634
3635 gfc_expr *
3636 gfc_convert_constant (gfc_expr * e, bt type, int kind)
3637 {
3638   gfc_expr *g, *result, *(*f) (gfc_expr *, int);
3639   gfc_constructor *head, *c, *tail = NULL;
3640
3641   switch (e->ts.type)
3642     {
3643     case BT_INTEGER:
3644       switch (type)
3645         {
3646         case BT_INTEGER:
3647           f = gfc_int2int;
3648           break;
3649         case BT_REAL:
3650           f = gfc_int2real;
3651           break;
3652         case BT_COMPLEX:
3653           f = gfc_int2complex;
3654           break;
3655         default:
3656           goto oops;
3657         }
3658       break;
3659
3660     case BT_REAL:
3661       switch (type)
3662         {
3663         case BT_INTEGER:
3664           f = gfc_real2int;
3665           break;
3666         case BT_REAL:
3667           f = gfc_real2real;
3668           break;
3669         case BT_COMPLEX:
3670           f = gfc_real2complex;
3671           break;
3672         default:
3673           goto oops;
3674         }
3675       break;
3676
3677     case BT_COMPLEX:
3678       switch (type)
3679         {
3680         case BT_INTEGER:
3681           f = gfc_complex2int;
3682           break;
3683         case BT_REAL:
3684           f = gfc_complex2real;
3685           break;
3686         case BT_COMPLEX:
3687           f = gfc_complex2complex;
3688           break;
3689
3690         default:
3691           goto oops;
3692         }
3693       break;
3694
3695     case BT_LOGICAL:
3696       if (type != BT_LOGICAL)
3697         goto oops;
3698       f = gfc_log2log;
3699       break;
3700
3701     default:
3702     oops:
3703       gfc_internal_error ("gfc_convert_constant(): Unexpected type");
3704     }
3705
3706   result = NULL;
3707
3708   switch (e->expr_type)
3709     {
3710     case EXPR_CONSTANT:
3711       result = f (e, kind);
3712       if (result == NULL)
3713         return &gfc_bad_expr;
3714       break;
3715
3716     case EXPR_ARRAY:
3717       if (!gfc_is_constant_expr (e))
3718         break;
3719
3720       head = NULL;
3721
3722       for (c = e->value.constructor; c; c = c->next)
3723         {
3724           if (head == NULL)
3725             head = tail = gfc_get_constructor ();
3726           else
3727             {
3728               tail->next = gfc_get_constructor ();
3729               tail = tail->next;
3730             }
3731
3732           tail->where = c->where;
3733
3734           if (c->iterator == NULL)
3735             tail->expr = f (c->expr, kind);
3736           else
3737             {
3738               g = gfc_convert_constant (c->expr, type, kind);
3739               if (g == &gfc_bad_expr)
3740                 return g;
3741               tail->expr = g;
3742             }
3743
3744           if (tail->expr == NULL)
3745             {
3746               gfc_free_constructor (head);
3747               return NULL;
3748             }
3749         }
3750
3751       result = gfc_get_expr ();
3752       result->ts.type = type;
3753       result->ts.kind = kind;
3754       result->expr_type = EXPR_ARRAY;
3755       result->value.constructor = head;
3756       result->shape = gfc_copy_shape (e->shape, e->rank);
3757       result->where = e->where;
3758       result->rank = e->rank;
3759       break;
3760
3761     default:
3762       break;
3763     }
3764
3765   return result;
3766 }
3767
3768
3769 /****************** Helper functions ***********************/
3770
3771 /* Given a collating table, create the inverse table.  */
3772
3773 static void
3774 invert_table (const int *table, int *xtable)
3775 {
3776   int i;
3777
3778   for (i = 0; i < 256; i++)
3779     xtable[i] = 0;
3780
3781   for (i = 0; i < 256; i++)
3782     xtable[table[i]] = i;
3783 }
3784
3785
3786 void
3787 gfc_simplify_init_1 (void)
3788 {
3789
3790   invert_table (ascii_table, xascii_table);
3791 }