OSDN Git Service

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