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