gcc/fortran/arith.c

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