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   result = gfc_constant_result (BT_CHARACTER, gfc_default_character_kind,
1073                                 &op1->where);
1074
1075   len = op1->value.character.length + op2->value.character.length;
1076
1077   result->value.character.string = gfc_get_wide_string (len + 1);
1078   result->value.character.length = len;
1079
1080   memcpy (result->value.character.string, op1->value.character.string,
1081           op1->value.character.length * sizeof (gfc_char_t));
1082
1083   memcpy (&result->value.character.string[op1->value.character.length],
1084           op2->value.character.string,
1085           op2->value.character.length * sizeof (gfc_char_t));
1086
1087   result->value.character.string[len] = '\0';
1088
1089   *resultp = result;
1090
1091   return ARITH_OK;
1092 }
1093
1094 /* Comparison between real values; returns 0 if (op1 .op. op2) is true.
1095    This function mimics mpfr_cmp but takes NaN into account.  */
1096
1097 static int
1098 compare_real (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1099 {
1100   int rc;
1101   switch (op)
1102     {
1103       case INTRINSIC_EQ:
1104         rc = mpfr_equal_p (op1->value.real, op2->value.real) ? 0 : 1;
1105         break;
1106       case INTRINSIC_GT:
1107         rc = mpfr_greater_p (op1->value.real, op2->value.real) ? 1 : -1;
1108         break;
1109       case INTRINSIC_GE:
1110         rc = mpfr_greaterequal_p (op1->value.real, op2->value.real) ? 1 : -1;
1111         break;
1112       case INTRINSIC_LT:
1113         rc = mpfr_less_p (op1->value.real, op2->value.real) ? -1 : 1;
1114         break;
1115       case INTRINSIC_LE:
1116         rc = mpfr_lessequal_p (op1->value.real, op2->value.real) ? -1 : 1;
1117         break;
1118       default:
1119         gfc_internal_error ("compare_real(): Bad operator");
1120     }
1121
1122   return rc;
1123 }
1124
1125 /* Comparison operators.  Assumes that the two expression nodes
1126    contain two constants of the same type. The op argument is
1127    needed to handle NaN correctly.  */
1128
1129 int
1130 gfc_compare_expr (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1131 {
1132   int rc;
1133
1134   switch (op1->ts.type)
1135     {
1136     case BT_INTEGER:
1137       rc = mpz_cmp (op1->value.integer, op2->value.integer);
1138       break;
1139
1140     case BT_REAL:
1141       rc = compare_real (op1, op2, op);
1142       break;
1143
1144     case BT_CHARACTER:
1145       rc = gfc_compare_string (op1, op2);
1146       break;
1147
1148     case BT_LOGICAL:
1149       rc = ((!op1->value.logical && op2->value.logical)
1150             || (op1->value.logical && !op2->value.logical));
1151       break;
1152
1153     default:
1154       gfc_internal_error ("gfc_compare_expr(): Bad basic type");
1155     }
1156
1157   return rc;
1158 }
1159
1160
1161 /* Compare a pair of complex numbers.  Naturally, this is only for
1162    equality and inequality.  */
1163
1164 static int
1165 compare_complex (gfc_expr *op1, gfc_expr *op2)
1166 {
1167   return (mpfr_equal_p (op1->value.complex.r, op2->value.complex.r)
1168           && mpfr_equal_p (op1->value.complex.i, op2->value.complex.i));
1169 }
1170
1171
1172 /* Given two constant strings and the inverse collating sequence, compare the
1173    strings.  We return -1 for a < b, 0 for a == b and 1 for a > b.
1174    We use the processor's default collating sequence.  */
1175
1176 int
1177 gfc_compare_string (gfc_expr *a, gfc_expr *b)
1178 {
1179   int len, alen, blen, i;
1180   gfc_char_t ac, bc;
1181
1182   alen = a->value.character.length;
1183   blen = b->value.character.length;
1184
1185   len = MAX(alen, blen);
1186
1187   for (i = 0; i < len; i++)
1188     {
1189       ac = ((i < alen) ? a->value.character.string[i] : ' ');
1190       bc = ((i < blen) ? b->value.character.string[i] : ' ');
1191
1192       if (ac < bc)
1193         return -1;
1194       if (ac > bc)
1195         return 1;
1196     }
1197
1198   /* Strings are equal */
1199   return 0;
1200 }
1201
1202
1203 int
1204 gfc_compare_with_Cstring (gfc_expr *a, const char *b, bool case_sensitive)
1205 {
1206   int len, alen, blen, i;
1207   gfc_char_t ac, bc;
1208
1209   alen = a->value.character.length;
1210   blen = strlen (b);
1211
1212   len = MAX(alen, blen);
1213
1214   for (i = 0; i < len; i++)
1215     {
1216       ac = ((i < alen) ? a->value.character.string[i] : ' ');
1217       bc = ((i < blen) ? b[i] : ' ');
1218
1219       if (!case_sensitive)
1220         {
1221           ac = TOLOWER (ac);
1222           bc = TOLOWER (bc);
1223         }
1224
1225       if (ac < bc)
1226         return -1;
1227       if (ac > bc)
1228         return 1;
1229     }
1230
1231   /* Strings are equal */
1232   return 0;
1233 }
1234
1235
1236 /* Specific comparison subroutines.  */
1237
1238 static arith
1239 gfc_arith_eq (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1240 {
1241   gfc_expr *result;
1242
1243   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1244                                 &op1->where);
1245   result->value.logical = (op1->ts.type == BT_COMPLEX)
1246                         ? compare_complex (op1, op2)
1247                         : (gfc_compare_expr (op1, op2, INTRINSIC_EQ) == 0);
1248
1249   *resultp = result;
1250   return ARITH_OK;
1251 }
1252
1253
1254 static arith
1255 gfc_arith_ne (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1256 {
1257   gfc_expr *result;
1258
1259   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1260                                 &op1->where);
1261   result->value.logical = (op1->ts.type == BT_COMPLEX)
1262                         ? !compare_complex (op1, op2)
1263                         : (gfc_compare_expr (op1, op2, INTRINSIC_EQ) != 0);
1264
1265   *resultp = result;
1266   return ARITH_OK;
1267 }
1268
1269
1270 static arith
1271 gfc_arith_gt (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1272 {
1273   gfc_expr *result;
1274
1275   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1276                                 &op1->where);
1277   result->value.logical = (gfc_compare_expr (op1, op2, INTRINSIC_GT) > 0);
1278   *resultp = result;
1279
1280   return ARITH_OK;
1281 }
1282
1283
1284 static arith
1285 gfc_arith_ge (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1286 {
1287   gfc_expr *result;
1288
1289   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1290                                 &op1->where);
1291   result->value.logical = (gfc_compare_expr (op1, op2, INTRINSIC_GE) >= 0);
1292   *resultp = result;
1293
1294   return ARITH_OK;
1295 }
1296
1297
1298 static arith
1299 gfc_arith_lt (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1300 {
1301   gfc_expr *result;
1302
1303   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1304                                 &op1->where);
1305   result->value.logical = (gfc_compare_expr (op1, op2, INTRINSIC_LT) < 0);
1306   *resultp = result;
1307
1308   return ARITH_OK;
1309 }
1310
1311
1312 static arith
1313 gfc_arith_le (gfc_expr *op1, gfc_expr *op2, gfc_expr **resultp)
1314 {
1315   gfc_expr *result;
1316
1317   result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind,
1318                                 &op1->where);
1319   result->value.logical = (gfc_compare_expr (op1, op2, INTRINSIC_LE) <= 0);
1320   *resultp = result;
1321
1322   return ARITH_OK;
1323 }
1324
1325
1326 static arith
1327 reduce_unary (arith (*eval) (gfc_expr *, gfc_expr **), gfc_expr *op,
1328               gfc_expr **result)
1329 {
1330   gfc_constructor *c, *head;
1331   gfc_expr *r;
1332   arith rc;
1333
1334   if (op->expr_type == EXPR_CONSTANT)
1335     return eval (op, result);
1336
1337   rc = ARITH_OK;
1338   head = gfc_copy_constructor (op->value.constructor);
1339
1340   for (c = head; c; c = c->next)
1341     {
1342       rc = reduce_unary (eval, c->expr, &r);
1343
1344       if (rc != ARITH_OK)
1345         break;
1346
1347       gfc_replace_expr (c->expr, r);
1348     }
1349
1350   if (rc != ARITH_OK)
1351     gfc_free_constructor (head);
1352   else
1353     {
1354       r = gfc_get_expr ();
1355       r->expr_type = EXPR_ARRAY;
1356       r->value.constructor = head;
1357       r->shape = gfc_copy_shape (op->shape, op->rank);
1358
1359       r->ts = head->expr->ts;
1360       r->where = op->where;
1361       r->rank = op->rank;
1362
1363       *result = r;
1364     }
1365
1366   return rc;
1367 }
1368
1369
1370 static arith
1371 reduce_binary_ac (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1372                   gfc_expr *op1, gfc_expr *op2, gfc_expr **result)
1373 {
1374   gfc_constructor *c, *head;
1375   gfc_expr *r;
1376   arith rc;
1377
1378   head = gfc_copy_constructor (op1->value.constructor);
1379   rc = ARITH_OK;
1380
1381   for (c = head; c; c = c->next)
1382     {
1383       if (c->expr->expr_type == EXPR_CONSTANT)
1384         rc = eval (c->expr, op2, &r);
1385       else
1386         rc = reduce_binary_ac (eval, c->expr, op2, &r);
1387
1388       if (rc != ARITH_OK)
1389         break;
1390
1391       gfc_replace_expr (c->expr, r);
1392     }
1393
1394   if (rc != ARITH_OK)
1395     gfc_free_constructor (head);
1396   else
1397     {
1398       r = gfc_get_expr ();
1399       r->expr_type = EXPR_ARRAY;
1400       r->value.constructor = head;
1401       r->shape = gfc_copy_shape (op1->shape, op1->rank);
1402
1403       r->ts = head->expr->ts;
1404       r->where = op1->where;
1405       r->rank = op1->rank;
1406
1407       *result = r;
1408     }
1409
1410   return rc;
1411 }
1412
1413
1414 static arith
1415 reduce_binary_ca (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1416                   gfc_expr *op1, gfc_expr *op2, gfc_expr **result)
1417 {
1418   gfc_constructor *c, *head;
1419   gfc_expr *r;
1420   arith rc;
1421
1422   head = gfc_copy_constructor (op2->value.constructor);
1423   rc = ARITH_OK;
1424
1425   for (c = head; c; c = c->next)
1426     {
1427       if (c->expr->expr_type == EXPR_CONSTANT)
1428         rc = eval (op1, c->expr, &r);
1429       else
1430         rc = reduce_binary_ca (eval, op1, c->expr, &r);
1431
1432       if (rc != ARITH_OK)
1433         break;
1434
1435       gfc_replace_expr (c->expr, r);
1436     }
1437
1438   if (rc != ARITH_OK)
1439     gfc_free_constructor (head);
1440   else
1441     {
1442       r = gfc_get_expr ();
1443       r->expr_type = EXPR_ARRAY;
1444       r->value.constructor = head;
1445       r->shape = gfc_copy_shape (op2->shape, op2->rank);
1446
1447       r->ts = head->expr->ts;
1448       r->where = op2->where;
1449       r->rank = op2->rank;
1450
1451       *result = r;
1452     }
1453
1454   return rc;
1455 }
1456
1457
1458 /* We need a forward declaration of reduce_binary.  */
1459 static arith reduce_binary (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1460                             gfc_expr *op1, gfc_expr *op2, gfc_expr **result);
1461
1462
1463 static arith
1464 reduce_binary_aa (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1465                   gfc_expr *op1, gfc_expr *op2, gfc_expr **result)
1466 {
1467   gfc_constructor *c, *d, *head;
1468   gfc_expr *r;
1469   arith rc;
1470
1471   head = gfc_copy_constructor (op1->value.constructor);
1472
1473   rc = ARITH_OK;
1474   d = op2->value.constructor;
1475
1476   if (gfc_check_conformance ("elemental binary operation", op1, op2)
1477       != SUCCESS)
1478     rc = ARITH_INCOMMENSURATE;
1479   else
1480     {
1481       for (c = head; c; c = c->next, d = d->next)
1482         {
1483           if (d == NULL)
1484             {
1485               rc = ARITH_INCOMMENSURATE;
1486               break;
1487             }
1488
1489           rc = reduce_binary (eval, c->expr, d->expr, &r);
1490           if (rc != ARITH_OK)
1491             break;
1492
1493           gfc_replace_expr (c->expr, r);
1494         }
1495
1496       if (d != NULL)
1497         rc = ARITH_INCOMMENSURATE;
1498     }
1499
1500   if (rc != ARITH_OK)
1501     gfc_free_constructor (head);
1502   else
1503     {
1504       r = gfc_get_expr ();
1505       r->expr_type = EXPR_ARRAY;
1506       r->value.constructor = head;
1507       r->shape = gfc_copy_shape (op1->shape, op1->rank);
1508
1509       r->ts = head->expr->ts;
1510       r->where = op1->where;
1511       r->rank = op1->rank;
1512
1513       *result = r;
1514     }
1515
1516   return rc;
1517 }
1518
1519
1520 static arith
1521 reduce_binary (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1522                gfc_expr *op1, gfc_expr *op2, gfc_expr **result)
1523 {
1524   if (op1->expr_type == EXPR_CONSTANT && op2->expr_type == EXPR_CONSTANT)
1525     return eval (op1, op2, result);
1526
1527   if (op1->expr_type == EXPR_CONSTANT && op2->expr_type == EXPR_ARRAY)
1528     return reduce_binary_ca (eval, op1, op2, result);
1529
1530   if (op1->expr_type == EXPR_ARRAY && op2->expr_type == EXPR_CONSTANT)
1531     return reduce_binary_ac (eval, op1, op2, result);
1532
1533   return reduce_binary_aa (eval, op1, op2, result);
1534 }
1535
1536
1537 typedef union
1538 {
1539   arith (*f2)(gfc_expr *, gfc_expr **);
1540   arith (*f3)(gfc_expr *, gfc_expr *, gfc_expr **);
1541 }
1542 eval_f;
1543
1544 /* High level arithmetic subroutines.  These subroutines go into
1545    eval_intrinsic(), which can do one of several things to its
1546    operands.  If the operands are incompatible with the intrinsic
1547    operation, we return a node pointing to the operands and hope that
1548    an operator interface is found during resolution.
1549
1550    If the operands are compatible and are constants, then we try doing
1551    the arithmetic.  We also handle the cases where either or both
1552    operands are array constructors.  */
1553
1554 static gfc_expr *
1555 eval_intrinsic (gfc_intrinsic_op op,
1556                 eval_f eval, gfc_expr *op1, gfc_expr *op2)
1557 {
1558   gfc_expr temp, *result;
1559   int unary;
1560   arith rc;
1561
1562   gfc_clear_ts (&temp.ts);
1563
1564   switch (op)
1565     {
1566     /* Logical unary  */
1567     case INTRINSIC_NOT:
1568       if (op1->ts.type != BT_LOGICAL)
1569         goto runtime;
1570
1571       temp.ts.type = BT_LOGICAL;
1572       temp.ts.kind = gfc_default_logical_kind;
1573       unary = 1;
1574       break;
1575
1576     /* Logical binary operators  */
1577     case INTRINSIC_OR:
1578     case INTRINSIC_AND:
1579     case INTRINSIC_NEQV:
1580     case INTRINSIC_EQV:
1581       if (op1->ts.type != BT_LOGICAL || op2->ts.type != BT_LOGICAL)
1582         goto runtime;
1583
1584       temp.ts.type = BT_LOGICAL;
1585       temp.ts.kind = gfc_default_logical_kind;
1586       unary = 0;
1587       break;
1588
1589     /* Numeric unary  */
1590     case INTRINSIC_UPLUS:
1591     case INTRINSIC_UMINUS:
1592       if (!gfc_numeric_ts (&op1->ts))
1593         goto runtime;
1594
1595       temp.ts = op1->ts;
1596       unary = 1;
1597       break;
1598
1599     case INTRINSIC_PARENTHESES:
1600       temp.ts = op1->ts;
1601       unary = 1;
1602       break;
1603
1604     /* Additional restrictions for ordering relations.  */
1605     case INTRINSIC_GE:
1606     case INTRINSIC_GE_OS:
1607     case INTRINSIC_LT:
1608     case INTRINSIC_LT_OS:
1609     case INTRINSIC_LE:
1610     case INTRINSIC_LE_OS:
1611     case INTRINSIC_GT:
1612     case INTRINSIC_GT_OS:
1613       if (op1->ts.type == BT_COMPLEX || op2->ts.type == BT_COMPLEX)
1614         {
1615           temp.ts.type = BT_LOGICAL;
1616           temp.ts.kind = gfc_default_logical_kind;
1617           goto runtime;
1618         }
1619
1620     /* Fall through  */
1621     case INTRINSIC_EQ:
1622     case INTRINSIC_EQ_OS:
1623     case INTRINSIC_NE:
1624     case INTRINSIC_NE_OS:
1625       if (op1->ts.type == BT_CHARACTER && op2->ts.type == BT_CHARACTER)
1626         {
1627           unary = 0;
1628           temp.ts.type = BT_LOGICAL;
1629           temp.ts.kind = gfc_default_logical_kind;
1630
1631           /* If kind mismatch, exit and we'll error out later.  */
1632           if (op1->ts.kind != op2->ts.kind)
1633             goto runtime;
1634
1635           break;
1636         }
1637
1638     /* Fall through  */
1639     /* Numeric binary  */
1640     case INTRINSIC_PLUS:
1641     case INTRINSIC_MINUS:
1642     case INTRINSIC_TIMES:
1643     case INTRINSIC_DIVIDE:
1644     case INTRINSIC_POWER:
1645       if (!gfc_numeric_ts (&op1->ts) || !gfc_numeric_ts (&op2->ts))
1646         goto runtime;
1647
1648       /* Insert any necessary type conversions to make the operands
1649          compatible.  */
1650
1651       temp.expr_type = EXPR_OP;
1652       gfc_clear_ts (&temp.ts);
1653       temp.value.op.op = op;
1654
1655       temp.value.op.op1 = op1;
1656       temp.value.op.op2 = op2;
1657
1658       gfc_type_convert_binary (&temp);
1659
1660       if (op == INTRINSIC_EQ || op == INTRINSIC_NE
1661           || op == INTRINSIC_GE || op == INTRINSIC_GT
1662           || op == INTRINSIC_LE || op == INTRINSIC_LT
1663           || op == INTRINSIC_EQ_OS || op == INTRINSIC_NE_OS
1664           || op == INTRINSIC_GE_OS || op == INTRINSIC_GT_OS
1665           || op == INTRINSIC_LE_OS || op == INTRINSIC_LT_OS)
1666         {
1667           temp.ts.type = BT_LOGICAL;
1668           temp.ts.kind = gfc_default_logical_kind;
1669         }
1670
1671       unary = 0;
1672       break;
1673
1674     /* Character binary  */
1675     case INTRINSIC_CONCAT:
1676       if (op1->ts.type != BT_CHARACTER || op2->ts.type != BT_CHARACTER
1677           || op1->ts.kind != op2->ts.kind)
1678         goto runtime;
1679
1680       temp.ts.type = BT_CHARACTER;
1681       temp.ts.kind = op1->ts.kind;
1682       unary = 0;
1683       break;
1684
1685     case INTRINSIC_USER:
1686       goto runtime;
1687
1688     default:
1689       gfc_internal_error ("eval_intrinsic(): Bad operator");
1690     }
1691
1692   /* Try to combine the operators.  */
1693   if (op == INTRINSIC_POWER && op2->ts.type != BT_INTEGER)
1694     goto runtime;
1695
1696   if (op1->expr_type != EXPR_CONSTANT
1697       && (op1->expr_type != EXPR_ARRAY
1698           || !gfc_is_constant_expr (op1) || !gfc_expanded_ac (op1)))
1699     goto runtime;
1700
1701   if (op2 != NULL
1702       && op2->expr_type != EXPR_CONSTANT
1703          && (op2->expr_type != EXPR_ARRAY
1704              || !gfc_is_constant_expr (op2) || !gfc_expanded_ac (op2)))
1705     goto runtime;
1706
1707   if (unary)
1708     rc = reduce_unary (eval.f2, op1, &result);
1709   else
1710     rc = reduce_binary (eval.f3, op1, op2, &result);
1711
1712   if (rc != ARITH_OK)
1713     { /* Something went wrong.  */
1714       gfc_error (gfc_arith_error (rc), &op1->where);
1715       return NULL;
1716     }
1717
1718   gfc_free_expr (op1);
1719   gfc_free_expr (op2);
1720   return result;
1721
1722 runtime:
1723   /* Create a run-time expression.  */
1724   result = gfc_get_expr ();
1725   result->ts = temp.ts;
1726
1727   result->expr_type = EXPR_OP;
1728   result->value.op.op = op;
1729
1730   result->value.op.op1 = op1;
1731   result->value.op.op2 = op2;
1732
1733   result->where = op1->where;
1734
1735   return result;
1736 }
1737
1738
1739 /* Modify type of expression for zero size array.  */
1740
1741 static gfc_expr *
1742 eval_type_intrinsic0 (gfc_intrinsic_op iop, gfc_expr *op)
1743 {
1744   if (op == NULL)
1745     gfc_internal_error ("eval_type_intrinsic0(): op NULL");
1746
1747   switch (iop)
1748     {
1749     case INTRINSIC_GE:
1750     case INTRINSIC_GE_OS:
1751     case INTRINSIC_LT:
1752     case INTRINSIC_LT_OS:
1753     case INTRINSIC_LE:
1754     case INTRINSIC_LE_OS:
1755     case INTRINSIC_GT:
1756     case INTRINSIC_GT_OS:
1757     case INTRINSIC_EQ:
1758     case INTRINSIC_EQ_OS:
1759     case INTRINSIC_NE:
1760     case INTRINSIC_NE_OS:
1761       op->ts.type = BT_LOGICAL;
1762       op->ts.kind = gfc_default_logical_kind;
1763       break;
1764
1765     default:
1766       break;
1767     }
1768
1769   return op;
1770 }
1771
1772
1773 /* Return nonzero if the expression is a zero size array.  */
1774
1775 static int
1776 gfc_zero_size_array (gfc_expr *e)
1777 {
1778   if (e->expr_type != EXPR_ARRAY)
1779     return 0;
1780
1781   return e->value.constructor == NULL;
1782 }
1783
1784
1785 /* Reduce a binary expression where at least one of the operands
1786    involves a zero-length array.  Returns NULL if neither of the
1787    operands is a zero-length array.  */
1788
1789 static gfc_expr *
1790 reduce_binary0 (gfc_expr *op1, gfc_expr *op2)
1791 {
1792   if (gfc_zero_size_array (op1))
1793     {
1794       gfc_free_expr (op2);
1795       return op1;
1796     }
1797
1798   if (gfc_zero_size_array (op2))
1799     {
1800       gfc_free_expr (op1);
1801       return op2;
1802     }
1803
1804   return NULL;
1805 }
1806
1807
1808 static gfc_expr *
1809 eval_intrinsic_f2 (gfc_intrinsic_op op,
1810                    arith (*eval) (gfc_expr *, gfc_expr **),
1811                    gfc_expr *op1, gfc_expr *op2)
1812 {
1813   gfc_expr *result;
1814   eval_f f;
1815
1816   if (op2 == NULL)
1817     {
1818       if (gfc_zero_size_array (op1))
1819         return eval_type_intrinsic0 (op, op1);
1820     }
1821   else
1822     {
1823       result = reduce_binary0 (op1, op2);
1824       if (result != NULL)
1825         return eval_type_intrinsic0 (op, result);
1826     }
1827
1828   f.f2 = eval;
1829   return eval_intrinsic (op, f, op1, op2);
1830 }
1831
1832
1833 static gfc_expr *
1834 eval_intrinsic_f3 (gfc_intrinsic_op op,
1835                    arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1836                    gfc_expr *op1, gfc_expr *op2)
1837 {
1838   gfc_expr *result;
1839   eval_f f;
1840
1841   result = reduce_binary0 (op1, op2);
1842   if (result != NULL)
1843     return eval_type_intrinsic0(op, result);
1844
1845   f.f3 = eval;
1846   return eval_intrinsic (op, f, op1, op2);
1847 }
1848
1849
1850 gfc_expr *
1851 gfc_parentheses (gfc_expr *op)
1852 {
1853   if (gfc_is_constant_expr (op))
1854     return op;
1855
1856   return eval_intrinsic_f2 (INTRINSIC_PARENTHESES, gfc_arith_identity,
1857                             op, NULL);
1858 }
1859
1860 gfc_expr *
1861 gfc_uplus (gfc_expr *op)
1862 {
1863   return eval_intrinsic_f2 (INTRINSIC_UPLUS, gfc_arith_identity, op, NULL);
1864 }
1865
1866
1867 gfc_expr *
1868 gfc_uminus (gfc_expr *op)
1869 {
1870   return eval_intrinsic_f2 (INTRINSIC_UMINUS, gfc_arith_uminus, op, NULL);
1871 }
1872
1873
1874 gfc_expr *
1875 gfc_add (gfc_expr *op1, gfc_expr *op2)
1876 {
1877   return eval_intrinsic_f3 (INTRINSIC_PLUS, gfc_arith_plus, op1, op2);
1878 }
1879
1880
1881 gfc_expr *
1882 gfc_subtract (gfc_expr *op1, gfc_expr *op2)
1883 {
1884   return eval_intrinsic_f3 (INTRINSIC_MINUS, gfc_arith_minus, op1, op2);
1885 }
1886
1887
1888 gfc_expr *
1889 gfc_multiply (gfc_expr *op1, gfc_expr *op2)
1890 {
1891   return eval_intrinsic_f3 (INTRINSIC_TIMES, gfc_arith_times, op1, op2);
1892 }
1893
1894
1895 gfc_expr *
1896 gfc_divide (gfc_expr *op1, gfc_expr *op2)
1897 {
1898   return eval_intrinsic_f3 (INTRINSIC_DIVIDE, gfc_arith_divide, op1, op2);
1899 }
1900
1901
1902 gfc_expr *
1903 gfc_power (gfc_expr *op1, gfc_expr *op2)
1904 {
1905   return eval_intrinsic_f3 (INTRINSIC_POWER, gfc_arith_power, op1, op2);
1906 }
1907
1908
1909 gfc_expr *
1910 gfc_concat (gfc_expr *op1, gfc_expr *op2)
1911 {
1912   return eval_intrinsic_f3 (INTRINSIC_CONCAT, gfc_arith_concat, op1, op2);
1913 }
1914
1915
1916 gfc_expr *
1917 gfc_and (gfc_expr *op1, gfc_expr *op2)
1918 {
1919   return eval_intrinsic_f3 (INTRINSIC_AND, gfc_arith_and, op1, op2);
1920 }
1921
1922
1923 gfc_expr *
1924 gfc_or (gfc_expr *op1, gfc_expr *op2)
1925 {
1926   return eval_intrinsic_f3 (INTRINSIC_OR, gfc_arith_or, op1, op2);
1927 }
1928
1929
1930 gfc_expr *
1931 gfc_not (gfc_expr *op1)
1932 {
1933   return eval_intrinsic_f2 (INTRINSIC_NOT, gfc_arith_not, op1, NULL);
1934 }
1935
1936
1937 gfc_expr *
1938 gfc_eqv (gfc_expr *op1, gfc_expr *op2)
1939 {
1940   return eval_intrinsic_f3 (INTRINSIC_EQV, gfc_arith_eqv, op1, op2);
1941 }
1942
1943
1944 gfc_expr *
1945 gfc_neqv (gfc_expr *op1, gfc_expr *op2)
1946 {
1947   return eval_intrinsic_f3 (INTRINSIC_NEQV, gfc_arith_neqv, op1, op2);
1948 }
1949
1950
1951 gfc_expr *
1952 gfc_eq (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1953 {
1954   return eval_intrinsic_f3 (op, gfc_arith_eq, op1, op2);
1955 }
1956
1957
1958 gfc_expr *
1959 gfc_ne (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1960 {
1961   return eval_intrinsic_f3 (op, gfc_arith_ne, op1, op2);
1962 }
1963
1964
1965 gfc_expr *
1966 gfc_gt (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1967 {
1968   return eval_intrinsic_f3 (op, gfc_arith_gt, op1, op2);
1969 }
1970
1971
1972 gfc_expr *
1973 gfc_ge (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1974 {
1975   return eval_intrinsic_f3 (op, gfc_arith_ge, op1, op2);
1976 }
1977
1978
1979 gfc_expr *
1980 gfc_lt (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1981 {
1982   return eval_intrinsic_f3 (op, gfc_arith_lt, op1, op2);
1983 }
1984
1985
1986 gfc_expr *
1987 gfc_le (gfc_expr *op1, gfc_expr *op2, gfc_intrinsic_op op)
1988 {
1989   return eval_intrinsic_f3 (op, gfc_arith_le, op1, op2);
1990 }
1991
1992
1993 /* Convert an integer string to an expression node.  */
1994
1995 gfc_expr *
1996 gfc_convert_integer (const char *buffer, int kind, int radix, locus *where)
1997 {
1998   gfc_expr *e;
1999   const char *t;
2000
2001   e = gfc_constant_result (BT_INTEGER, kind, where);
2002   /* A leading plus is allowed, but not by mpz_set_str.  */
2003   if (buffer[0] == '+')
2004     t = buffer + 1;
2005   else
2006     t = buffer;
2007   mpz_set_str (e->value.integer, t, radix);
2008
2009   return e;
2010 }
2011
2012
2013 /* Convert a real string to an expression node.  */
2014
2015 gfc_expr *
2016 gfc_convert_real (const char *buffer, int kind, locus *where)
2017 {
2018   gfc_expr *e;
2019
2020   e = gfc_constant_result (BT_REAL, kind, where);
2021   mpfr_set_str (e->value.real, buffer, 10, GFC_RND_MODE);
2022
2023   return e;
2024 }
2025
2026
2027 /* Convert a pair of real, constant expression nodes to a single
2028    complex expression node.  */
2029
2030 gfc_expr *
2031 gfc_convert_complex (gfc_expr *real, gfc_expr *imag, int kind)
2032 {
2033   gfc_expr *e;
2034
2035   e = gfc_constant_result (BT_COMPLEX, kind, &real->where);
2036   mpfr_set (e->value.complex.r, real->value.real, GFC_RND_MODE);
2037   mpfr_set (e->value.complex.i, imag->value.real, GFC_RND_MODE);
2038
2039   return e;
2040 }
2041
2042
2043 /******* Simplification of intrinsic functions with constant arguments *****/
2044
2045
2046 /* Deal with an arithmetic error.  */
2047
2048 static void
2049 arith_error (arith rc, gfc_typespec *from, gfc_typespec *to, locus *where)
2050 {
2051   switch (rc)
2052     {
2053     case ARITH_OK:
2054       gfc_error ("Arithmetic OK converting %s to %s at %L",
2055                  gfc_typename (from), gfc_typename (to), where);
2056       break;
2057     case ARITH_OVERFLOW:
2058       gfc_error ("Arithmetic overflow converting %s to %s at %L. This check "
2059                  "can be disabled with the option -fno-range-check",
2060                  gfc_typename (from), gfc_typename (to), where);
2061       break;
2062     case ARITH_UNDERFLOW:
2063       gfc_error ("Arithmetic underflow converting %s to %s at %L",
2064                  gfc_typename (from), gfc_typename (to), where);
2065       break;
2066     case ARITH_NAN:
2067       gfc_error ("Arithmetic NaN converting %s to %s at %L",
2068                  gfc_typename (from), gfc_typename (to), where);
2069       break;
2070     case ARITH_DIV0:
2071       gfc_error ("Division by zero converting %s to %s at %L",
2072                  gfc_typename (from), gfc_typename (to), where);
2073       break;
2074     case ARITH_INCOMMENSURATE:
2075       gfc_error ("Array operands are incommensurate converting %s to %s at %L",
2076                  gfc_typename (from), gfc_typename (to), where);
2077       break;
2078     case ARITH_ASYMMETRIC:
2079       gfc_error ("Integer outside symmetric range implied by Standard Fortran"
2080                  " converting %s to %s at %L",
2081                  gfc_typename (from), gfc_typename (to), where);
2082       break;
2083     default:
2084       gfc_internal_error ("gfc_arith_error(): Bad error code");
2085     }
2086
2087   /* TODO: Do something about the error, i.e., throw exception, return
2088      NaN, etc.  */
2089 }
2090
2091
2092 /* Convert integers to integers.  */
2093
2094 gfc_expr *
2095 gfc_int2int (gfc_expr *src, int kind)
2096 {
2097   gfc_expr *result;
2098   arith rc;
2099
2100   result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2101
2102   mpz_set (result->value.integer, src->value.integer);
2103
2104   if ((rc = gfc_check_integer_range (result->value.integer, kind)) != ARITH_OK)
2105     {
2106       if (rc == ARITH_ASYMMETRIC)
2107         {
2108           gfc_warning (gfc_arith_error (rc), &src->where);
2109         }
2110       else
2111         {
2112           arith_error (rc, &src->ts, &result->ts, &src->where);
2113           gfc_free_expr (result);
2114           return NULL;
2115         }
2116     }
2117
2118   return result;
2119 }
2120
2121
2122 /* Convert integers to reals.  */
2123
2124 gfc_expr *
2125 gfc_int2real (gfc_expr *src, int kind)
2126 {
2127   gfc_expr *result;
2128   arith rc;
2129
2130   result = gfc_constant_result (BT_REAL, kind, &src->where);
2131
2132   mpfr_set_z (result->value.real, src->value.integer, GFC_RND_MODE);
2133
2134   if ((rc = gfc_check_real_range (result->value.real, kind)) != ARITH_OK)
2135     {
2136       arith_error (rc, &src->ts, &result->ts, &src->where);
2137       gfc_free_expr (result);
2138       return NULL;
2139     }
2140
2141   return result;
2142 }
2143
2144
2145 /* Convert default integer to default complex.  */
2146
2147 gfc_expr *
2148 gfc_int2complex (gfc_expr *src, int kind)
2149 {
2150   gfc_expr *result;
2151   arith rc;
2152
2153   result = gfc_constant_result (BT_COMPLEX, kind, &src->where);
2154
2155   mpfr_set_z (result->value.complex.r, src->value.integer, GFC_RND_MODE);
2156   mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
2157
2158   if ((rc = gfc_check_real_range (result->value.complex.r, kind)) != ARITH_OK)
2159     {
2160       arith_error (rc, &src->ts, &result->ts, &src->where);
2161       gfc_free_expr (result);
2162       return NULL;
2163     }
2164
2165   return result;
2166 }
2167
2168
2169 /* Convert default real to default integer.  */
2170
2171 gfc_expr *
2172 gfc_real2int (gfc_expr *src, int kind)
2173 {
2174   gfc_expr *result;
2175   arith rc;
2176
2177   result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2178
2179   gfc_mpfr_to_mpz (result->value.integer, src->value.real);
2180
2181   if ((rc = gfc_check_integer_range (result->value.integer, kind)) != ARITH_OK)
2182     {
2183       arith_error (rc, &src->ts, &result->ts, &src->where);
2184       gfc_free_expr (result);
2185       return NULL;
2186     }
2187
2188   return result;
2189 }
2190
2191
2192 /* Convert real to real.  */
2193
2194 gfc_expr *
2195 gfc_real2real (gfc_expr *src, int kind)
2196 {
2197   gfc_expr *result;
2198   arith rc;
2199
2200   result = gfc_constant_result (BT_REAL, kind, &src->where);
2201
2202   mpfr_set (result->value.real, src->value.real, GFC_RND_MODE);
2203
2204   rc = gfc_check_real_range (result->value.real, kind);
2205
2206   if (rc == ARITH_UNDERFLOW)
2207     {
2208       if (gfc_option.warn_underflow)
2209         gfc_warning (gfc_arith_error (rc), &src->where);
2210       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
2211     }
2212   else if (rc != ARITH_OK)
2213     {
2214       arith_error (rc, &src->ts, &result->ts, &src->where);
2215       gfc_free_expr (result);
2216       return NULL;
2217     }
2218
2219   return result;
2220 }
2221
2222
2223 /* Convert real to complex.  */
2224
2225 gfc_expr *
2226 gfc_real2complex (gfc_expr *src, int kind)
2227 {
2228   gfc_expr *result;
2229   arith rc;
2230
2231   result = gfc_constant_result (BT_COMPLEX, kind, &src->where);
2232
2233   mpfr_set (result->value.complex.r, src->value.real, GFC_RND_MODE);
2234   mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
2235
2236   rc = gfc_check_real_range (result->value.complex.r, kind);
2237
2238   if (rc == ARITH_UNDERFLOW)
2239     {
2240       if (gfc_option.warn_underflow)
2241         gfc_warning (gfc_arith_error (rc), &src->where);
2242       mpfr_set_ui (result->value.complex.r, 0, GFC_RND_MODE);
2243     }
2244   else if (rc != ARITH_OK)
2245     {
2246       arith_error (rc, &src->ts, &result->ts, &src->where);
2247       gfc_free_expr (result);
2248       return NULL;
2249     }
2250
2251   return result;
2252 }
2253
2254
2255 /* Convert complex to integer.  */
2256
2257 gfc_expr *
2258 gfc_complex2int (gfc_expr *src, int kind)
2259 {
2260   gfc_expr *result;
2261   arith rc;
2262
2263   result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2264
2265   gfc_mpfr_to_mpz (result->value.integer, src->value.complex.r);
2266
2267   if ((rc = gfc_check_integer_range (result->value.integer, kind)) != ARITH_OK)
2268     {
2269       arith_error (rc, &src->ts, &result->ts, &src->where);
2270       gfc_free_expr (result);
2271       return NULL;
2272     }
2273
2274   return result;
2275 }
2276
2277
2278 /* Convert complex to real.  */
2279
2280 gfc_expr *
2281 gfc_complex2real (gfc_expr *src, int kind)
2282 {
2283   gfc_expr *result;
2284   arith rc;
2285
2286   result = gfc_constant_result (BT_REAL, kind, &src->where);
2287
2288   mpfr_set (result->value.real, src->value.complex.r, GFC_RND_MODE);
2289
2290   rc = gfc_check_real_range (result->value.real, kind);
2291
2292   if (rc == ARITH_UNDERFLOW)
2293     {
2294       if (gfc_option.warn_underflow)
2295         gfc_warning (gfc_arith_error (rc), &src->where);
2296       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
2297     }
2298   if (rc != ARITH_OK)
2299     {
2300       arith_error (rc, &src->ts, &result->ts, &src->where);
2301       gfc_free_expr (result);
2302       return NULL;
2303     }
2304
2305   return result;
2306 }
2307
2308
2309 /* Convert complex to complex.  */
2310
2311 gfc_expr *
2312 gfc_complex2complex (gfc_expr *src, int kind)
2313 {
2314   gfc_expr *result;
2315   arith rc;
2316
2317   result = gfc_constant_result (BT_COMPLEX, kind, &src->where);
2318
2319   mpfr_set (result->value.complex.r, src->value.complex.r, GFC_RND_MODE);
2320   mpfr_set (result->value.complex.i, src->value.complex.i, GFC_RND_MODE);
2321
2322   rc = gfc_check_real_range (result->value.complex.r, kind);
2323
2324   if (rc == ARITH_UNDERFLOW)
2325     {
2326       if (gfc_option.warn_underflow)
2327         gfc_warning (gfc_arith_error (rc), &src->where);
2328       mpfr_set_ui (result->value.complex.r, 0, GFC_RND_MODE);
2329     }
2330   else if (rc != ARITH_OK)
2331     {
2332       arith_error (rc, &src->ts, &result->ts, &src->where);
2333       gfc_free_expr (result);
2334       return NULL;
2335     }
2336
2337   rc = gfc_check_real_range (result->value.complex.i, kind);
2338
2339   if (rc == ARITH_UNDERFLOW)
2340     {
2341       if (gfc_option.warn_underflow)
2342         gfc_warning (gfc_arith_error (rc), &src->where);
2343       mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
2344     }
2345   else if (rc != ARITH_OK)
2346     {
2347       arith_error (rc, &src->ts, &result->ts, &src->where);
2348       gfc_free_expr (result);
2349       return NULL;
2350     }
2351
2352   return result;
2353 }
2354
2355
2356 /* Logical kind conversion.  */
2357
2358 gfc_expr *
2359 gfc_log2log (gfc_expr *src, int kind)
2360 {
2361   gfc_expr *result;
2362
2363   result = gfc_constant_result (BT_LOGICAL, kind, &src->where);
2364   result->value.logical = src->value.logical;
2365
2366   return result;
2367 }
2368
2369
2370 /* Convert logical to integer.  */
2371
2372 gfc_expr *
2373 gfc_log2int (gfc_expr *src, int kind)
2374 {
2375   gfc_expr *result;
2376
2377   result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2378   mpz_set_si (result->value.integer, src->value.logical);
2379
2380   return result;
2381 }
2382
2383
2384 /* Convert integer to logical.  */
2385
2386 gfc_expr *
2387 gfc_int2log (gfc_expr *src, int kind)
2388 {
2389   gfc_expr *result;
2390
2391   result = gfc_constant_result (BT_LOGICAL, kind, &src->where);
2392   result->value.logical = (mpz_cmp_si (src->value.integer, 0) != 0);
2393
2394   return result;
2395 }
2396
2397
2398 /* Helper function to set the representation in a Hollerith conversion.
2399    This assumes that the ts.type and ts.kind of the result have already
2400    been set.  */
2401
2402 static void
2403 hollerith2representation (gfc_expr *result, gfc_expr *src)
2404 {
2405   int src_len, result_len;
2406
2407   src_len = src->representation.length;
2408   result_len = gfc_target_expr_size (result);
2409
2410   if (src_len > result_len)
2411     {
2412       gfc_warning ("The Hollerith constant at %L is too long to convert to %s",
2413                    &src->where, gfc_typename(&result->ts));
2414     }
2415
2416   result->representation.string = XCNEWVEC (char, result_len + 1);
2417   memcpy (result->representation.string, src->representation.string,
2418           MIN (result_len, src_len));
2419
2420   if (src_len < result_len)
2421     memset (&result->representation.string[src_len], ' ', result_len - src_len);
2422
2423   result->representation.string[result_len] = '\0'; /* For debugger  */
2424   result->representation.length = result_len;
2425 }
2426
2427
2428 /* Convert Hollerith to integer. The constant will be padded or truncated.  */
2429
2430 gfc_expr *
2431 gfc_hollerith2int (gfc_expr *src, int kind)
2432 {
2433   gfc_expr *result;
2434
2435   result = gfc_get_expr ();
2436   result->expr_type = EXPR_CONSTANT;
2437   result->ts.type = BT_INTEGER;
2438   result->ts.kind = kind;
2439   result->where = src->where;
2440
2441   hollerith2representation (result, src);
2442   gfc_interpret_integer (kind, (unsigned char *) result->representation.string,
2443                          result->representation.length, result->value.integer);
2444
2445   return result;
2446 }
2447
2448
2449 /* Convert Hollerith to real. The constant will be padded or truncated.  */
2450
2451 gfc_expr *
2452 gfc_hollerith2real (gfc_expr *src, int kind)
2453 {
2454   gfc_expr *result;
2455   int len;
2456
2457   len = src->value.character.length;
2458
2459   result = gfc_get_expr ();
2460   result->expr_type = EXPR_CONSTANT;
2461   result->ts.type = BT_REAL;
2462   result->ts.kind = kind;
2463   result->where = src->where;
2464
2465   hollerith2representation (result, src);
2466   gfc_interpret_float (kind, (unsigned char *) result->representation.string,
2467                        result->representation.length, result->value.real);
2468
2469   return result;
2470 }
2471
2472
2473 /* Convert Hollerith to complex. The constant will be padded or truncated.  */
2474
2475 gfc_expr *
2476 gfc_hollerith2complex (gfc_expr *src, int kind)
2477 {
2478   gfc_expr *result;
2479   int len;
2480
2481   len = src->value.character.length;
2482
2483   result = gfc_get_expr ();
2484   result->expr_type = EXPR_CONSTANT;
2485   result->ts.type = BT_COMPLEX;
2486   result->ts.kind = kind;
2487   result->where = src->where;
2488
2489   hollerith2representation (result, src);
2490   gfc_interpret_complex (kind, (unsigned char *) result->representation.string,
2491                          result->representation.length, result->value.complex.r,
2492                          result->value.complex.i);
2493
2494   return result;
2495 }
2496
2497
2498 /* Convert Hollerith to character. */
2499
2500 gfc_expr *
2501 gfc_hollerith2character (gfc_expr *src, int kind)
2502 {
2503   gfc_expr *result;
2504
2505   result = gfc_copy_expr (src);
2506   result->ts.type = BT_CHARACTER;
2507   result->ts.kind = kind;
2508
2509   result->value.character.length = result->representation.length;
2510   result->value.character.string
2511     = gfc_char_to_widechar (result->representation.string);
2512
2513   return result;
2514 }
2515
2516
2517 /* Convert Hollerith to logical. The constant will be padded or truncated.  */
2518
2519 gfc_expr *
2520 gfc_hollerith2logical (gfc_expr *src, int kind)
2521 {
2522   gfc_expr *result;
2523   int len;
2524
2525   len = src->value.character.length;
2526
2527   result = gfc_get_expr ();
2528   result->expr_type = EXPR_CONSTANT;
2529   result->ts.type = BT_LOGICAL;
2530   result->ts.kind = kind;
2531   result->where = src->where;
2532
2533   hollerith2representation (result, src);
2534   gfc_interpret_logical (kind, (unsigned char *) result->representation.string,
2535                          result->representation.length, &result->value.logical);
2536
2537   return result;
2538 }
2539
2540
2541 /* Returns an initializer whose value is one higher than the value of the
2542    LAST_INITIALIZER argument.  If the argument is NULL, the
2543    initializers value will be set to zero.  The initializer's kind
2544    will be set to gfc_c_int_kind.
2545
2546    If -fshort-enums is given, the appropriate kind will be selected
2547    later after all enumerators have been parsed.  A warning is issued
2548    here if an initializer exceeds gfc_c_int_kind.  */
2549
2550 gfc_expr *
2551 gfc_enum_initializer (gfc_expr *last_initializer, locus where)
2552 {
2553   gfc_expr *result;
2554
2555   result = gfc_get_expr ();
2556   result->expr_type = EXPR_CONSTANT;
2557   result->ts.type = BT_INTEGER;
2558   result->ts.kind = gfc_c_int_kind;
2559   result->where = where;
2560
2561   mpz_init (result->value.integer);
2562
2563   if (last_initializer != NULL)
2564     {
2565       mpz_add_ui (result->value.integer, last_initializer->value.integer, 1);
2566       result->where = last_initializer->where;
2567
2568       if (gfc_check_integer_range (result->value.integer,
2569              gfc_c_int_kind) != ARITH_OK)
2570         {
2571           gfc_error ("Enumerator exceeds the C integer type at %C");
2572           return NULL;
2573         }
2574     }
2575   else
2576     {
2577       /* Control comes here, if it's the very first enumerator and no
2578          initializer has been given.  It will be initialized to zero.  */
2579       mpz_set_si (result->value.integer, 0);
2580     }
2581
2582   return result;
2583 }