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