OSDN Git Service

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