OSDN Git Service

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