OSDN Git Service

2008-11-24 Jerry DeLisle <jvdelisle@gcc.gnu.org>
[pf3gnuchains/gcc-fork.git] / gcc / fortran / arith.c
1 /* Compiler arithmetic
2    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008
3    Free Software Foundation, Inc.
4    Contributed by Andy Vaught
5
6 This file is part of GCC.
7
8 GCC is free software; you can redistribute it and/or modify it under
9 the terms of the GNU General Public License as published by the Free
10 Software Foundation; either version 3, or (at your option) any later
11 version.
12
13 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
14 WARRANTY; without even the implied warranty of MERCHANTABILITY or
15 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
16 for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with GCC; see the file COPYING3.  If not see
20 <http://www.gnu.org/licenses/>.  */
21
22 /* Since target arithmetic must be done on the host, there has to
23    be some way of evaluating arithmetic expressions as the host
24    would evaluate them.  We use the GNU MP library and the MPFR
25    library to do arithmetic, and this file provides the interface.  */
26
27 #include "config.h"
28 #include "system.h"
29 #include "flags.h"
30 #include "gfortran.h"
31 #include "arith.h"
32 #include "target-memory.h"
33
34 /* MPFR does not have a direct replacement for mpz_set_f() from GMP.
35    It's easily implemented with a few calls though.  */
36
37 void
38 gfc_mpfr_to_mpz (mpz_t z, mpfr_t x, locus *where)
39 {
40   mp_exp_t e;
41
42   if (mpfr_inf_p (x) || mpfr_nan_p (x))
43     {
44       gfc_error ("Conversion of an Infinity or Not-a-Number at %L "
45                  "to INTEGER", where);
46       mpz_set_ui (z, 0);
47       return;
48     }
49
50   e = mpfr_get_z_exp (z, x);
51
52   if (e > 0)
53     mpz_mul_2exp (z, z, e);
54   else
55     mpz_tdiv_q_2exp (z, z, -e);
56 }
57
58
59 /* Set the model number precision by the requested KIND.  */
60
61 void
62 gfc_set_model_kind (int kind)
63 {
64   int index = gfc_validate_kind (BT_REAL, kind, false);
65   int base2prec;
66
67   base2prec = gfc_real_kinds[index].digits;
68   if (gfc_real_kinds[index].radix != 2)
69     base2prec *= gfc_real_kinds[index].radix / 2;
70   mpfr_set_default_prec (base2prec);
71 }
72
73
74 /* Set the model number precision from mpfr_t x.  */
75
76 void
77 gfc_set_model (mpfr_t x)
78 {
79   mpfr_set_default_prec (mpfr_get_prec (x));
80 }
81
82
83 /* Given an arithmetic error code, return a pointer to a string that
84    explains the error.  */
85
86 static const char *
87 gfc_arith_error (arith code)
88 {
89   const char *p;
90
91   switch (code)
92     {
93     case ARITH_OK:
94       p = _("Arithmetic OK at %L");
95       break;
96     case ARITH_OVERFLOW:
97       p = _("Arithmetic overflow at %L");
98       break;
99     case ARITH_UNDERFLOW:
100       p = _("Arithmetic underflow at %L");
101       break;
102     case ARITH_NAN:
103       p = _("Arithmetic NaN at %L");
104       break;
105     case ARITH_DIV0:
106       p = _("Division by zero at %L");
107       break;
108     case ARITH_INCOMMENSURATE:
109       p = _("Array operands are incommensurate at %L");
110       break;
111     case ARITH_ASYMMETRIC:
112       p =
113         _("Integer outside symmetric range implied by Standard Fortran at %L");
114       break;
115     default:
116       gfc_internal_error ("gfc_arith_error(): Bad error code");
117     }
118
119   return p;
120 }
121
122
123 /* Get things ready to do math.  */
124
125 void
126 gfc_arith_init_1 (void)
127 {
128   gfc_integer_info *int_info;
129   gfc_real_info *real_info;
130   mpfr_t a, b;
131   int i;
132
133   mpfr_set_default_prec (128);
134   mpfr_init (a);
135
136   /* Convert the minimum and maximum values for each kind into their
137      GNU MP representation.  */
138   for (int_info = gfc_integer_kinds; int_info->kind != 0; int_info++)
139     {
140       /* Huge  */
141       mpz_init (int_info->huge);
142       mpz_set_ui (int_info->huge, int_info->radix);
143       mpz_pow_ui (int_info->huge, int_info->huge, int_info->digits);
144       mpz_sub_ui (int_info->huge, int_info->huge, 1);
145
146       /* These are the numbers that are actually representable by the
147          target.  For bases other than two, this needs to be changed.  */
148       if (int_info->radix != 2)
149         gfc_internal_error ("Fix min_int calculation");
150
151       /* See PRs 13490 and 17912, related to integer ranges.
152          The pedantic_min_int exists for range checking when a program
153          is compiled with -pedantic, and reflects the belief that
154          Standard Fortran requires integers to be symmetrical, i.e.
155          every negative integer must have a representable positive
156          absolute value, and vice versa.  */
157
158       mpz_init (int_info->pedantic_min_int);
159       mpz_neg (int_info->pedantic_min_int, int_info->huge);
160
161       mpz_init (int_info->min_int);
162       mpz_sub_ui (int_info->min_int, int_info->pedantic_min_int, 1);
163
164       /* Range  */
165       mpfr_set_z (a, int_info->huge, GFC_RND_MODE);
166       mpfr_log10 (a, a, GFC_RND_MODE);
167       mpfr_trunc (a, a);
168       int_info->range = (int) mpfr_get_si (a, GFC_RND_MODE);
169     }
170
171   mpfr_clear (a);
172
173   for (real_info = gfc_real_kinds; real_info->kind != 0; real_info++)
174     {
175       gfc_set_model_kind (real_info->kind);
176
177       mpfr_init (a);
178       mpfr_init (b);
179
180       /* huge(x) = (1 - b**(-p)) * b**(emax-1) * b  */
181       /* 1 - b**(-p)  */
182       mpfr_init (real_info->huge);
183       mpfr_set_ui (real_info->huge, 1, GFC_RND_MODE);
184       mpfr_set_ui (a, real_info->radix, GFC_RND_MODE);
185       mpfr_pow_si (a, a, -real_info->digits, GFC_RND_MODE);
186       mpfr_sub (real_info->huge, real_info->huge, a, GFC_RND_MODE);
187
188       /* b**(emax-1)  */
189       mpfr_set_ui (a, real_info->radix, GFC_RND_MODE);
190       mpfr_pow_ui (a, a, real_info->max_exponent - 1, GFC_RND_MODE);
191
192       /* (1 - b**(-p)) * b**(emax-1)  */
193       mpfr_mul (real_info->huge, real_info->huge, a, GFC_RND_MODE);
194
195       /* (1 - b**(-p)) * b**(emax-1) * b  */
196       mpfr_mul_ui (real_info->huge, real_info->huge, real_info->radix,
197                    GFC_RND_MODE);
198
199       /* tiny(x) = b**(emin-1)  */
200       mpfr_init (real_info->tiny);
201       mpfr_set_ui (real_info->tiny, real_info->radix, GFC_RND_MODE);
202       mpfr_pow_si (real_info->tiny, real_info->tiny,
203                    real_info->min_exponent - 1, GFC_RND_MODE);
204
205       /* subnormal (x) = b**(emin - digit)  */
206       mpfr_init (real_info->subnormal);
207       mpfr_set_ui (real_info->subnormal, real_info->radix, GFC_RND_MODE);
208       mpfr_pow_si (real_info->subnormal, real_info->subnormal,
209                    real_info->min_exponent - real_info->digits, GFC_RND_MODE);
210
211       /* epsilon(x) = b**(1-p)  */
212       mpfr_init (real_info->epsilon);
213       mpfr_set_ui (real_info->epsilon, real_info->radix, GFC_RND_MODE);
214       mpfr_pow_si (real_info->epsilon, real_info->epsilon,
215                    1 - real_info->digits, GFC_RND_MODE);
216
217       /* range(x) = int(min(log10(huge(x)), -log10(tiny))  */
218       mpfr_log10 (a, real_info->huge, GFC_RND_MODE);
219       mpfr_log10 (b, real_info->tiny, GFC_RND_MODE);
220       mpfr_neg (b, b, GFC_RND_MODE);
221
222       /* a = min(a, b)  */
223       mpfr_min (a, a, b, GFC_RND_MODE);
224       mpfr_trunc (a, a);
225       real_info->range = (int) mpfr_get_si (a, GFC_RND_MODE);
226
227       /* precision(x) = int((p - 1) * log10(b)) + k  */
228       mpfr_set_ui (a, real_info->radix, GFC_RND_MODE);
229       mpfr_log10 (a, a, GFC_RND_MODE);
230       mpfr_mul_ui (a, a, real_info->digits - 1, GFC_RND_MODE);
231       mpfr_trunc (a, a);
232       real_info->precision = (int) mpfr_get_si (a, GFC_RND_MODE);
233
234       /* If the radix is an integral power of 10, add one to the precision.  */
235       for (i = 10; i <= real_info->radix; i *= 10)
236         if (i == real_info->radix)
237           real_info->precision++;
238
239       mpfr_clears (a, b, NULL);
240     }
241 }
242
243
244 /* Clean up, get rid of numeric constants.  */
245
246 void
247 gfc_arith_done_1 (void)
248 {
249   gfc_integer_info *ip;
250   gfc_real_info *rp;
251
252   for (ip = gfc_integer_kinds; ip->kind; ip++)
253     {
254       mpz_clear (ip->min_int);
255       mpz_clear (ip->pedantic_min_int);
256       mpz_clear (ip->huge);
257     }
258
259   for (rp = gfc_real_kinds; rp->kind; rp++)
260     mpfr_clears (rp->epsilon, rp->huge, rp->tiny, rp->subnormal, NULL);
261 }
262
263
264 /* Given a wide character value and a character kind, determine whether
265    the character is representable for that kind.  */
266 bool
267 gfc_check_character_range (gfc_char_t c, int kind)
268 {
269   /* As wide characters are stored as 32-bit values, they're all
270      representable in UCS=4.  */
271   if (kind == 4)
272     return true;
273
274   if (kind == 1)
275     return c <= 255 ? true : false;
276
277   gcc_unreachable ();
278 }
279
280
281 /* Given an integer and a kind, make sure that the integer lies within
282    the range of the kind.  Returns ARITH_OK, ARITH_ASYMMETRIC or
283    ARITH_OVERFLOW.  */
284
285 arith
286 gfc_check_integer_range (mpz_t p, int kind)
287 {
288   arith result;
289   int i;
290
291   i = gfc_validate_kind (BT_INTEGER, kind, false);
292   result = ARITH_OK;
293
294   if (pedantic)
295     {
296       if (mpz_cmp (p, gfc_integer_kinds[i].pedantic_min_int) < 0)
297         result = ARITH_ASYMMETRIC;
298     }
299
300
301   if (gfc_option.flag_range_check == 0)
302     return result;
303
304   if (mpz_cmp (p, gfc_integer_kinds[i].min_int) < 0
305       || mpz_cmp (p, gfc_integer_kinds[i].huge) > 0)
306     result = ARITH_OVERFLOW;
307
308   return result;
309 }
310
311
312 /* Given a real and a kind, make sure that the real lies within the
313    range of the kind.  Returns ARITH_OK, ARITH_OVERFLOW or
314    ARITH_UNDERFLOW.  */
315
316 static arith
317 gfc_check_real_range (mpfr_t p, int kind)
318 {
319   arith retval;
320   mpfr_t q;
321   int i;
322
323   i = gfc_validate_kind (BT_REAL, kind, false);
324
325   gfc_set_model (p);
326   mpfr_init (q);
327   mpfr_abs (q, p, GFC_RND_MODE);
328
329   retval = ARITH_OK;
330
331   if (mpfr_inf_p (p))
332     {
333       if (gfc_option.flag_range_check != 0)
334         retval = ARITH_OVERFLOW;
335     }
336   else if (mpfr_nan_p (p))
337     {
338       if (gfc_option.flag_range_check != 0)
339         retval = ARITH_NAN;
340     }
341   else if (mpfr_sgn (q) == 0)
342     {
343       mpfr_clear (q);
344       return retval;
345     }
346   else if (mpfr_cmp (q, gfc_real_kinds[i].huge) > 0)
347     {
348       if (gfc_option.flag_range_check == 0)
349         mpfr_set_inf (p, mpfr_sgn (p));
350       else
351         retval = ARITH_OVERFLOW;
352     }
353   else if (mpfr_cmp (q, gfc_real_kinds[i].subnormal) < 0)
354     {
355       if (gfc_option.flag_range_check == 0)
356         {
357           if (mpfr_sgn (p) < 0)
358             {
359               mpfr_set_ui (p, 0, GFC_RND_MODE);
360               mpfr_set_si (q, -1, GFC_RND_MODE);
361               mpfr_copysign (p, p, q, GFC_RND_MODE);
362             }
363           else
364             mpfr_set_ui (p, 0, GFC_RND_MODE);
365         }
366       else
367         retval = ARITH_UNDERFLOW;
368     }
369   else if (mpfr_cmp (q, gfc_real_kinds[i].tiny) < 0)
370     {
371       mp_exp_t emin, emax;
372       int en;
373
374       /* Save current values of emin and emax.  */
375       emin = mpfr_get_emin ();
376       emax = mpfr_get_emax ();
377
378       /* Set emin and emax for the current model number.  */
379       en = gfc_real_kinds[i].min_exponent - gfc_real_kinds[i].digits + 1;
380       mpfr_set_emin ((mp_exp_t) en);
381       mpfr_set_emax ((mp_exp_t) gfc_real_kinds[i].max_exponent);
382       mpfr_check_range (q, 0, GFC_RND_MODE);
383       mpfr_subnormalize (q, 0, GFC_RND_MODE);
384
385       /* Reset emin and emax.  */
386       mpfr_set_emin (emin);
387       mpfr_set_emax (emax);
388
389       /* Copy sign if needed.  */
390       if (mpfr_sgn (p) < 0)
391         mpfr_neg (p, q, GMP_RNDN);
392       else
393         mpfr_set (p, q, GMP_RNDN);
394     }
395
396   mpfr_clear (q);
397
398   return retval;
399 }
400
401
402 /* Function to return a constant expression node of a given type and kind.  */
403
404 gfc_expr *
405 gfc_constant_result (bt type, int kind, locus *where)
406 {
407   gfc_expr *result;
408
409   if (!where)
410     gfc_internal_error ("gfc_constant_result(): locus 'where' cannot be NULL");
411
412   result = gfc_get_expr ();
413
414   result->expr_type = EXPR_CONSTANT;
415   result->ts.type = type;
416   result->ts.kind = kind;
417   result->where = *where;
418
419   switch (type)
420     {
421     case BT_INTEGER:
422       mpz_init (result->value.integer);
423       break;
424
425     case BT_REAL:
426       gfc_set_model_kind (kind);
427       mpfr_init (result->value.real);
428       break;
429
430     case BT_COMPLEX:
431       gfc_set_model_kind (kind);
432       mpfr_init (result->value.complex.r);
433       mpfr_init (result->value.complex.i);
434       break;
435
436     default:
437       break;
438     }
439
440   return result;
441 }
442
443
444 /* Low-level arithmetic functions.  All of these subroutines assume
445    that all operands are of the same type and return an operand of the
446    same type.  The other thing about these subroutines is that they
447    can fail in various ways -- overflow, underflow, division by zero,
448    zero raised to the zero, etc.  */
449
450 static arith
451 gfc_arith_not (gfc_expr *op1, gfc_expr **resultp)
452 {
453   gfc_expr *result;
454
455   result = gfc_constant_result (BT_LOGICAL, op1->ts.kind, &op1->where);
456   result->value.logical = !op1->value.logical;
457   *resultp = result;
458
459   return ARITH_OK;
460 }
461
462
463 static arith
464 gfc_arith_and (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
465 {
466   gfc_expr *result;
467
468   result = gfc_constant_result (BT_LOGICAL, gfc_kind_max (op1, op2),
469                                 &op1->where);
470   result->value.logical = op1->value.logical && op2->value.logical;
471   *resultp = result;
472
473   return ARITH_OK;
474 }
475
476
477 static arith
478 gfc_arith_or (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
479 {
480   gfc_expr *result;
481
482   result = gfc_constant_result (BT_LOGICAL, gfc_kind_max (op1, op2),
483                                 &op1->where);
484   result->value.logical = op1->value.logical || op2->value.logical;
485   *resultp = result;
486
487   return ARITH_OK;
488 }
489
490
491 static arith
492 gfc_arith_eqv (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
493 {
494   gfc_expr *result;
495
496   result = gfc_constant_result (BT_LOGICAL, gfc_kind_max (op1, op2),
497                                 &op1->where);
498   result->value.logical = op1->value.logical == op2->value.logical;
499   *resultp = result;
500
501   return ARITH_OK;
502 }
503
504
505 static arith
506 gfc_arith_neqv (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
507 {
508   gfc_expr *result;
509
510   result = gfc_constant_result (BT_LOGICAL, gfc_kind_max (op1, op2),
511                                 &op1->where);
512   result->value.logical = op1->value.logical != op2->value.logical;
513   *resultp = result;
514
515   return ARITH_OK;
516 }
517
518
519 /* Make sure a constant numeric expression is within the range for
520    its type and kind.  Note that there's also a gfc_check_range(),
521    but that one deals with the intrinsic RANGE function.  */
522
523 arith
524 gfc_range_check (gfc_expr *e)
525 {
526   arith rc;
527   arith rc2;
528
529   switch (e->ts.type)
530     {
531     case BT_INTEGER:
532       rc = gfc_check_integer_range (e->value.integer, e->ts.kind);
533       break;
534
535     case BT_REAL:
536       rc = gfc_check_real_range (e->value.real, e->ts.kind);
537       if (rc == ARITH_UNDERFLOW)
538         mpfr_set_ui (e->value.real, 0, GFC_RND_MODE);
539       if (rc == ARITH_OVERFLOW)
540         mpfr_set_inf (e->value.real, mpfr_sgn (e->value.real));
541       if (rc == ARITH_NAN)
542         mpfr_set_nan (e->value.real);
543       break;
544
545     case BT_COMPLEX:
546       rc = gfc_check_real_range (e->value.complex.r, e->ts.kind);
547       if (rc == ARITH_UNDERFLOW)
548         mpfr_set_ui (e->value.complex.r, 0, GFC_RND_MODE);
549       if (rc == ARITH_OVERFLOW)
550         mpfr_set_inf (e->value.complex.r, mpfr_sgn (e->value.complex.r));
551       if (rc == ARITH_NAN)
552         mpfr_set_nan (e->value.complex.r);
553
554       rc2 = gfc_check_real_range (e->value.complex.i, e->ts.kind);
555       if (rc == ARITH_UNDERFLOW)
556         mpfr_set_ui (e->value.complex.i, 0, GFC_RND_MODE);
557       if (rc == ARITH_OVERFLOW)
558         mpfr_set_inf (e->value.complex.i, mpfr_sgn (e->value.complex.i));
559       if (rc == ARITH_NAN)
560         mpfr_set_nan (e->value.complex.i);
561
562       if (rc == ARITH_OK)
563         rc = rc2;
564       break;
565
566     default:
567       gfc_internal_error ("gfc_range_check(): Bad type");
568     }
569
570   return rc;
571 }
572
573
574 /* Several of the following routines use the same set of statements to
575    check the validity of the result.  Encapsulate the checking here.  */
576
577 static arith
578 check_result (arith rc, gfc_expr *x, gfc_expr *r, gfc_expr **rp)
579 {
580   arith val = rc;
581
582   if (val == ARITH_UNDERFLOW)
583     {
584       if (gfc_option.warn_underflow)
585         gfc_warning (gfc_arith_error (val), &x->where);
586       val = ARITH_OK;
587     }
588
589   if (val == ARITH_ASYMMETRIC)
590     {
591       gfc_warning (gfc_arith_error (val), &x->where);
592       val = ARITH_OK;
593     }
594
595   if (val != ARITH_OK)
596     gfc_free_expr (r);
597   else
598     *rp = r;
599
600   return val;
601 }
602
603
604 /* It may seem silly to have a subroutine that actually computes the
605    unary plus of a constant, but it prevents us from making exceptions
606    in the code elsewhere.  Used for unary plus and parenthesized
607    expressions.  */
608
609 static arith
610 gfc_arith_identity (gfc_expr *op1, gfc_expr **resultp)
611 {
612   *resultp = gfc_copy_expr (op1);
613   return ARITH_OK;
614 }
615
616
617 static arith
618 gfc_arith_uminus (gfc_expr *op1, gfc_expr **resultp)
619 {
620   gfc_expr *result;
621   arith rc;
622
623   result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
624
625   switch (op1->ts.type)
626     {
627     case BT_INTEGER:
628       mpz_neg (result->value.integer, op1->value.integer);
629       break;
630
631     case BT_REAL:
632       mpfr_neg (result->value.real, op1->value.real, GFC_RND_MODE);
633       break;
634
635     case BT_COMPLEX:
636       mpfr_neg (result->value.complex.r, op1->value.complex.r, GFC_RND_MODE);
637       mpfr_neg (result->value.complex.i, op1->value.complex.i, GFC_RND_MODE);
638       break;
639
640     default:
641       gfc_internal_error ("gfc_arith_uminus(): Bad basic type");
642     }
643
644   rc = gfc_range_check (result);
645
646   return check_result (rc, op1, result, resultp);
647 }
648
649
650 static arith
651 gfc_arith_plus (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
652 {
653   gfc_expr *result;
654   arith rc;
655
656   result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
657
658   switch (op1->ts.type)
659     {
660     case BT_INTEGER:
661       mpz_add (result->value.integer, op1->value.integer, op2->value.integer);
662       break;
663
664     case BT_REAL:
665       mpfr_add (result->value.real, op1->value.real, op2->value.real,
666                GFC_RND_MODE);
667       break;
668
669     case BT_COMPLEX:
670       mpfr_add (result->value.complex.r, op1->value.complex.r,
671                 op2->value.complex.r, GFC_RND_MODE);
672
673       mpfr_add (result->value.complex.i, op1->value.complex.i,
674                 op2->value.complex.i, GFC_RND_MODE);
675       break;
676
677     default:
678       gfc_internal_error ("gfc_arith_plus(): Bad basic type");
679     }
680
681   rc = gfc_range_check (result);
682
683   return check_result (rc, op1, result, resultp);
684 }
685
686
687 static arith
688 gfc_arith_minus (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
689 {
690   gfc_expr *result;
691   arith rc;
692
693   result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
694
695   switch (op1->ts.type)
696     {
697     case BT_INTEGER:
698       mpz_sub (result->value.integer, op1->value.integer, op2->value.integer);
699       break;
700
701     case BT_REAL:
702       mpfr_sub (result->value.real, op1->value.real, op2->value.real,
703                 GFC_RND_MODE);
704       break;
705
706     case BT_COMPLEX:
707       mpfr_sub (result->value.complex.r, op1->value.complex.r,
708                 op2->value.complex.r, GFC_RND_MODE);
709
710       mpfr_sub (result->value.complex.i, op1->value.complex.i,
711                 op2->value.complex.i, GFC_RND_MODE);
712       break;
713
714     default:
715       gfc_internal_error ("gfc_arith_minus(): Bad basic type");
716     }
717
718   rc = gfc_range_check (result);
719
720   return check_result (rc, op1, result, resultp);
721 }
722
723
724 static arith
725 gfc_arith_times (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
726 {
727   gfc_expr *result;
728   mpfr_t x, y;
729   arith rc;
730
731   result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
732
733   switch (op1->ts.type)
734     {
735     case BT_INTEGER:
736       mpz_mul (result->value.integer, op1->value.integer, op2->value.integer);
737       break;
738
739     case BT_REAL:
740       mpfr_mul (result->value.real, op1->value.real, op2->value.real,
741                GFC_RND_MODE);
742       break;
743
744     case BT_COMPLEX:
745       gfc_set_model (op1->value.complex.r);
746       mpfr_init (x);
747       mpfr_init (y);
748
749       mpfr_mul (x, op1->value.complex.r, op2->value.complex.r, GFC_RND_MODE);
750       mpfr_mul (y, op1->value.complex.i, op2->value.complex.i, GFC_RND_MODE);
751       mpfr_sub (result->value.complex.r, x, y, GFC_RND_MODE);
752
753       mpfr_mul (x, op1->value.complex.r, op2->value.complex.i, GFC_RND_MODE);
754       mpfr_mul (y, op1->value.complex.i, op2->value.complex.r, GFC_RND_MODE);
755       mpfr_add (result->value.complex.i, x, y, GFC_RND_MODE);
756
757       mpfr_clears (x, y, NULL);
758       break;
759
760     default:
761       gfc_internal_error ("gfc_arith_times(): Bad basic type");
762     }
763
764   rc = gfc_range_check (result);
765
766   return check_result (rc, op1, result, resultp);
767 }
768
769
770 static arith
771 gfc_arith_divide (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
772 {
773   gfc_expr *result;
774   mpfr_t x, y, div;
775   arith rc;
776
777   rc = ARITH_OK;
778
779   result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
780
781   switch (op1->ts.type)
782     {
783     case BT_INTEGER:
784       if (mpz_sgn (op2->value.integer) == 0)
785         {
786           rc = ARITH_DIV0;
787           break;
788         }
789
790       mpz_tdiv_q (result->value.integer, op1->value.integer,
791                   op2->value.integer);
792       break;
793
794     case BT_REAL:
795       if (mpfr_sgn (op2->value.real) == 0 && gfc_option.flag_range_check == 1)
796         {
797           rc = ARITH_DIV0;
798           break;
799         }
800
801       mpfr_div (result->value.real, op1->value.real, op2->value.real,
802                GFC_RND_MODE);
803       break;
804
805     case BT_COMPLEX:
806       if (mpfr_sgn (op2->value.complex.r) == 0
807           && mpfr_sgn (op2->value.complex.i) == 0
808           && gfc_option.flag_range_check == 1)
809         {
810           rc = ARITH_DIV0;
811           break;
812         }
813
814       gfc_set_model (op1->value.complex.r);
815       mpfr_init (x);
816       mpfr_init (y);
817       mpfr_init (div);
818
819       mpfr_mul (x, op2->value.complex.r, op2->value.complex.r, GFC_RND_MODE);
820       mpfr_mul (y, op2->value.complex.i, op2->value.complex.i, GFC_RND_MODE);
821       mpfr_add (div, x, y, GFC_RND_MODE);
822
823       mpfr_mul (x, op1->value.complex.r, op2->value.complex.r, GFC_RND_MODE);
824       mpfr_mul (y, op1->value.complex.i, op2->value.complex.i, GFC_RND_MODE);
825       mpfr_add (result->value.complex.r, x, y, GFC_RND_MODE);
826       mpfr_div (result->value.complex.r, result->value.complex.r, div,
827                 GFC_RND_MODE);
828
829       mpfr_mul (x, op1->value.complex.i, op2->value.complex.r, GFC_RND_MODE);
830       mpfr_mul (y, op1->value.complex.r, op2->value.complex.i, GFC_RND_MODE);
831       mpfr_sub (result->value.complex.i, x, y, GFC_RND_MODE);
832       mpfr_div (result->value.complex.i, result->value.complex.i, div,
833                 GFC_RND_MODE);
834
835       mpfr_clears (x, y, div, NULL);
836       break;
837
838     default:
839       gfc_internal_error ("gfc_arith_divide(): Bad basic type");
840     }
841
842   if (rc == ARITH_OK)
843     rc = gfc_range_check (result);
844
845   return check_result (rc, op1, result, resultp);
846 }
847
848
849 /* Compute the reciprocal of a complex number (guaranteed nonzero).  */
850
851 static void
852 complex_reciprocal (gfc_expr *op)
853 {
854   mpfr_t mod, tmp;
855
856   gfc_set_model (op->value.complex.r);
857   mpfr_init (mod);
858   mpfr_init (tmp);
859
860   mpfr_mul (mod, op->value.complex.r, op->value.complex.r, GFC_RND_MODE);
861   mpfr_mul (tmp, op->value.complex.i, op->value.complex.i, GFC_RND_MODE);
862   mpfr_add (mod, mod, tmp, GFC_RND_MODE);
863
864   mpfr_div (op->value.complex.r, op->value.complex.r, mod, GFC_RND_MODE);
865
866   mpfr_neg (op->value.complex.i, op->value.complex.i, GFC_RND_MODE);
867   mpfr_div (op->value.complex.i, op->value.complex.i, mod, GFC_RND_MODE);
868
869   mpfr_clears (tmp, mod, NULL);
870 }
871
872
873 /* Raise a complex number to positive power (power > 0).
874    This function will modify the content of power.
875
876    Use Binary Method, which is not an optimal but a simple and reasonable
877    arithmetic. See section 4.6.3, "Evaluation of Powers" of Donald E. Knuth,
878    "Seminumerical Algorithms", Vol. 2, "The Art of Computer Programming",
879    3rd Edition, 1998.  */
880
881 static void
882 complex_pow (gfc_expr *result, gfc_expr *base, mpz_t power)
883 {
884   mpfr_t x_r, x_i, tmp, re, im;
885
886   gfc_set_model (base->value.complex.r);
887   mpfr_init (x_r);
888   mpfr_init (x_i);
889   mpfr_init (tmp);
890   mpfr_init (re);
891   mpfr_init (im);
892
893   /* res = 1 */
894   mpfr_set_ui (result->value.complex.r, 1, GFC_RND_MODE);
895   mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
896
897   /* x = base */
898   mpfr_set (x_r, base->value.complex.r, GFC_RND_MODE);
899   mpfr_set (x_i, base->value.complex.i, GFC_RND_MODE);
900
901   /* Macro for complex multiplication. We have to take care that
902      res_r/res_i and a_r/a_i can (and will) be the same variable.  */
903 #define CMULT(res_r,res_i,a_r,a_i,b_r,b_i) \
904     mpfr_mul (re, a_r, b_r, GFC_RND_MODE), \
905     mpfr_mul (tmp, a_i, b_i, GFC_RND_MODE), \
906     mpfr_sub (re, re, tmp, GFC_RND_MODE), \
907     \
908     mpfr_mul (im, a_r, b_i, GFC_RND_MODE), \
909     mpfr_mul (tmp, a_i, b_r, GFC_RND_MODE), \
910     mpfr_add (res_i, im, tmp, GFC_RND_MODE), \
911     mpfr_set (res_r, re, GFC_RND_MODE)
912   
913 #define res_r result->value.complex.r
914 #define res_i result->value.complex.i
915
916   /* for (; power > 0; x *= x) */
917   for (; mpz_cmp_si (power, 0) > 0; CMULT(x_r,x_i,x_r,x_i,x_r,x_i))
918     {
919       /* if (power & 1) res = res * x; */
920       if (mpz_congruent_ui_p (power, 1, 2))
921         CMULT(res_r,res_i,res_r,res_i,x_r,x_i);
922
923       /* power /= 2; */
924       mpz_fdiv_q_ui (power, power, 2);
925     }
926
927 #undef res_r
928 #undef res_i
929 #undef CMULT
930
931   mpfr_clears (x_r, x_i, tmp, re, im, NULL);
932 }
933
934
935 /* Raise a number to an integer power.  */
936
937 static arith
938 gfc_arith_power (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
939 {
940   int power_sign;
941   gfc_expr *result;
942   arith rc;
943
944   gcc_assert (op2->expr_type == EXPR_CONSTANT && op2->ts.type == BT_INTEGER);
945
946   rc = ARITH_OK;
947   result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
948   power_sign = mpz_sgn (op2->value.integer);
949
950   if (power_sign == 0)
951     {
952       /* Handle something to the zeroth power.  Since we're dealing
953          with integral exponents, there is no ambiguity in the
954          limiting procedure used to determine the value of 0**0.  */
955       switch (op1->ts.type)
956         {
957         case BT_INTEGER:
958           mpz_set_ui (result->value.integer, 1);
959           break;
960
961         case BT_REAL:
962           mpfr_set_ui (result->value.real, 1, GFC_RND_MODE);
963           break;
964
965         case BT_COMPLEX:
966           mpfr_set_ui (result->value.complex.r, 1, GFC_RND_MODE);
967           mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
968           break;
969
970         default:
971           gfc_internal_error ("gfc_arith_power(): Bad base");
972         }
973     }
974   else
975     {
976       switch (op1->ts.type)
977         {
978         case BT_INTEGER:
979           {
980             int power;
981
982             /* First, we simplify the cases of op1 == 1, 0 or -1.  */
983             if (mpz_cmp_si (op1->value.integer, 1) == 0)
984               {
985                 /* 1**op2 == 1 */
986                 mpz_set_si (result->value.integer, 1);
987               }
988             else if (mpz_cmp_si (op1->value.integer, 0) == 0)
989               {
990                 /* 0**op2 == 0, if op2 > 0
991                    0**op2 overflow, if op2 < 0 ; in that case, we
992                    set the result to 0 and return ARITH_DIV0.  */
993                 mpz_set_si (result->value.integer, 0);
994                 if (mpz_cmp_si (op2->value.integer, 0) < 0)
995                   rc = ARITH_DIV0;
996               }
997             else if (mpz_cmp_si (op1->value.integer, -1) == 0)
998               {
999                 /* (-1)**op2 == (-1)**(mod(op2,2)) */
1000                 unsigned int odd = mpz_fdiv_ui (op2->value.integer, 2);
1001                 if (odd)
1002                   mpz_set_si (result->value.integer, -1);
1003                 else
1004                   mpz_set_si (result->value.integer, 1);
1005               }
1006             /* Then, we take care of op2 < 0.  */
1007             else if (mpz_cmp_si (op2->value.integer, 0) < 0)
1008               {
1009                 /* if op2 < 0, op1**op2 == 0  because abs(op1) > 1.  */
1010                 mpz_set_si (result->value.integer, 0);
1011               }
1012             else if (gfc_extract_int (op2, &power) != NULL)
1013               {
1014                 /* If op2 doesn't fit in an int, the exponentiation will
1015                    overflow, because op2 > 0 and abs(op1) > 1.  */
1016                 mpz_t max;
1017                 int i = gfc_validate_kind (BT_INTEGER, result->ts.kind, false);
1018
1019                 if (gfc_option.flag_range_check)
1020                   rc = ARITH_OVERFLOW;
1021
1022                 /* Still, we want to give the same value as the processor.  */
1023                 mpz_init (max);
1024                 mpz_add_ui (max, gfc_integer_kinds[i].huge, 1);
1025                 mpz_mul_ui (max, max, 2);
1026                 mpz_powm (result->value.integer, op1->value.integer,
1027                           op2->value.integer, max);
1028                 mpz_clear (max);
1029               }
1030             else
1031               mpz_pow_ui (result->value.integer, op1->value.integer, power);
1032           }
1033           break;
1034
1035         case BT_REAL:
1036           mpfr_pow_z (result->value.real, op1->value.real, op2->value.integer,
1037                       GFC_RND_MODE);
1038           break;
1039
1040         case BT_COMPLEX:
1041           {
1042             mpz_t apower;
1043
1044             /* Compute op1**abs(op2)  */
1045             mpz_init (apower);
1046             mpz_abs (apower, op2->value.integer);
1047             complex_pow (result, op1, apower);
1048             mpz_clear (apower);
1049
1050             /* If (op2 < 0), compute the inverse.  */
1051             if (power_sign < 0)
1052               complex_reciprocal (result);
1053
1054             break;
1055           }
1056
1057         default:
1058           break;
1059         }
1060     }
1061
1062   if (rc == ARITH_OK)
1063     rc = gfc_range_check (result);
1064
1065   return check_result (rc, op1, result, resultp);
1066 }
1067
1068
1069 /* Concatenate two string constants.  */
1070
1071 static arith
1072 gfc_arith_concat (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1073 {
1074   gfc_expr *result;
1075   int len;
1076
1077   gcc_assert (op1->ts.kind == op2->ts.kind);
1078   result = gfc_constant_result (BT_CHARACTER, op1->ts.kind,
1079                                 &op1->where);
1080
1081   len = op1->value.character.length + op2->value.character.length;
1082
1083   result->value.character.string = gfc_get_wide_string (len + 1);
1084   result->value.character.length = len;
1085
1086   memcpy (result->value.character.string, op1->value.character.string,
1087           op1->value.character.length * sizeof (gfc_char_t));
1088
1089   memcpy (&result->value.character.string[op1->value.character.length],
1090           op2->value.character.string,
1091           op2->value.character.length * sizeof (gfc_char_t));
1092
1093   result->value.character.string[len] = '\0';
1094
1095   *resultp = result;
1096
1097   return ARITH_OK;
1098 }
1099
1100 /* Comparison between real values; returns 0 if (op1 .op. op2) is true.
1101    This function mimics mpfr_cmp but takes NaN into account.  */
1102
1103 static int
1104 compare_real (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1105 {
1106   int rc;
1107   switch (op)
1108     {
1109       case INTRINSIC_EQ:
1110         rc = mpfr_equal_p (op1->value.real, op2->value.real) ? 0 : 1;
1111         break;
1112       case INTRINSIC_GT:
1113         rc = mpfr_greater_p (op1->value.real, op2->value.real) ? 1 : -1;
1114         break;
1115       case INTRINSIC_GE:
1116         rc = mpfr_greaterequal_p (op1->value.real, op2->value.real) ? 1 : -1;
1117         break;
1118       case INTRINSIC_LT:
1119         rc = mpfr_less_p (op1->value.real, op2->value.real) ? -1 : 1;
1120         break;
1121       case INTRINSIC_LE:
1122         rc = mpfr_lessequal_p (op1->value.real, op2->value.real) ? -1 : 1;
1123         break;
1124       default:
1125         gfc_internal_error ("compare_real(): Bad operator");
1126     }
1127
1128   return rc;
1129 }
1130
1131 /* Comparison operators.  Assumes that the two expression nodes
1132    contain two constants of the same type. The op argument is
1133    needed to handle NaN correctly.  */
1134
1135 int
1136 gfc_compare_expr (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1137 {
1138   int rc;
1139
1140   switch (op1->ts.type)
1141     {
1142     case BT_INTEGER:
1143       rc = mpz_cmp (op1->value.integer, op2->value.integer);
1144       break;
1145
1146     case BT_REAL:
1147       rc = compare_real (op1, op2, op);
1148       break;
1149
1150     case BT_CHARACTER:
1151       rc = gfc_compare_string (op1, op2);
1152       break;
1153
1154     case BT_LOGICAL:
1155       rc = ((!op1->value.logical && op2->value.logical)
1156             || (op1->value.logical && !op2->value.logical));
1157       break;
1158
1159     default:
1160       gfc_internal_error ("gfc_compare_expr(): Bad basic type");
1161     }
1162
1163   return rc;
1164 }
1165
1166
1167 /* Compare a pair of complex numbers.  Naturally, this is only for
1168    equality and inequality.  */
1169
1170 static int
1171 compare_complex (gfc_expr *op1, gfc_expr *op2)
1172 {
1173   return (mpfr_equal_p (op1->value.complex.r, op2->value.complex.r)
1174           && mpfr_equal_p (op1->value.complex.i, op2->value.complex.i));
1175 }
1176
1177
1178 /* Given two constant strings and the inverse collating sequence, compare the
1179    strings.  We return -1 for a < b, 0 for a == b and 1 for a > b. 
1180    We use the processor's default collating sequence.  */
1181
1182 int
1183 gfc_compare_string (gfc_expr *a, gfc_expr *b)
1184 {
1185   int len, alen, blen, i;
1186   gfc_char_t ac, bc;
1187
1188   alen = a->value.character.length;
1189   blen = b->value.character.length;
1190
1191   len = MAX(alen, blen);
1192
1193   for (i = 0; i < len; i++)
1194     {
1195       ac = ((i < alen) ? a->value.character.string[i] : ' ');
1196       bc = ((i < blen) ? b->value.character.string[i] : ' ');
1197
1198       if (ac < bc)
1199         return -1;
1200       if (ac > bc)
1201         return 1;
1202     }
1203
1204   /* Strings are equal */
1205   return 0;
1206 }
1207
1208
1209 int
1210 gfc_compare_with_Cstring (gfc_expr *a, const char *b, bool case_sensitive)
1211 {
1212   int len, alen, blen, i;
1213   gfc_char_t ac, bc;
1214
1215   alen = a->value.character.length;
1216   blen = strlen (b);
1217
1218   len = MAX(alen, blen);
1219
1220   for (i = 0; i < len; i++)
1221     {
1222       ac = ((i < alen) ? a->value.character.string[i] : ' ');
1223       bc = ((i < blen) ? b[i] : ' ');
1224
1225       if (!case_sensitive)
1226         {
1227           ac = TOLOWER (ac);
1228           bc = TOLOWER (bc);
1229         }
1230
1231       if (ac < bc)
1232         return -1;
1233       if (ac > bc)
1234         return 1;
1235     }
1236
1237   /* Strings are equal */
1238   return 0;
1239 }
1240
1241
1242 /* Specific comparison subroutines.  */
1243
1244 static arith
1245 gfc_arith_eq (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1246 {
1247   gfc_expr *result;
1248
1249   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1250                                 &op1->where);
1251   result->value.logical = (op1->ts.type == BT_COMPLEX)
1252                         ? compare_complex (op1, op2)
1253                         : (gfc_compare_expr (op1, op2, INTRINSIC_EQ) == 0);
1254
1255   *resultp = result;
1256   return ARITH_OK;
1257 }
1258
1259
1260 static arith
1261 gfc_arith_ne (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1262 {
1263   gfc_expr *result;
1264
1265   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1266                                 &op1->where);
1267   result->value.logical = (op1->ts.type == BT_COMPLEX)
1268                         ? !compare_complex (op1, op2)
1269                         : (gfc_compare_expr (op1, op2, INTRINSIC_EQ) != 0);
1270
1271   *resultp = result;
1272   return ARITH_OK;
1273 }
1274
1275
1276 static arith
1277 gfc_arith_gt (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1278 {
1279   gfc_expr *result;
1280
1281   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1282                                 &op1->where);
1283   result->value.logical = (gfc_compare_expr (op1, op2, INTRINSIC_GT) > 0);
1284   *resultp = result;
1285
1286   return ARITH_OK;
1287 }
1288
1289
1290 static arith
1291 gfc_arith_ge (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1292 {
1293   gfc_expr *result;
1294
1295   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1296                                 &op1->where);
1297   result->value.logical = (gfc_compare_expr (op1, op2, INTRINSIC_GE) >= 0);
1298   *resultp = result;
1299
1300   return ARITH_OK;
1301 }
1302
1303
1304 static arith
1305 gfc_arith_lt (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1306 {
1307   gfc_expr *result;
1308
1309   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1310                                 &op1->where);
1311   result->value.logical = (gfc_compare_expr (op1, op2, INTRINSIC_LT) < 0);
1312   *resultp = result;
1313
1314   return ARITH_OK;
1315 }
1316
1317
1318 static arith
1319 gfc_arith_le (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1320 {
1321   gfc_expr *result;
1322
1323   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1324                                 &op1->where);
1325   result->value.logical = (gfc_compare_expr (op1, op2, INTRINSIC_LE) <= 0);
1326   *resultp = result;
1327
1328   return ARITH_OK;
1329 }
1330
1331
1332 static arith
1333 reduce_unary (arith (*eval) (gfc_expr *, gfc_expr **), gfc_expr *op,
1334               gfc_expr **result)
1335 {
1336   gfc_constructor *c, *head;
1337   gfc_expr *r;
1338   arith rc;
1339
1340   if (op->expr_type == EXPR_CONSTANT)
1341     return eval (op, result);
1342
1343   rc = ARITH_OK;
1344   head = gfc_copy_constructor (op->value.constructor);
1345
1346   for (c = head; c; c = c->next)
1347     {
1348       rc = reduce_unary (eval, c->expr, &r);
1349
1350       if (rc != ARITH_OK)
1351         break;
1352
1353       gfc_replace_expr (c->expr, r);
1354     }
1355
1356   if (rc != ARITH_OK)
1357     gfc_free_constructor (head);
1358   else
1359     {
1360       r = gfc_get_expr ();
1361       r->expr_type = EXPR_ARRAY;
1362       r->value.constructor = head;
1363       r->shape = gfc_copy_shape (op->shape, op->rank);
1364
1365       r->ts = head->expr->ts;
1366       r->where = op->where;
1367       r->rank = op->rank;
1368
1369       *result = r;
1370     }
1371
1372   return rc;
1373 }
1374
1375
1376 static arith
1377 reduce_binary_ac (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1378                   gfc_expr *op1, gfc_expr *op2, gfc_expr **result)
1379 {
1380   gfc_constructor *c, *head;
1381   gfc_expr *r;
1382   arith rc;
1383
1384   head = gfc_copy_constructor (op1->value.constructor);
1385   rc = ARITH_OK;
1386
1387   for (c = head; c; c = c->next)
1388     {
1389       if (c->expr->expr_type == EXPR_CONSTANT)
1390         rc = eval (c->expr, op2, &r);
1391       else
1392         rc = reduce_binary_ac (eval, c->expr, op2, &r);
1393
1394       if (rc != ARITH_OK)
1395         break;
1396
1397       gfc_replace_expr (c->expr, r);
1398     }
1399
1400   if (rc != ARITH_OK)
1401     gfc_free_constructor (head);
1402   else
1403     {
1404       r = gfc_get_expr ();
1405       r->expr_type = EXPR_ARRAY;
1406       r->value.constructor = head;
1407       r->shape = gfc_copy_shape (op1->shape, op1->rank);
1408
1409       r->ts = head->expr->ts;
1410       r->where = op1->where;
1411       r->rank = op1->rank;
1412
1413       *result = r;
1414     }
1415
1416   return rc;
1417 }
1418
1419
1420 static arith
1421 reduce_binary_ca (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1422                   gfc_expr *op1, gfc_expr *op2, gfc_expr **result)
1423 {
1424   gfc_constructor *c, *head;
1425   gfc_expr *r;
1426   arith rc;
1427
1428   head = gfc_copy_constructor (op2->value.constructor);
1429   rc = ARITH_OK;
1430
1431   for (c = head; c; c = c->next)
1432     {
1433       if (c->expr->expr_type == EXPR_CONSTANT)
1434         rc = eval (op1, c->expr, &r);
1435       else
1436         rc = reduce_binary_ca (eval, op1, c->expr, &r);
1437
1438       if (rc != ARITH_OK)
1439         break;
1440
1441       gfc_replace_expr (c->expr, r);
1442     }
1443
1444   if (rc != ARITH_OK)
1445     gfc_free_constructor (head);
1446   else
1447     {
1448       r = gfc_get_expr ();
1449       r->expr_type = EXPR_ARRAY;
1450       r->value.constructor = head;
1451       r->shape = gfc_copy_shape (op2->shape, op2->rank);
1452
1453       r->ts = head->expr->ts;
1454       r->where = op2->where;
1455       r->rank = op2->rank;
1456
1457       *result = r;
1458     }
1459
1460   return rc;
1461 }
1462
1463
1464 /* We need a forward declaration of reduce_binary.  */
1465 static arith reduce_binary (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1466                             gfc_expr *op1, gfc_expr *op2, gfc_expr **result);
1467
1468
1469 static arith
1470 reduce_binary_aa (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1471                   gfc_expr *op1, gfc_expr *op2, gfc_expr **result)
1472 {
1473   gfc_constructor *c, *d, *head;
1474   gfc_expr *r;
1475   arith rc;
1476
1477   head = gfc_copy_constructor (op1->value.constructor);
1478
1479   rc = ARITH_OK;
1480   d = op2->value.constructor;
1481
1482   if (gfc_check_conformance ("elemental binary operation", op1, op2)
1483       != SUCCESS)
1484     rc = ARITH_INCOMMENSURATE;
1485   else
1486     {
1487       for (c = head; c; c = c->next, d = d->next)
1488         {
1489           if (d == NULL)
1490             {
1491               rc = ARITH_INCOMMENSURATE;
1492               break;
1493             }
1494
1495           rc = reduce_binary (eval, c->expr, d->expr, &r);
1496           if (rc != ARITH_OK)
1497             break;
1498
1499           gfc_replace_expr (c->expr, r);
1500         }
1501
1502       if (d != NULL)
1503         rc = ARITH_INCOMMENSURATE;
1504     }
1505
1506   if (rc != ARITH_OK)
1507     gfc_free_constructor (head);
1508   else
1509     {
1510       r = gfc_get_expr ();
1511       r->expr_type = EXPR_ARRAY;
1512       r->value.constructor = head;
1513       r->shape = gfc_copy_shape (op1->shape, op1->rank);
1514
1515       r->ts = head->expr->ts;
1516       r->where = op1->where;
1517       r->rank = op1->rank;
1518
1519       *result = r;
1520     }
1521
1522   return rc;
1523 }
1524
1525
1526 static arith
1527 reduce_binary (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1528                gfc_expr *op1, gfc_expr *op2, gfc_expr **result)
1529 {
1530   if (op1->expr_type == EXPR_CONSTANT && op2->expr_type == EXPR_CONSTANT)
1531     return eval (op1, op2, result);
1532
1533   if (op1->expr_type == EXPR_CONSTANT && op2->expr_type == EXPR_ARRAY)
1534     return reduce_binary_ca (eval, op1, op2, result);
1535
1536   if (op1->expr_type == EXPR_ARRAY && op2->expr_type == EXPR_CONSTANT)
1537     return reduce_binary_ac (eval, op1, op2, result);
1538
1539   return reduce_binary_aa (eval, op1, op2, result);
1540 }
1541
1542
1543 typedef union
1544 {
1545   arith (*f2)(gfc_expr *, gfc_expr **);
1546   arith (*f3)(gfc_expr *, gfc_expr *, gfc_expr **);
1547 }
1548 eval_f;
1549
1550 /* High level arithmetic subroutines.  These subroutines go into
1551    eval_intrinsic(), which can do one of several things to its
1552    operands.  If the operands are incompatible with the intrinsic
1553    operation, we return a node pointing to the operands and hope that
1554    an operator interface is found during resolution.
1555
1556    If the operands are compatible and are constants, then we try doing
1557    the arithmetic.  We also handle the cases where either or both
1558    operands are array constructors.  */
1559
1560 static gfc_expr *
1561 eval_intrinsic (gfc_intrinsic_op op,
1562                 eval_f eval, gfc_expr *op1, gfc_expr *op2)
1563 {
1564   gfc_expr temp, *result;
1565   int unary;
1566   arith rc;
1567
1568   gfc_clear_ts (&temp.ts);
1569
1570   switch (op)
1571     {
1572     /* Logical unary  */
1573     case INTRINSIC_NOT:
1574       if (op1->ts.type != BT_LOGICAL)
1575         goto runtime;
1576
1577       temp.ts.type = BT_LOGICAL;
1578       temp.ts.kind = gfc_default_logical_kind;
1579       unary = 1;
1580       break;
1581
1582     /* Logical binary operators  */
1583     case INTRINSIC_OR:
1584     case INTRINSIC_AND:
1585     case INTRINSIC_NEQV:
1586     case INTRINSIC_EQV:
1587       if (op1->ts.type != BT_LOGICAL || op2->ts.type != BT_LOGICAL)
1588         goto runtime;
1589
1590       temp.ts.type = BT_LOGICAL;
1591       temp.ts.kind = gfc_default_logical_kind;
1592       unary = 0;
1593       break;
1594
1595     /* Numeric unary  */
1596     case INTRINSIC_UPLUS:
1597     case INTRINSIC_UMINUS:
1598       if (!gfc_numeric_ts (&op1->ts))
1599         goto runtime;
1600
1601       temp.ts = op1->ts;
1602       unary = 1;
1603       break;
1604
1605     case INTRINSIC_PARENTHESES:
1606       temp.ts = op1->ts;
1607       unary = 1;
1608       break;
1609
1610     /* Additional restrictions for ordering relations.  */
1611     case INTRINSIC_GE:
1612     case INTRINSIC_GE_OS:
1613     case INTRINSIC_LT:
1614     case INTRINSIC_LT_OS:
1615     case INTRINSIC_LE:
1616     case INTRINSIC_LE_OS:
1617     case INTRINSIC_GT:
1618     case INTRINSIC_GT_OS:
1619       if (op1->ts.type == BT_COMPLEX || op2->ts.type == BT_COMPLEX)
1620         {
1621           temp.ts.type = BT_LOGICAL;
1622           temp.ts.kind = gfc_default_logical_kind;
1623           goto runtime;
1624         }
1625
1626     /* Fall through  */
1627     case INTRINSIC_EQ:
1628     case INTRINSIC_EQ_OS:
1629     case INTRINSIC_NE:
1630     case INTRINSIC_NE_OS:
1631       if (op1->ts.type == BT_CHARACTER && op2->ts.type == BT_CHARACTER)
1632         {
1633           unary = 0;
1634           temp.ts.type = BT_LOGICAL;
1635           temp.ts.kind = gfc_default_logical_kind;
1636
1637           /* If kind mismatch, exit and we'll error out later.  */
1638           if (op1->ts.kind != op2->ts.kind)
1639             goto runtime;
1640
1641           break;
1642         }
1643
1644     /* Fall through  */
1645     /* Numeric binary  */
1646     case INTRINSIC_PLUS:
1647     case INTRINSIC_MINUS:
1648     case INTRINSIC_TIMES:
1649     case INTRINSIC_DIVIDE:
1650     case INTRINSIC_POWER:
1651       if (!gfc_numeric_ts (&op1->ts) || !gfc_numeric_ts (&op2->ts))
1652         goto runtime;
1653
1654       /* Insert any necessary type conversions to make the operands
1655          compatible.  */
1656
1657       temp.expr_type = EXPR_OP;
1658       gfc_clear_ts (&temp.ts);
1659       temp.value.op.op = op;
1660
1661       temp.value.op.op1 = op1;
1662       temp.value.op.op2 = op2;
1663
1664       gfc_type_convert_binary (&temp);
1665
1666       if (op == INTRINSIC_EQ || op == INTRINSIC_NE
1667           || op == INTRINSIC_GE || op == INTRINSIC_GT
1668           || op == INTRINSIC_LE || op == INTRINSIC_LT
1669           || op == INTRINSIC_EQ_OS || op == INTRINSIC_NE_OS
1670           || op == INTRINSIC_GE_OS || op == INTRINSIC_GT_OS
1671           || op == INTRINSIC_LE_OS || op == INTRINSIC_LT_OS)
1672         {
1673           temp.ts.type = BT_LOGICAL;
1674           temp.ts.kind = gfc_default_logical_kind;
1675         }
1676
1677       unary = 0;
1678       break;
1679
1680     /* Character binary  */
1681     case INTRINSIC_CONCAT:
1682       if (op1->ts.type != BT_CHARACTER || op2->ts.type != BT_CHARACTER
1683           || op1->ts.kind != op2->ts.kind)
1684         goto runtime;
1685
1686       temp.ts.type = BT_CHARACTER;
1687       temp.ts.kind = op1->ts.kind;
1688       unary = 0;
1689       break;
1690
1691     case INTRINSIC_USER:
1692       goto runtime;
1693
1694     default:
1695       gfc_internal_error ("eval_intrinsic(): Bad operator");
1696     }
1697
1698   /* Try to combine the operators.  */
1699   if (op == INTRINSIC_POWER && op2->ts.type != BT_INTEGER)
1700     goto runtime;
1701
1702   if (op1->expr_type != EXPR_CONSTANT
1703       && (op1->expr_type != EXPR_ARRAY
1704           || !gfc_is_constant_expr (op1) || !gfc_expanded_ac (op1)))
1705     goto runtime;
1706
1707   if (op2 != NULL
1708       && op2->expr_type != EXPR_CONSTANT
1709          && (op2->expr_type != EXPR_ARRAY
1710              || !gfc_is_constant_expr (op2) || !gfc_expanded_ac (op2)))
1711     goto runtime;
1712
1713   if (unary)
1714     rc = reduce_unary (eval.f2, op1, &result);
1715   else
1716     rc = reduce_binary (eval.f3, op1, op2, &result);
1717
1718   if (rc != ARITH_OK)
1719     { /* Something went wrong.  */
1720       gfc_error (gfc_arith_error (rc), &op1->where);
1721       return NULL;
1722     }
1723
1724   gfc_free_expr (op1);
1725   gfc_free_expr (op2);
1726   return result;
1727
1728 runtime:
1729   /* Create a run-time expression.  */
1730   result = gfc_get_expr ();
1731   result->ts = temp.ts;
1732
1733   result->expr_type = EXPR_OP;
1734   result->value.op.op = op;
1735
1736   result->value.op.op1 = op1;
1737   result->value.op.op2 = op2;
1738
1739   result->where = op1->where;
1740
1741   return result;
1742 }
1743
1744
1745 /* Modify type of expression for zero size array.  */
1746
1747 static gfc_expr *
1748 eval_type_intrinsic0 (gfc_intrinsic_op iop, gfc_expr *op)
1749 {
1750   if (op == NULL)
1751     gfc_internal_error ("eval_type_intrinsic0(): op NULL");
1752
1753   switch (iop)
1754     {
1755     case INTRINSIC_GE:
1756     case INTRINSIC_GE_OS:
1757     case INTRINSIC_LT:
1758     case INTRINSIC_LT_OS:
1759     case INTRINSIC_LE:
1760     case INTRINSIC_LE_OS:
1761     case INTRINSIC_GT:
1762     case INTRINSIC_GT_OS:
1763     case INTRINSIC_EQ:
1764     case INTRINSIC_EQ_OS:
1765     case INTRINSIC_NE:
1766     case INTRINSIC_NE_OS:
1767       op->ts.type = BT_LOGICAL;
1768       op->ts.kind = gfc_default_logical_kind;
1769       break;
1770
1771     default:
1772       break;
1773     }
1774
1775   return op;
1776 }
1777
1778
1779 /* Return nonzero if the expression is a zero size array.  */
1780
1781 static int
1782 gfc_zero_size_array (gfc_expr *e)
1783 {
1784   if (e->expr_type != EXPR_ARRAY)
1785     return 0;
1786
1787   return e->value.constructor == NULL;
1788 }
1789
1790
1791 /* Reduce a binary expression where at least one of the operands
1792    involves a zero-length array.  Returns NULL if neither of the
1793    operands is a zero-length array.  */
1794
1795 static gfc_expr *
1796 reduce_binary0 (gfc_expr *op1, gfc_expr *op2)
1797 {
1798   if (gfc_zero_size_array (op1))
1799     {
1800       gfc_free_expr (op2);
1801       return op1;
1802     }
1803
1804   if (gfc_zero_size_array (op2))
1805     {
1806       gfc_free_expr (op1);
1807       return op2;
1808     }
1809
1810   return NULL;
1811 }
1812
1813
1814 static gfc_expr *
1815 eval_intrinsic_f2 (gfc_intrinsic_op op,
1816                    arith (*eval) (gfc_expr *, gfc_expr **),
1817                    gfc_expr *op1, gfc_expr *op2)
1818 {
1819   gfc_expr *result;
1820   eval_f f;
1821
1822   if (op2 == NULL)
1823     {
1824       if (gfc_zero_size_array (op1))
1825         return eval_type_intrinsic0 (op, op1);
1826     }
1827   else
1828     {
1829       result = reduce_binary0 (op1, op2);
1830       if (result != NULL)
1831         return eval_type_intrinsic0 (op, result);
1832     }
1833
1834   f.f2 = eval;
1835   return eval_intrinsic (op, f, op1, op2);
1836 }
1837
1838
1839 static gfc_expr *
1840 eval_intrinsic_f3 (gfc_intrinsic_op op,
1841                    arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1842                    gfc_expr *op1, gfc_expr *op2)
1843 {
1844   gfc_expr *result;
1845   eval_f f;
1846
1847   result = reduce_binary0 (op1, op2);
1848   if (result != NULL)
1849     return eval_type_intrinsic0(op, result);
1850
1851   f.f3 = eval;
1852   return eval_intrinsic (op, f, op1, op2);
1853 }
1854
1855
1856 gfc_expr *
1857 gfc_parentheses (gfc_expr *op)
1858 {
1859   if (gfc_is_constant_expr (op))
1860     return op;
1861
1862   return eval_intrinsic_f2 (INTRINSIC_PARENTHESES, gfc_arith_identity,
1863                             op, NULL);
1864 }
1865
1866 gfc_expr *
1867 gfc_uplus (gfc_expr *op)
1868 {
1869   return eval_intrinsic_f2 (INTRINSIC_UPLUS, gfc_arith_identity, op, NULL);
1870 }
1871
1872
1873 gfc_expr *
1874 gfc_uminus (gfc_expr *op)
1875 {
1876   return eval_intrinsic_f2 (INTRINSIC_UMINUS, gfc_arith_uminus, op, NULL);
1877 }
1878
1879
1880 gfc_expr *
1881 gfc_add (gfc_expr *op1, gfc_expr *op2)
1882 {
1883   return eval_intrinsic_f3 (INTRINSIC_PLUS, gfc_arith_plus, op1, op2);
1884 }
1885
1886
1887 gfc_expr *
1888 gfc_subtract (gfc_expr *op1, gfc_expr *op2)
1889 {
1890   return eval_intrinsic_f3 (INTRINSIC_MINUS, gfc_arith_minus, op1, op2);
1891 }
1892
1893
1894 gfc_expr *
1895 gfc_multiply (gfc_expr *op1, gfc_expr *op2)
1896 {
1897   return eval_intrinsic_f3 (INTRINSIC_TIMES, gfc_arith_times, op1, op2);
1898 }
1899
1900
1901 gfc_expr *
1902 gfc_divide (gfc_expr *op1, gfc_expr *op2)
1903 {
1904   return eval_intrinsic_f3 (INTRINSIC_DIVIDE, gfc_arith_divide, op1, op2);
1905 }
1906
1907
1908 gfc_expr *
1909 gfc_power (gfc_expr *op1, gfc_expr *op2)
1910 {
1911   return eval_intrinsic_f3 (INTRINSIC_POWER, gfc_arith_power, op1, op2);
1912 }
1913
1914
1915 gfc_expr *
1916 gfc_concat (gfc_expr *op1, gfc_expr *op2)
1917 {
1918   return eval_intrinsic_f3 (INTRINSIC_CONCAT, gfc_arith_concat, op1, op2);
1919 }
1920
1921
1922 gfc_expr *
1923 gfc_and (gfc_expr *op1, gfc_expr *op2)
1924 {
1925   return eval_intrinsic_f3 (INTRINSIC_AND, gfc_arith_and, op1, op2);
1926 }
1927
1928
1929 gfc_expr *
1930 gfc_or (gfc_expr *op1, gfc_expr *op2)
1931 {
1932   return eval_intrinsic_f3 (INTRINSIC_OR, gfc_arith_or, op1, op2);
1933 }
1934
1935
1936 gfc_expr *
1937 gfc_not (gfc_expr *op1)
1938 {
1939   return eval_intrinsic_f2 (INTRINSIC_NOT, gfc_arith_not, op1, NULL);
1940 }
1941
1942
1943 gfc_expr *
1944 gfc_eqv (gfc_expr *op1, gfc_expr *op2)
1945 {
1946   return eval_intrinsic_f3 (INTRINSIC_EQV, gfc_arith_eqv, op1, op2);
1947 }
1948
1949
1950 gfc_expr *
1951 gfc_neqv (gfc_expr *op1, gfc_expr *op2)
1952 {
1953   return eval_intrinsic_f3 (INTRINSIC_NEQV, gfc_arith_neqv, op1, op2);
1954 }
1955
1956
1957 gfc_expr *
1958 gfc_eq (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1959 {
1960   return eval_intrinsic_f3 (op, gfc_arith_eq, op1, op2);
1961 }
1962
1963
1964 gfc_expr *
1965 gfc_ne (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1966 {
1967   return eval_intrinsic_f3 (op, gfc_arith_ne, op1, op2);
1968 }
1969
1970
1971 gfc_expr *
1972 gfc_gt (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1973 {
1974   return eval_intrinsic_f3 (op, gfc_arith_gt, op1, op2);
1975 }
1976
1977
1978 gfc_expr *
1979 gfc_ge (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1980 {
1981   return eval_intrinsic_f3 (op, gfc_arith_ge, op1, op2);
1982 }
1983
1984
1985 gfc_expr *
1986 gfc_lt (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1987 {
1988   return eval_intrinsic_f3 (op, gfc_arith_lt, op1, op2);
1989 }
1990
1991
1992 gfc_expr *
1993 gfc_le (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1994 {
1995   return eval_intrinsic_f3 (op, gfc_arith_le, op1, op2);
1996 }
1997
1998
1999 /* Convert an integer string to an expression node.  */
2000
2001 gfc_expr *
2002 gfc_convert_integer (const char *buffer, int kind, int radix, locus *where)
2003 {
2004   gfc_expr *e;
2005   const char *t;
2006
2007   e = gfc_constant_result (BT_INTEGER, kind, where);
2008   /* A leading plus is allowed, but not by mpz_set_str.  */
2009   if (buffer[0] == '+')
2010     t = buffer + 1;
2011   else
2012     t = buffer;
2013   mpz_set_str (e->value.integer, t, radix);
2014
2015   return e;
2016 }
2017
2018
2019 /* Convert a real string to an expression node.  */
2020
2021 gfc_expr *
2022 gfc_convert_real (const char *buffer, int kind, locus *where)
2023 {
2024   gfc_expr *e;
2025
2026   e = gfc_constant_result (BT_REAL, kind, where);
2027   mpfr_set_str (e->value.real, buffer, 10, GFC_RND_MODE);
2028
2029   return e;
2030 }
2031
2032
2033 /* Convert a pair of real, constant expression nodes to a single
2034    complex expression node.  */
2035
2036 gfc_expr *
2037 gfc_convert_complex (gfc_expr *real, gfc_expr *imag, int kind)
2038 {
2039   gfc_expr *e;
2040
2041   e = gfc_constant_result (BT_COMPLEX, kind, &real->where);
2042   mpfr_set (e->value.complex.r, real->value.real, GFC_RND_MODE);
2043   mpfr_set (e->value.complex.i, imag->value.real, GFC_RND_MODE);
2044
2045   return e;
2046 }
2047
2048
2049 /******* Simplification of intrinsic functions with constant arguments *****/
2050
2051
2052 /* Deal with an arithmetic error.  */
2053
2054 static void
2055 arith_error (arith rc, gfc_typespec *from, gfc_typespec *to, locus *where)
2056 {
2057   switch (rc)
2058     {
2059     case ARITH_OK:
2060       gfc_error ("Arithmetic OK converting %s to %s at %L",
2061                  gfc_typename (from), gfc_typename (to), where);
2062       break;
2063     case ARITH_OVERFLOW:
2064       gfc_error ("Arithmetic overflow converting %s to %s at %L. This check "
2065                  "can be disabled with the option -fno-range-check",
2066                  gfc_typename (from), gfc_typename (to), where);
2067       break;
2068     case ARITH_UNDERFLOW:
2069       gfc_error ("Arithmetic underflow converting %s to %s at %L. This check "
2070                  "can be disabled with the option -fno-range-check",
2071                  gfc_typename (from), gfc_typename (to), where);
2072       break;
2073     case ARITH_NAN:
2074       gfc_error ("Arithmetic NaN converting %s to %s at %L. This check "
2075                  "can be disabled with the option -fno-range-check",
2076                  gfc_typename (from), gfc_typename (to), where);
2077       break;
2078     case ARITH_DIV0:
2079       gfc_error ("Division by zero converting %s to %s at %L",
2080                  gfc_typename (from), gfc_typename (to), where);
2081       break;
2082     case ARITH_INCOMMENSURATE:
2083       gfc_error ("Array operands are incommensurate converting %s to %s at %L",
2084                  gfc_typename (from), gfc_typename (to), where);
2085       break;
2086     case ARITH_ASYMMETRIC:
2087       gfc_error ("Integer outside symmetric range implied by Standard Fortran"
2088                  " converting %s to %s at %L",
2089                  gfc_typename (from), gfc_typename (to), where);
2090       break;
2091     default:
2092       gfc_internal_error ("gfc_arith_error(): Bad error code");
2093     }
2094
2095   /* TODO: Do something about the error, i.e., throw exception, return
2096      NaN, etc.  */
2097 }
2098
2099
2100 /* Convert integers to integers.  */
2101
2102 gfc_expr *
2103 gfc_int2int (gfc_expr *src, int kind)
2104 {
2105   gfc_expr *result;
2106   arith rc;
2107
2108   result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2109
2110   mpz_set (result->value.integer, src->value.integer);
2111
2112   if ((rc = gfc_check_integer_range (result->value.integer, kind)) != ARITH_OK)
2113     {
2114       if (rc == ARITH_ASYMMETRIC)
2115         {
2116           gfc_warning (gfc_arith_error (rc), &src->where);
2117         }
2118       else
2119         {
2120           arith_error (rc, &src->ts, &result->ts, &src->where);
2121           gfc_free_expr (result);
2122           return NULL;
2123         }
2124     }
2125
2126   return result;
2127 }
2128
2129
2130 /* Convert integers to reals.  */
2131
2132 gfc_expr *
2133 gfc_int2real (gfc_expr *src, int kind)
2134 {
2135   gfc_expr *result;
2136   arith rc;
2137
2138   result = gfc_constant_result (BT_REAL, kind, &src->where);
2139
2140   mpfr_set_z (result->value.real, src->value.integer, GFC_RND_MODE);
2141
2142   if ((rc = gfc_check_real_range (result->value.real, kind)) != ARITH_OK)
2143     {
2144       arith_error (rc, &src->ts, &result->ts, &src->where);
2145       gfc_free_expr (result);
2146       return NULL;
2147     }
2148
2149   return result;
2150 }
2151
2152
2153 /* Convert default integer to default complex.  */
2154
2155 gfc_expr *
2156 gfc_int2complex (gfc_expr *src, int kind)
2157 {
2158   gfc_expr *result;
2159   arith rc;
2160
2161   result = gfc_constant_result (BT_COMPLEX, kind, &src->where);
2162
2163   mpfr_set_z (result->value.complex.r, src->value.integer, GFC_RND_MODE);
2164   mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
2165
2166   if ((rc = gfc_check_real_range (result->value.complex.r, kind)) != ARITH_OK)
2167     {
2168       arith_error (rc, &src->ts, &result->ts, &src->where);
2169       gfc_free_expr (result);
2170       return NULL;
2171     }
2172
2173   return result;
2174 }
2175
2176
2177 /* Convert default real to default integer.  */
2178
2179 gfc_expr *
2180 gfc_real2int (gfc_expr *src, int kind)
2181 {
2182   gfc_expr *result;
2183   arith rc;
2184
2185   result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2186
2187   gfc_mpfr_to_mpz (result->value.integer, src->value.real, &src->where);
2188
2189   if ((rc = gfc_check_integer_range (result->value.integer, kind)) != ARITH_OK)
2190     {
2191       arith_error (rc, &src->ts, &result->ts, &src->where);
2192       gfc_free_expr (result);
2193       return NULL;
2194     }
2195
2196   return result;
2197 }
2198
2199
2200 /* Convert real to real.  */
2201
2202 gfc_expr *
2203 gfc_real2real (gfc_expr *src, int kind)
2204 {
2205   gfc_expr *result;
2206   arith rc;
2207
2208   result = gfc_constant_result (BT_REAL, kind, &src->where);
2209
2210   mpfr_set (result->value.real, src->value.real, GFC_RND_MODE);
2211
2212   rc = gfc_check_real_range (result->value.real, kind);
2213
2214   if (rc == ARITH_UNDERFLOW)
2215     {
2216       if (gfc_option.warn_underflow)
2217         gfc_warning (gfc_arith_error (rc), &src->where);
2218       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
2219     }
2220   else if (rc != ARITH_OK)
2221     {
2222       arith_error (rc, &src->ts, &result->ts, &src->where);
2223       gfc_free_expr (result);
2224       return NULL;
2225     }
2226
2227   return result;
2228 }
2229
2230
2231 /* Convert real to complex.  */
2232
2233 gfc_expr *
2234 gfc_real2complex (gfc_expr *src, int kind)
2235 {
2236   gfc_expr *result;
2237   arith rc;
2238
2239   result = gfc_constant_result (BT_COMPLEX, kind, &src->where);
2240
2241   mpfr_set (result->value.complex.r, src->value.real, GFC_RND_MODE);
2242   mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
2243
2244   rc = gfc_check_real_range (result->value.complex.r, kind);
2245
2246   if (rc == ARITH_UNDERFLOW)
2247     {
2248       if (gfc_option.warn_underflow)
2249         gfc_warning (gfc_arith_error (rc), &src->where);
2250       mpfr_set_ui (result->value.complex.r, 0, GFC_RND_MODE);
2251     }
2252   else if (rc != ARITH_OK)
2253     {
2254       arith_error (rc, &src->ts, &result->ts, &src->where);
2255       gfc_free_expr (result);
2256       return NULL;
2257     }
2258
2259   return result;
2260 }
2261
2262
2263 /* Convert complex to integer.  */
2264
2265 gfc_expr *
2266 gfc_complex2int (gfc_expr *src, int kind)
2267 {
2268   gfc_expr *result;
2269   arith rc;
2270
2271   result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2272
2273   gfc_mpfr_to_mpz (result->value.integer, src->value.complex.r, &src->where);
2274
2275   if ((rc = gfc_check_integer_range (result->value.integer, kind)) != ARITH_OK)
2276     {
2277       arith_error (rc, &src->ts, &result->ts, &src->where);
2278       gfc_free_expr (result);
2279       return NULL;
2280     }
2281
2282   return result;
2283 }
2284
2285
2286 /* Convert complex to real.  */
2287
2288 gfc_expr *
2289 gfc_complex2real (gfc_expr *src, int kind)
2290 {
2291   gfc_expr *result;
2292   arith rc;
2293
2294   result = gfc_constant_result (BT_REAL, kind, &src->where);
2295
2296   mpfr_set (result->value.real, src->value.complex.r, GFC_RND_MODE);
2297
2298   rc = gfc_check_real_range (result->value.real, kind);
2299
2300   if (rc == ARITH_UNDERFLOW)
2301     {
2302       if (gfc_option.warn_underflow)
2303         gfc_warning (gfc_arith_error (rc), &src->where);
2304       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
2305     }
2306   if (rc != ARITH_OK)
2307     {
2308       arith_error (rc, &src->ts, &result->ts, &src->where);
2309       gfc_free_expr (result);
2310       return NULL;
2311     }
2312
2313   return result;
2314 }
2315
2316
2317 /* Convert complex to complex.  */
2318
2319 gfc_expr *
2320 gfc_complex2complex (gfc_expr *src, int kind)
2321 {
2322   gfc_expr *result;
2323   arith rc;
2324
2325   result = gfc_constant_result (BT_COMPLEX, kind, &src->where);
2326
2327   mpfr_set (result->value.complex.r, src->value.complex.r, GFC_RND_MODE);
2328   mpfr_set (result->value.complex.i, src->value.complex.i, GFC_RND_MODE);
2329
2330   rc = gfc_check_real_range (result->value.complex.r, kind);
2331
2332   if (rc == ARITH_UNDERFLOW)
2333     {
2334       if (gfc_option.warn_underflow)
2335         gfc_warning (gfc_arith_error (rc), &src->where);
2336       mpfr_set_ui (result->value.complex.r, 0, GFC_RND_MODE);
2337     }
2338   else if (rc != ARITH_OK)
2339     {
2340       arith_error (rc, &src->ts, &result->ts, &src->where);
2341       gfc_free_expr (result);
2342       return NULL;
2343     }
2344
2345   rc = gfc_check_real_range (result->value.complex.i, kind);
2346
2347   if (rc == ARITH_UNDERFLOW)
2348     {
2349       if (gfc_option.warn_underflow)
2350         gfc_warning (gfc_arith_error (rc), &src->where);
2351       mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
2352     }
2353   else if (rc != ARITH_OK)
2354     {
2355       arith_error (rc, &src->ts, &result->ts, &src->where);
2356       gfc_free_expr (result);
2357       return NULL;
2358     }
2359
2360   return result;
2361 }
2362
2363
2364 /* Logical kind conversion.  */
2365
2366 gfc_expr *
2367 gfc_log2log (gfc_expr *src, int kind)
2368 {
2369   gfc_expr *result;
2370
2371   result = gfc_constant_result (BT_LOGICAL, kind, &src->where);
2372   result->value.logical = src->value.logical;
2373
2374   return result;
2375 }
2376
2377
2378 /* Convert logical to integer.  */
2379
2380 gfc_expr *
2381 gfc_log2int (gfc_expr *src, int kind)
2382 {
2383   gfc_expr *result;
2384
2385   result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2386   mpz_set_si (result->value.integer, src->value.logical);
2387
2388   return result;
2389 }
2390
2391
2392 /* Convert integer to logical.  */
2393
2394 gfc_expr *
2395 gfc_int2log (gfc_expr *src, int kind)
2396 {
2397   gfc_expr *result;
2398
2399   result = gfc_constant_result (BT_LOGICAL, kind, &src->where);
2400   result->value.logical = (mpz_cmp_si (src->value.integer, 0) != 0);
2401
2402   return result;
2403 }
2404
2405
2406 /* Helper function to set the representation in a Hollerith conversion.  
2407    This assumes that the ts.type and ts.kind of the result have already
2408    been set.  */
2409
2410 static void
2411 hollerith2representation (gfc_expr *result, gfc_expr *src)
2412 {
2413   int src_len, result_len;
2414
2415   src_len = src->representation.length;
2416   result_len = gfc_target_expr_size (result);
2417
2418   if (src_len > result_len)
2419     {
2420       gfc_warning ("The Hollerith constant at %L is too long to convert to %s",
2421                    &src->where, gfc_typename(&result->ts));
2422     }
2423
2424   result->representation.string = XCNEWVEC (char, result_len + 1);
2425   memcpy (result->representation.string, src->representation.string,
2426           MIN (result_len, src_len));
2427
2428   if (src_len < result_len)
2429     memset (&result->representation.string[src_len], ' ', result_len - src_len);
2430
2431   result->representation.string[result_len] = '\0'; /* For debugger  */
2432   result->representation.length = result_len;
2433 }
2434
2435
2436 /* Convert Hollerith to integer. The constant will be padded or truncated.  */
2437
2438 gfc_expr *
2439 gfc_hollerith2int (gfc_expr *src, int kind)
2440 {
2441   gfc_expr *result;
2442
2443   result = gfc_get_expr ();
2444   result->expr_type = EXPR_CONSTANT;
2445   result->ts.type = BT_INTEGER;
2446   result->ts.kind = kind;
2447   result->where = src->where;
2448
2449   hollerith2representation (result, src);
2450   gfc_interpret_integer (kind, (unsigned char *) result->representation.string,
2451                          result->representation.length, result->value.integer);
2452
2453   return result;
2454 }
2455
2456
2457 /* Convert Hollerith to real. The constant will be padded or truncated.  */
2458
2459 gfc_expr *
2460 gfc_hollerith2real (gfc_expr *src, int kind)
2461 {
2462   gfc_expr *result;
2463   int len;
2464
2465   len = src->value.character.length;
2466
2467   result = gfc_get_expr ();
2468   result->expr_type = EXPR_CONSTANT;
2469   result->ts.type = BT_REAL;
2470   result->ts.kind = kind;
2471   result->where = src->where;
2472
2473   hollerith2representation (result, src);
2474   gfc_interpret_float (kind, (unsigned char *) result->representation.string,
2475                        result->representation.length, result->value.real);
2476
2477   return result;
2478 }
2479
2480
2481 /* Convert Hollerith to complex. The constant will be padded or truncated.  */
2482
2483 gfc_expr *
2484 gfc_hollerith2complex (gfc_expr *src, int kind)
2485 {
2486   gfc_expr *result;
2487   int len;
2488
2489   len = src->value.character.length;
2490
2491   result = gfc_get_expr ();
2492   result->expr_type = EXPR_CONSTANT;
2493   result->ts.type = BT_COMPLEX;
2494   result->ts.kind = kind;
2495   result->where = src->where;
2496
2497   hollerith2representation (result, src);
2498   gfc_interpret_complex (kind, (unsigned char *) result->representation.string,
2499                          result->representation.length, result->value.complex.r,
2500                          result->value.complex.i);
2501
2502   return result;
2503 }
2504
2505
2506 /* Convert Hollerith to character. */
2507
2508 gfc_expr *
2509 gfc_hollerith2character (gfc_expr *src, int kind)
2510 {
2511   gfc_expr *result;
2512
2513   result = gfc_copy_expr (src);
2514   result->ts.type = BT_CHARACTER;
2515   result->ts.kind = kind;
2516
2517   result->value.character.length = result->representation.length;
2518   result->value.character.string
2519     = gfc_char_to_widechar (result->representation.string);
2520
2521   return result;
2522 }
2523
2524
2525 /* Convert Hollerith to logical. The constant will be padded or truncated.  */
2526
2527 gfc_expr *
2528 gfc_hollerith2logical (gfc_expr *src, int kind)
2529 {
2530   gfc_expr *result;
2531   int len;
2532
2533   len = src->value.character.length;
2534
2535   result = gfc_get_expr ();
2536   result->expr_type = EXPR_CONSTANT;
2537   result->ts.type = BT_LOGICAL;
2538   result->ts.kind = kind;
2539   result->where = src->where;
2540
2541   hollerith2representation (result, src);
2542   gfc_interpret_logical (kind, (unsigned char *) result->representation.string,
2543                          result->representation.length, &result->value.logical);
2544
2545   return result;
2546 }
2547
2548
2549 /* Returns an initializer whose value is one higher than the value of the
2550    LAST_INITIALIZER argument.  If the argument is NULL, the
2551    initializers value will be set to zero.  The initializer's kind
2552    will be set to gfc_c_int_kind.
2553
2554    If -fshort-enums is given, the appropriate kind will be selected
2555    later after all enumerators have been parsed.  A warning is issued
2556    here if an initializer exceeds gfc_c_int_kind.  */
2557
2558 gfc_expr *
2559 gfc_enum_initializer (gfc_expr *last_initializer, locus where)
2560 {
2561   gfc_expr *result;
2562
2563   result = gfc_get_expr ();
2564   result->expr_type = EXPR_CONSTANT;
2565   result->ts.type = BT_INTEGER;
2566   result->ts.kind = gfc_c_int_kind;
2567   result->where = where;
2568
2569   mpz_init (result->value.integer);
2570
2571   if (last_initializer != NULL)
2572     {
2573       mpz_add_ui (result->value.integer, last_initializer->value.integer, 1);
2574       result->where = last_initializer->where;
2575
2576       if (gfc_check_integer_range (result->value.integer,
2577              gfc_c_int_kind) != ARITH_OK)
2578         {
2579           gfc_error ("Enumerator exceeds the C integer type at %C");
2580           return NULL;
2581         }
2582     }
2583   else
2584     {
2585       /* Control comes here, if it's the very first enumerator and no
2586          initializer has been given.  It will be initialized to zero.  */
2587       mpz_set_si (result->value.integer, 0);
2588     }
2589
2590   return result;
2591 }