gcc/fortran/arith.c

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