gcc/fortran/simplify.c

   1 /* Simplify intrinsic functions at compile-time.
   2    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005 Free Software Foundation,
   3    Inc.
   4    Contributed by Andy Vaught & Katherine Holcomb
   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 #include "config.h"
  24 #include "system.h"
  25 #include "flags.h"
  26 #include "gfortran.h"
  27 #include "arith.h"
  28 #include "intrinsic.h"
  29
  30 gfc_expr gfc_bad_expr;
  31
  32
  33 /* Note that 'simplification' is not just transforming expressions.
  34    For functions that are not simplified at compile time, range
  35    checking is done if possible.
  36
  37    The return convention is that each simplification function returns:
  38
  39      A new expression node corresponding to the simplified arguments.
  40      The original arguments are destroyed by the caller, and must not
  41      be a part of the new expression.
  42
  43      NULL pointer indicating that no simplification was possible and
  44      the original expression should remain intact.  If the
  45      simplification function sets the type and/or the function name
  46      via the pointer gfc_simple_expression, then this type is
  47      retained.
  48
  49      An expression pointer to gfc_bad_expr (a static placeholder)
  50      indicating that some error has prevented simplification.  For
  51      example, sqrt(-1.0).  The error is generated within the function
  52      and should be propagated upwards
  53
  54    By the time a simplification function gets control, it has been
  55    decided that the function call is really supposed to be the
  56    intrinsic.  No type checking is strictly necessary, since only
  57    valid types will be passed on.  On the other hand, a simplification
  58    subroutine may have to look at the type of an argument as part of
  59    its processing.
  60
  61    Array arguments are never passed to these subroutines.
  62
  63    The functions in this file don't have much comment with them, but
  64    everything is reasonably straight-forward.  The Standard, chapter 13
  65    is the best comment you'll find for this file anyway.  */
  66
  67 /* Static table for converting non-ascii character sets to ascii.
  68    The xascii_table[] is the inverse table.  */
  69
  70 static int ascii_table[256] = {
  71   '\0', '\0', '\0', '\0', '\0', '\0', '\0', '\0',
  72   '\b', '\t', '\n', '\v', '\0', '\r', '\0', '\0',
  73   '\0', '\0', '\0', '\0', '\0', '\0', '\0', '\0',
  74   '\0', '\0', '\0', '\0', '\0', '\0', '\0', '\0',
  75   ' ', '!', '\'', '#', '$', '%', '&', '\'',
  76   '(', ')', '*', '+', ',', '-', '.', '/',
  77   '0', '1', '2', '3', '4', '5', '6', '7',
  78   '8', '9', ':', ';', '<', '=', '>', '?',
  79   '@', 'A', 'B', 'C', 'D', 'E', 'F', 'G',
  80   'H', 'I', 'J', 'K', 'L', 'M', 'N', 'O',
  81   'P', 'Q', 'R', 'S', 'T', 'U', 'V', 'W',
  82   'X', 'Y', 'Z', '[', '\\', ']', '^', '_',
  83   '`', 'a', 'b', 'c', 'd', 'e', 'f', 'g',
  84   'h', 'i', 'j', 'k', 'l', 'm', 'n', 'o',
  85   'p', 'q', 'r', 's', 't', 'u', 'v', 'w',
  86   'x', 'y', 'z', '{', '|', '}', '~', '\?'
  87 };
  88
  89 static int xascii_table[256];
  90
  91
  92 /* Range checks an expression node.  If all goes well, returns the
  93    node, otherwise returns &gfc_bad_expr and frees the node.  */
  94
  95 static gfc_expr *
  96 range_check (gfc_expr * result, const char *name)
  97 {
  98   if (gfc_range_check (result) == ARITH_OK)
  99     return result;
 100
 101   gfc_error ("Result of %s overflows its kind at %L", name, &result->where);
 102   gfc_free_expr (result);
 103   return &gfc_bad_expr;
 104 }
 105
 106
 107 /* A helper function that gets an optional and possibly missing
 108    kind parameter.  Returns the kind, -1 if something went wrong.  */
 109
 110 static int
 111 get_kind (bt type, gfc_expr * k, const char *name, int default_kind)
 112 {
 113   int kind;
 114
 115   if (k == NULL)
 116     return default_kind;
 117
 118   if (k->expr_type != EXPR_CONSTANT)
 119     {
 120       gfc_error ("KIND parameter of %s at %L must be an initialization "
 121                  "expression", name, &k->where);
 122
 123       return -1;
 124     }
 125
 126   if (gfc_extract_int (k, &kind) != NULL
 127       || gfc_validate_kind (type, kind, true) < 0)
 128     {
 129
 130       gfc_error ("Invalid KIND parameter of %s at %L", name, &k->where);
 131       return -1;
 132     }
 133
 134   return kind;
 135 }
 136
 137
 138 /* Checks if X, which is assumed to represent a two's complement
 139    integer of binary width BITSIZE, has the signbit set.  If so, makes
 140    X the corresponding negative number.  */
 141
 142 static void
 143 twos_complement (mpz_t x, int bitsize)
 144 {
 145   mpz_t mask;
 146
 147   if (mpz_tstbit (x, bitsize - 1) == 1)
 148     {
 149       mpz_init_set_ui(mask, 1);
 150       mpz_mul_2exp(mask, mask, bitsize);
 151       mpz_sub_ui(mask, mask, 1);
 152
 153       /* We negate the number by hand, zeroing the high bits, that is
 154         make it the corresponding positive number, and then have it
 155         negated by GMP, giving the correct representation of the
 156         negative number.  */
 157       mpz_com (x, x);
 158       mpz_add_ui (x, x, 1);
 159       mpz_and (x, x, mask);
 160
 161       mpz_neg (x, x);
 162
 163       mpz_clear (mask);
 164     }
 165 }
 166
 167
 168 /********************** Simplification functions *****************************/
 169
 170 gfc_expr *
 171 gfc_simplify_abs (gfc_expr * e)
 172 {
 173   gfc_expr *result;
 174
 175   if (e->expr_type != EXPR_CONSTANT)
 176     return NULL;
 177
 178   switch (e->ts.type)
 179     {
 180     case BT_INTEGER:
 181       result = gfc_constant_result (BT_INTEGER, e->ts.kind, &e->where);
 182
 183       mpz_abs (result->value.integer, e->value.integer);
 184
 185       result = range_check (result, "IABS");
 186       break;
 187
 188     case BT_REAL:
 189       result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
 190
 191       mpfr_abs (result->value.real, e->value.real, GFC_RND_MODE);
 192
 193       result = range_check (result, "ABS");
 194       break;
 195
 196     case BT_COMPLEX:
 197       result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
 198
 199       gfc_set_model_kind (e->ts.kind);
 200
 201       mpfr_hypot (result->value.real, e->value.complex.r,
 202                   e->value.complex.i, GFC_RND_MODE);
 203       result = range_check (result, "CABS");
 204       break;
 205
 206     default:
 207       gfc_internal_error ("gfc_simplify_abs(): Bad type");
 208     }
 209
 210   return result;
 211 }
 212
 213
 214 gfc_expr *
 215 gfc_simplify_achar (gfc_expr * e)
 216 {
 217   gfc_expr *result;
 218   int index;
 219
 220   if (e->expr_type != EXPR_CONSTANT)
 221     return NULL;
 222
 223   /* We cannot assume that the native character set is ASCII in this
 224      function.  */
 225   if (gfc_extract_int (e, &index) != NULL || index < 0 || index > 127)
 226     {
 227       gfc_error ("Extended ASCII not implemented: argument of ACHAR at %L "
 228                  "must be between 0 and 127", &e->where);
 229       return &gfc_bad_expr;
 230     }
 231
 232   result = gfc_constant_result (BT_CHARACTER, gfc_default_character_kind,
 233                                 &e->where);
 234
 235   result->value.character.string = gfc_getmem (2);
 236
 237   result->value.character.length = 1;
 238   result->value.character.string[0] = ascii_table[index];
 239   result->value.character.string[1] = '\0';     /* For debugger */
 240   return result;
 241 }
 242
 243
 244 gfc_expr *
 245 gfc_simplify_acos (gfc_expr * x)
 246 {
 247   gfc_expr *result;
 248
 249   if (x->expr_type != EXPR_CONSTANT)
 250     return NULL;
 251
 252   if (mpfr_cmp_si (x->value.real, 1) > 0 || mpfr_cmp_si (x->value.real, -1) < 0)
 253     {
 254       gfc_error ("Argument of ACOS at %L must be between -1 and 1",
 255                  &x->where);
 256       return &gfc_bad_expr;
 257     }
 258
 259   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 260
 261   mpfr_acos (result->value.real, x->value.real, GFC_RND_MODE);
 262
 263   return range_check (result, "ACOS");
 264 }
 265
 266 gfc_expr *
 267 gfc_simplify_acosh (gfc_expr * x)
 268 {
 269   gfc_expr *result;
 270
 271   if (x->expr_type != EXPR_CONSTANT)
 272     return NULL;
 273
 274   if (mpfr_cmp_si (x->value.real, 1) < 0)
 275     {
 276       gfc_error ("Argument of ACOSH at %L must not be less than 1",
 277                  &x->where);
 278       return &gfc_bad_expr;
 279     }
 280
 281   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 282
 283   mpfr_acosh (result->value.real, x->value.real, GFC_RND_MODE);
 284
 285   return range_check (result, "ACOSH");
 286 }
 287
 288 gfc_expr *
 289 gfc_simplify_adjustl (gfc_expr * e)
 290 {
 291   gfc_expr *result;
 292   int count, i, len;
 293   char ch;
 294
 295   if (e->expr_type != EXPR_CONSTANT)
 296     return NULL;
 297
 298   len = e->value.character.length;
 299
 300   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
 301
 302   result->value.character.length = len;
 303   result->value.character.string = gfc_getmem (len + 1);
 304
 305   for (count = 0, i = 0; i < len; ++i)
 306     {
 307       ch = e->value.character.string[i];
 308       if (ch != ' ')
 309         break;
 310       ++count;
 311     }
 312
 313   for (i = 0; i < len - count; ++i)
 314     {
 315       result->value.character.string[i] =
 316         e->value.character.string[count + i];
 317     }
 318
 319   for (i = len - count; i < len; ++i)
 320     {
 321       result->value.character.string[i] = ' ';
 322     }
 323
 324   result->value.character.string[len] = '\0';   /* For debugger */
 325
 326   return result;
 327 }
 328
 329
 330 gfc_expr *
 331 gfc_simplify_adjustr (gfc_expr * e)
 332 {
 333   gfc_expr *result;
 334   int count, i, len;
 335   char ch;
 336
 337   if (e->expr_type != EXPR_CONSTANT)
 338     return NULL;
 339
 340   len = e->value.character.length;
 341
 342   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
 343
 344   result->value.character.length = len;
 345   result->value.character.string = gfc_getmem (len + 1);
 346
 347   for (count = 0, i = len - 1; i >= 0; --i)
 348     {
 349       ch = e->value.character.string[i];
 350       if (ch != ' ')
 351         break;
 352       ++count;
 353     }
 354
 355   for (i = 0; i < count; ++i)
 356     {
 357       result->value.character.string[i] = ' ';
 358     }
 359
 360   for (i = count; i < len; ++i)
 361     {
 362       result->value.character.string[i] =
 363         e->value.character.string[i - count];
 364     }
 365
 366   result->value.character.string[len] = '\0';   /* For debugger */
 367
 368   return result;
 369 }
 370
 371
 372 gfc_expr *
 373 gfc_simplify_aimag (gfc_expr * e)
 374 {
 375   gfc_expr *result;
 376
 377   if (e->expr_type != EXPR_CONSTANT)
 378     return NULL;
 379
 380   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
 381   mpfr_set (result->value.real, e->value.complex.i, GFC_RND_MODE);
 382
 383   return range_check (result, "AIMAG");
 384 }
 385
 386
 387 gfc_expr *
 388 gfc_simplify_aint (gfc_expr * e, gfc_expr * k)
 389 {
 390   gfc_expr *rtrunc, *result;
 391   int kind;
 392
 393   kind = get_kind (BT_REAL, k, "AINT", e->ts.kind);
 394   if (kind == -1)
 395     return &gfc_bad_expr;
 396
 397   if (e->expr_type != EXPR_CONSTANT)
 398     return NULL;
 399
 400   rtrunc = gfc_copy_expr (e);
 401
 402   mpfr_trunc (rtrunc->value.real, e->value.real);
 403
 404   result = gfc_real2real (rtrunc, kind);
 405   gfc_free_expr (rtrunc);
 406
 407   return range_check (result, "AINT");
 408 }
 409
 410
 411 gfc_expr *
 412 gfc_simplify_dint (gfc_expr * e)
 413 {
 414   gfc_expr *rtrunc, *result;
 415
 416   if (e->expr_type != EXPR_CONSTANT)
 417     return NULL;
 418
 419   rtrunc = gfc_copy_expr (e);
 420
 421   mpfr_trunc (rtrunc->value.real, e->value.real);
 422
 423   result = gfc_real2real (rtrunc, gfc_default_double_kind);
 424   gfc_free_expr (rtrunc);
 425
 426   return range_check (result, "DINT");
 427 }
 428
 429
 430 gfc_expr *
 431 gfc_simplify_anint (gfc_expr * e, gfc_expr * k)
 432 {
 433   gfc_expr *result;
 434   int kind;
 435
 436   kind = get_kind (BT_REAL, k, "ANINT", e->ts.kind);
 437   if (kind == -1)
 438     return &gfc_bad_expr;
 439
 440   if (e->expr_type != EXPR_CONSTANT)
 441     return NULL;
 442
 443   result = gfc_constant_result (e->ts.type, kind, &e->where);
 444
 445   mpfr_round (result->value.real, e->value.real);
 446
 447   return range_check (result, "ANINT");
 448 }
 449
 450
 451 gfc_expr *
 452 gfc_simplify_dnint (gfc_expr * e)
 453 {
 454   gfc_expr *result;
 455
 456   if (e->expr_type != EXPR_CONSTANT)
 457     return NULL;
 458
 459   result = gfc_constant_result (BT_REAL, gfc_default_double_kind, &e->where);
 460
 461   mpfr_round (result->value.real, e->value.real);
 462
 463   return range_check (result, "DNINT");
 464 }
 465
 466
 467 gfc_expr *
 468 gfc_simplify_asin (gfc_expr * x)
 469 {
 470   gfc_expr *result;
 471
 472   if (x->expr_type != EXPR_CONSTANT)
 473     return NULL;
 474
 475   if (mpfr_cmp_si (x->value.real, 1) > 0 || mpfr_cmp_si (x->value.real, -1) < 0)
 476     {
 477       gfc_error ("Argument of ASIN at %L must be between -1 and 1",
 478                  &x->where);
 479       return &gfc_bad_expr;
 480     }
 481
 482   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 483
 484   mpfr_asin(result->value.real, x->value.real, GFC_RND_MODE);
 485
 486   return range_check (result, "ASIN");
 487 }
 488
 489
 490 gfc_expr *
 491 gfc_simplify_asinh (gfc_expr * x)
 492 {
 493   gfc_expr *result;
 494
 495   if (x->expr_type != EXPR_CONSTANT)
 496     return NULL;
 497
 498   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 499
 500   mpfr_asinh(result->value.real, x->value.real, GFC_RND_MODE);
 501
 502   return range_check (result, "ASINH");
 503 }
 504
 505
 506 gfc_expr *
 507 gfc_simplify_atan (gfc_expr * x)
 508 {
 509   gfc_expr *result;
 510
 511   if (x->expr_type != EXPR_CONSTANT)
 512     return NULL;
 513
 514   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 515
 516   mpfr_atan(result->value.real, x->value.real, GFC_RND_MODE);
 517
 518   return range_check (result, "ATAN");
 519 }
 520
 521
 522 gfc_expr *
 523 gfc_simplify_atanh (gfc_expr * x)
 524 {
 525   gfc_expr *result;
 526
 527   if (x->expr_type != EXPR_CONSTANT)
 528     return NULL;
 529
 530   if (mpfr_cmp_si (x->value.real, 1) >= 0 ||
 531       mpfr_cmp_si (x->value.real, -1) <= 0)
 532     {
 533       gfc_error ("Argument of ATANH at %L must be inside the range -1 to 1",
 534                  &x->where);
 535       return &gfc_bad_expr;
 536     }
 537
 538   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 539
 540   mpfr_atanh(result->value.real, x->value.real, GFC_RND_MODE);
 541
 542   return range_check (result, "ATANH");
 543 }
 544
 545
 546 gfc_expr *
 547 gfc_simplify_atan2 (gfc_expr * y, gfc_expr * x)
 548 {
 549   gfc_expr *result;
 550
 551   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
 552     return NULL;
 553
 554   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 555
 556   if (mpfr_sgn (y->value.real) == 0 && mpfr_sgn (x->value.real) == 0)
 557     {
 558       gfc_error
 559         ("If first argument of ATAN2 %L is zero, then the second argument "
 560           "must not be zero", &x->where);
 561       gfc_free_expr (result);
 562       return &gfc_bad_expr;
 563     }
 564
 565   arctangent2 (y->value.real, x->value.real, result->value.real);
 566
 567   return range_check (result, "ATAN2");
 568 }
 569
 570
 571 gfc_expr *
 572 gfc_simplify_bit_size (gfc_expr * e)
 573 {
 574   gfc_expr *result;
 575   int i;
 576
 577   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
 578   result = gfc_constant_result (BT_INTEGER, e->ts.kind, &e->where);
 579   mpz_set_ui (result->value.integer, gfc_integer_kinds[i].bit_size);
 580
 581   return result;
 582 }
 583
 584
 585 gfc_expr *
 586 gfc_simplify_btest (gfc_expr * e, gfc_expr * bit)
 587 {
 588   int b;
 589
 590   if (e->expr_type != EXPR_CONSTANT || bit->expr_type != EXPR_CONSTANT)
 591     return NULL;
 592
 593   if (gfc_extract_int (bit, &b) != NULL || b < 0)
 594     return gfc_logical_expr (0, &e->where);
 595
 596   return gfc_logical_expr (mpz_tstbit (e->value.integer, b), &e->where);
 597 }
 598
 599
 600 gfc_expr *
 601 gfc_simplify_ceiling (gfc_expr * e, gfc_expr * k)
 602 {
 603   gfc_expr *ceil, *result;
 604   int kind;
 605
 606   kind = get_kind (BT_INTEGER, k, "CEILING", gfc_default_integer_kind);
 607   if (kind == -1)
 608     return &gfc_bad_expr;
 609
 610   if (e->expr_type != EXPR_CONSTANT)
 611     return NULL;
 612
 613   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
 614
 615   ceil = gfc_copy_expr (e);
 616
 617   mpfr_ceil (ceil->value.real, e->value.real);
 618   gfc_mpfr_to_mpz(result->value.integer, ceil->value.real);
 619
 620   gfc_free_expr (ceil);
 621
 622   return range_check (result, "CEILING");
 623 }
 624
 625
 626 gfc_expr *
 627 gfc_simplify_char (gfc_expr * e, gfc_expr * k)
 628 {
 629   gfc_expr *result;
 630   int c, kind;
 631
 632   kind = get_kind (BT_CHARACTER, k, "CHAR", gfc_default_character_kind);
 633   if (kind == -1)
 634     return &gfc_bad_expr;
 635
 636   if (e->expr_type != EXPR_CONSTANT)
 637     return NULL;
 638
 639   if (gfc_extract_int (e, &c) != NULL || c < 0 || c > 255)
 640     {
 641       gfc_error ("Bad character in CHAR function at %L", &e->where);
 642       return &gfc_bad_expr;
 643     }
 644
 645   result = gfc_constant_result (BT_CHARACTER, kind, &e->where);
 646
 647   result->value.character.length = 1;
 648   result->value.character.string = gfc_getmem (2);
 649
 650   result->value.character.string[0] = c;
 651   result->value.character.string[1] = '\0';     /* For debugger */
 652
 653   return result;
 654 }
 655
 656
 657 /* Common subroutine for simplifying CMPLX and DCMPLX.  */
 658
 659 static gfc_expr *
 660 simplify_cmplx (const char *name, gfc_expr * x, gfc_expr * y, int kind)
 661 {
 662   gfc_expr *result;
 663
 664   result = gfc_constant_result (BT_COMPLEX, kind, &x->where);
 665
 666   mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
 667
 668   switch (x->ts.type)
 669     {
 670     case BT_INTEGER:
 671       mpfr_set_z (result->value.complex.r, x->value.integer, GFC_RND_MODE);
 672       break;
 673
 674     case BT_REAL:
 675       mpfr_set (result->value.complex.r, x->value.real, GFC_RND_MODE);
 676       break;
 677
 678     case BT_COMPLEX:
 679       mpfr_set (result->value.complex.r, x->value.complex.r, GFC_RND_MODE);
 680       mpfr_set (result->value.complex.i, x->value.complex.i, GFC_RND_MODE);
 681       break;
 682
 683     default:
 684       gfc_internal_error ("gfc_simplify_dcmplx(): Bad type (x)");
 685     }
 686
 687   if (y != NULL)
 688     {
 689       switch (y->ts.type)
 690         {
 691         case BT_INTEGER:
 692           mpfr_set_z (result->value.complex.i, y->value.integer, GFC_RND_MODE);
 693           break;
 694
 695         case BT_REAL:
 696           mpfr_set (result->value.complex.i, y->value.real, GFC_RND_MODE);
 697           break;
 698
 699         default:
 700           gfc_internal_error ("gfc_simplify_dcmplx(): Bad type (y)");
 701         }
 702     }
 703
 704   return range_check (result, name);
 705 }
 706
 707
 708 gfc_expr *
 709 gfc_simplify_cmplx (gfc_expr * x, gfc_expr * y, gfc_expr * k)
 710 {
 711   int kind;
 712
 713   if (x->expr_type != EXPR_CONSTANT
 714       || (y != NULL && y->expr_type != EXPR_CONSTANT))
 715     return NULL;
 716
 717   kind = get_kind (BT_REAL, k, "CMPLX", gfc_default_real_kind);
 718   if (kind == -1)
 719     return &gfc_bad_expr;
 720
 721   return simplify_cmplx ("CMPLX", x, y, kind);
 722 }
 723
 724
 725 gfc_expr *
 726 gfc_simplify_conjg (gfc_expr * e)
 727 {
 728   gfc_expr *result;
 729
 730   if (e->expr_type != EXPR_CONSTANT)
 731     return NULL;
 732
 733   result = gfc_copy_expr (e);
 734   mpfr_neg (result->value.complex.i, result->value.complex.i, GFC_RND_MODE);
 735
 736   return range_check (result, "CONJG");
 737 }
 738
 739
 740 gfc_expr *
 741 gfc_simplify_cos (gfc_expr * x)
 742 {
 743   gfc_expr *result;
 744   mpfr_t xp, xq;
 745
 746   if (x->expr_type != EXPR_CONSTANT)
 747     return NULL;
 748
 749   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 750
 751   switch (x->ts.type)
 752     {
 753     case BT_REAL:
 754       mpfr_cos (result->value.real, x->value.real, GFC_RND_MODE);
 755       break;
 756     case BT_COMPLEX:
 757       gfc_set_model_kind (x->ts.kind);
 758       mpfr_init (xp);
 759       mpfr_init (xq);
 760
 761       mpfr_cos  (xp, x->value.complex.r, GFC_RND_MODE);
 762       mpfr_cosh (xq, x->value.complex.i, GFC_RND_MODE);
 763       mpfr_mul(result->value.complex.r, xp, xq, GFC_RND_MODE);
 764
 765       mpfr_sin  (xp, x->value.complex.r, GFC_RND_MODE);
 766       mpfr_sinh (xq, x->value.complex.i, GFC_RND_MODE);
 767       mpfr_mul (xp, xp, xq, GFC_RND_MODE);
 768       mpfr_neg (result->value.complex.i, xp, GFC_RND_MODE );
 769
 770       mpfr_clear (xp);
 771       mpfr_clear (xq);
 772       break;
 773     default:
 774       gfc_internal_error ("in gfc_simplify_cos(): Bad type");
 775     }
 776
 777   return range_check (result, "COS");
 778
 779 }
 780
 781
 782 gfc_expr *
 783 gfc_simplify_cosh (gfc_expr * x)
 784 {
 785   gfc_expr *result;
 786
 787   if (x->expr_type != EXPR_CONSTANT)
 788     return NULL;
 789
 790   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 791
 792   mpfr_cosh (result->value.real, x->value.real, GFC_RND_MODE);
 793
 794   return range_check (result, "COSH");
 795 }
 796
 797
 798 gfc_expr *
 799 gfc_simplify_dcmplx (gfc_expr * x, gfc_expr * y)
 800 {
 801
 802   if (x->expr_type != EXPR_CONSTANT
 803       || (y != NULL && y->expr_type != EXPR_CONSTANT))
 804     return NULL;
 805
 806   return simplify_cmplx ("DCMPLX", x, y, gfc_default_double_kind);
 807 }
 808
 809
 810 gfc_expr *
 811 gfc_simplify_dble (gfc_expr * e)
 812 {
 813   gfc_expr *result;
 814
 815   if (e->expr_type != EXPR_CONSTANT)
 816     return NULL;
 817
 818   switch (e->ts.type)
 819     {
 820     case BT_INTEGER:
 821       result = gfc_int2real (e, gfc_default_double_kind);
 822       break;
 823
 824     case BT_REAL:
 825       result = gfc_real2real (e, gfc_default_double_kind);
 826       break;
 827
 828     case BT_COMPLEX:
 829       result = gfc_complex2real (e, gfc_default_double_kind);
 830       break;
 831
 832     default:
 833       gfc_internal_error ("gfc_simplify_dble(): bad type at %L", &e->where);
 834     }
 835
 836   return range_check (result, "DBLE");
 837 }
 838
 839
 840 gfc_expr *
 841 gfc_simplify_digits (gfc_expr * x)
 842 {
 843   int i, digits;
 844
 845   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
 846   switch (x->ts.type)
 847     {
 848     case BT_INTEGER:
 849       digits = gfc_integer_kinds[i].digits;
 850       break;
 851
 852     case BT_REAL:
 853     case BT_COMPLEX:
 854       digits = gfc_real_kinds[i].digits;
 855       break;
 856
 857     default:
 858       gcc_unreachable ();
 859     }
 860
 861   return gfc_int_expr (digits);
 862 }
 863
 864
 865 gfc_expr *
 866 gfc_simplify_dim (gfc_expr * x, gfc_expr * y)
 867 {
 868   gfc_expr *result;
 869
 870   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
 871     return NULL;
 872
 873   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 874
 875   switch (x->ts.type)
 876     {
 877     case BT_INTEGER:
 878       if (mpz_cmp (x->value.integer, y->value.integer) > 0)
 879         mpz_sub (result->value.integer, x->value.integer, y->value.integer);
 880       else
 881         mpz_set_ui (result->value.integer, 0);
 882
 883       break;
 884
 885     case BT_REAL:
 886       if (mpfr_cmp (x->value.real, y->value.real) > 0)
 887         mpfr_sub (result->value.real, x->value.real, y->value.real, GFC_RND_MODE);
 888       else
 889         mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
 890
 891       break;
 892
 893     default:
 894       gfc_internal_error ("gfc_simplify_dim(): Bad type");
 895     }
 896
 897   return range_check (result, "DIM");
 898 }
 899
 900
 901 gfc_expr *
 902 gfc_simplify_dprod (gfc_expr * x, gfc_expr * y)
 903 {
 904   gfc_expr *a1, *a2, *result;
 905
 906   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
 907     return NULL;
 908
 909   result =
 910     gfc_constant_result (BT_REAL, gfc_default_double_kind, &x->where);
 911
 912   a1 = gfc_real2real (x, gfc_default_double_kind);
 913   a2 = gfc_real2real (y, gfc_default_double_kind);
 914
 915   mpfr_mul (result->value.real, a1->value.real, a2->value.real, GFC_RND_MODE);
 916
 917   gfc_free_expr (a1);
 918   gfc_free_expr (a2);
 919
 920   return range_check (result, "DPROD");
 921 }
 922
 923
 924 gfc_expr *
 925 gfc_simplify_epsilon (gfc_expr * e)
 926 {
 927   gfc_expr *result;
 928   int i;
 929
 930   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
 931
 932   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
 933
 934   mpfr_set (result->value.real, gfc_real_kinds[i].epsilon, GFC_RND_MODE);
 935
 936   return range_check (result, "EPSILON");
 937 }
 938
 939
 940 gfc_expr *
 941 gfc_simplify_exp (gfc_expr * x)
 942 {
 943   gfc_expr *result;
 944   mpfr_t xp, xq;
 945
 946   if (x->expr_type != EXPR_CONSTANT)
 947     return NULL;
 948
 949   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
 950
 951   switch (x->ts.type)
 952     {
 953     case BT_REAL:
 954       mpfr_exp(result->value.real, x->value.real, GFC_RND_MODE);
 955       break;
 956
 957     case BT_COMPLEX:
 958       gfc_set_model_kind (x->ts.kind);
 959       mpfr_init (xp);
 960       mpfr_init (xq);
 961       mpfr_exp (xq, x->value.complex.r, GFC_RND_MODE);
 962       mpfr_cos (xp, x->value.complex.i, GFC_RND_MODE);
 963       mpfr_mul (result->value.complex.r, xq, xp, GFC_RND_MODE);
 964       mpfr_sin (xp, x->value.complex.i, GFC_RND_MODE);
 965       mpfr_mul (result->value.complex.i, xq, xp, GFC_RND_MODE);
 966       mpfr_clear (xp);
 967       mpfr_clear (xq);
 968       break;
 969
 970     default:
 971       gfc_internal_error ("in gfc_simplify_exp(): Bad type");
 972     }
 973
 974   return range_check (result, "EXP");
 975 }
 976
 977 /* FIXME:  MPFR should be able to do this better */
 978 gfc_expr *
 979 gfc_simplify_exponent (gfc_expr * x)
 980 {
 981   int i;
 982   mpfr_t tmp;
 983   gfc_expr *result;
 984
 985   if (x->expr_type != EXPR_CONSTANT)
 986     return NULL;
 987
 988   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
 989                                 &x->where);
 990
 991   gfc_set_model (x->value.real);
 992
 993   if (mpfr_sgn (x->value.real) == 0)
 994     {
 995       mpz_set_ui (result->value.integer, 0);
 996       return result;
 997     }
 998
 999   mpfr_init (tmp);
1000
1001   mpfr_abs (tmp, x->value.real, GFC_RND_MODE);
1002   mpfr_log2 (tmp, tmp, GFC_RND_MODE);
1003
1004   gfc_mpfr_to_mpz (result->value.integer, tmp);
1005
1006   /* The model number for tiny(x) is b**(emin - 1) where b is the base and emin
1007      is the smallest exponent value.  So, we need to add 1 if x is tiny(x).  */
1008   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
1009   if (mpfr_cmp (x->value.real, gfc_real_kinds[i].tiny) == 0)
1010     mpz_add_ui (result->value.integer,result->value.integer, 1);
1011
1012   mpfr_clear (tmp);
1013
1014   return range_check (result, "EXPONENT");
1015 }
1016
1017
1018 gfc_expr *
1019 gfc_simplify_float (gfc_expr * a)
1020 {
1021   gfc_expr *result;
1022
1023   if (a->expr_type != EXPR_CONSTANT)
1024     return NULL;
1025
1026   result = gfc_int2real (a, gfc_default_real_kind);
1027   return range_check (result, "FLOAT");
1028 }
1029
1030
1031 gfc_expr *
1032 gfc_simplify_floor (gfc_expr * e, gfc_expr * k)
1033 {
1034   gfc_expr *result;
1035   mpfr_t floor;
1036   int kind;
1037
1038   kind = get_kind (BT_INTEGER, k, "FLOOR", gfc_default_integer_kind);
1039   if (kind == -1)
1040     gfc_internal_error ("gfc_simplify_floor(): Bad kind");
1041
1042   if (e->expr_type != EXPR_CONSTANT)
1043     return NULL;
1044
1045   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
1046
1047   gfc_set_model_kind (kind);
1048   mpfr_init (floor);
1049   mpfr_floor (floor, e->value.real);
1050
1051   gfc_mpfr_to_mpz (result->value.integer, floor);
1052
1053   mpfr_clear (floor);
1054
1055   return range_check (result, "FLOOR");
1056 }
1057
1058
1059 gfc_expr *
1060 gfc_simplify_fraction (gfc_expr * x)
1061 {
1062   gfc_expr *result;
1063   mpfr_t absv, exp, pow2;
1064
1065   if (x->expr_type != EXPR_CONSTANT)
1066     return NULL;
1067
1068   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
1069
1070   gfc_set_model_kind (x->ts.kind);
1071
1072   if (mpfr_sgn (x->value.real) == 0)
1073     {
1074       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
1075       return result;
1076     }
1077
1078   mpfr_init (exp);
1079   mpfr_init (absv);
1080   mpfr_init (pow2);
1081
1082   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
1083   mpfr_log2 (exp, absv, GFC_RND_MODE);
1084
1085   mpfr_trunc (exp, exp);
1086   mpfr_add_ui (exp, exp, 1, GFC_RND_MODE);
1087
1088   mpfr_ui_pow (pow2, 2, exp, GFC_RND_MODE);
1089
1090   mpfr_div (result->value.real, absv, pow2, GFC_RND_MODE);
1091
1092   mpfr_clear (exp);
1093   mpfr_clear (absv);
1094   mpfr_clear (pow2);
1095
1096   return range_check (result, "FRACTION");
1097 }
1098
1099
1100 gfc_expr *
1101 gfc_simplify_huge (gfc_expr * e)
1102 {
1103   gfc_expr *result;
1104   int i;
1105
1106   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
1107
1108   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
1109
1110   switch (e->ts.type)
1111     {
1112     case BT_INTEGER:
1113       mpz_set (result->value.integer, gfc_integer_kinds[i].huge);
1114       break;
1115
1116     case BT_REAL:
1117       mpfr_set (result->value.real, gfc_real_kinds[i].huge, GFC_RND_MODE);
1118       break;
1119
1120     default:
1121       gcc_unreachable ();
1122     }
1123
1124   return result;
1125 }
1126
1127
1128 gfc_expr *
1129 gfc_simplify_iachar (gfc_expr * e)
1130 {
1131   gfc_expr *result;
1132   int index;
1133
1134   if (e->expr_type != EXPR_CONSTANT)
1135     return NULL;
1136
1137   if (e->value.character.length != 1)
1138     {
1139       gfc_error ("Argument of IACHAR at %L must be of length one", &e->where);
1140       return &gfc_bad_expr;
1141     }
1142
1143   index = xascii_table[(int) e->value.character.string[0] & 0xFF];
1144
1145   result = gfc_int_expr (index);
1146   result->where = e->where;
1147
1148   return range_check (result, "IACHAR");
1149 }
1150
1151
1152 gfc_expr *
1153 gfc_simplify_iand (gfc_expr * x, gfc_expr * y)
1154 {
1155   gfc_expr *result;
1156
1157   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1158     return NULL;
1159
1160   result = gfc_constant_result (BT_INTEGER, x->ts.kind, &x->where);
1161
1162   mpz_and (result->value.integer, x->value.integer, y->value.integer);
1163
1164   return range_check (result, "IAND");
1165 }
1166
1167
1168 gfc_expr *
1169 gfc_simplify_ibclr (gfc_expr * x, gfc_expr * y)
1170 {
1171   gfc_expr *result;
1172   int k, pos;
1173
1174   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1175     return NULL;
1176
1177   if (gfc_extract_int (y, &pos) != NULL || pos < 0)
1178     {
1179       gfc_error ("Invalid second argument of IBCLR at %L", &y->where);
1180       return &gfc_bad_expr;
1181     }
1182
1183   k = gfc_validate_kind (x->ts.type, x->ts.kind, false);
1184
1185   if (pos > gfc_integer_kinds[k].bit_size)
1186     {
1187       gfc_error ("Second argument of IBCLR exceeds bit size at %L",
1188                  &y->where);
1189       return &gfc_bad_expr;
1190     }
1191
1192   result = gfc_copy_expr (x);
1193
1194   mpz_clrbit (result->value.integer, pos);
1195   return range_check (result, "IBCLR");
1196 }
1197
1198
1199 gfc_expr *
1200 gfc_simplify_ibits (gfc_expr * x, gfc_expr * y, gfc_expr * z)
1201 {
1202   gfc_expr *result;
1203   int pos, len;
1204   int i, k, bitsize;
1205   int *bits;
1206
1207   if (x->expr_type != EXPR_CONSTANT
1208       || y->expr_type != EXPR_CONSTANT
1209       || z->expr_type != EXPR_CONSTANT)
1210     return NULL;
1211
1212   if (gfc_extract_int (y, &pos) != NULL || pos < 0)
1213     {
1214       gfc_error ("Invalid second argument of IBITS at %L", &y->where);
1215       return &gfc_bad_expr;
1216     }
1217
1218   if (gfc_extract_int (z, &len) != NULL || len < 0)
1219     {
1220       gfc_error ("Invalid third argument of IBITS at %L", &z->where);
1221       return &gfc_bad_expr;
1222     }
1223
1224   k = gfc_validate_kind (BT_INTEGER, x->ts.kind, false);
1225
1226   bitsize = gfc_integer_kinds[k].bit_size;
1227
1228   if (pos + len > bitsize)
1229     {
1230       gfc_error
1231         ("Sum of second and third arguments of IBITS exceeds bit size "
1232          "at %L", &y->where);
1233       return &gfc_bad_expr;
1234     }
1235
1236   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1237
1238   bits = gfc_getmem (bitsize * sizeof (int));
1239
1240   for (i = 0; i < bitsize; i++)
1241     bits[i] = 0;
1242
1243   for (i = 0; i < len; i++)
1244     bits[i] = mpz_tstbit (x->value.integer, i + pos);
1245
1246   for (i = 0; i < bitsize; i++)
1247     {
1248       if (bits[i] == 0)
1249         {
1250           mpz_clrbit (result->value.integer, i);
1251         }
1252       else if (bits[i] == 1)
1253         {
1254           mpz_setbit (result->value.integer, i);
1255         }
1256       else
1257         {
1258           gfc_internal_error ("IBITS: Bad bit");
1259         }
1260     }
1261
1262   gfc_free (bits);
1263
1264   return range_check (result, "IBITS");
1265 }
1266
1267
1268 gfc_expr *
1269 gfc_simplify_ibset (gfc_expr * x, gfc_expr * y)
1270 {
1271   gfc_expr *result;
1272   int k, pos;
1273
1274   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1275     return NULL;
1276
1277   if (gfc_extract_int (y, &pos) != NULL || pos < 0)
1278     {
1279       gfc_error ("Invalid second argument of IBSET at %L", &y->where);
1280       return &gfc_bad_expr;
1281     }
1282
1283   k = gfc_validate_kind (x->ts.type, x->ts.kind, false);
1284
1285   if (pos > gfc_integer_kinds[k].bit_size)
1286     {
1287       gfc_error ("Second argument of IBSET exceeds bit size at %L",
1288                  &y->where);
1289       return &gfc_bad_expr;
1290     }
1291
1292   result = gfc_copy_expr (x);
1293
1294   mpz_setbit (result->value.integer, pos);
1295   return range_check (result, "IBSET");
1296 }
1297
1298
1299 gfc_expr *
1300 gfc_simplify_ichar (gfc_expr * e)
1301 {
1302   gfc_expr *result;
1303   int index;
1304
1305   if (e->expr_type != EXPR_CONSTANT)
1306     return NULL;
1307
1308   if (e->value.character.length != 1)
1309     {
1310       gfc_error ("Argument of ICHAR at %L must be of length one", &e->where);
1311       return &gfc_bad_expr;
1312     }
1313
1314   index = (int) e->value.character.string[0];
1315
1316   if (index < CHAR_MIN || index > CHAR_MAX)
1317     {
1318       gfc_error ("Argument of ICHAR at %L out of range of this processor",
1319                  &e->where);
1320       return &gfc_bad_expr;
1321     }
1322
1323   result = gfc_int_expr (index);
1324   result->where = e->where;
1325   return range_check (result, "ICHAR");
1326 }
1327
1328
1329 gfc_expr *
1330 gfc_simplify_ieor (gfc_expr * x, gfc_expr * y)
1331 {
1332   gfc_expr *result;
1333
1334   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1335     return NULL;
1336
1337   result = gfc_constant_result (BT_INTEGER, x->ts.kind, &x->where);
1338
1339   mpz_xor (result->value.integer, x->value.integer, y->value.integer);
1340
1341   return range_check (result, "IEOR");
1342 }
1343
1344
1345 gfc_expr *
1346 gfc_simplify_index (gfc_expr * x, gfc_expr * y, gfc_expr * b)
1347 {
1348   gfc_expr *result;
1349   int back, len, lensub;
1350   int i, j, k, count, index = 0, start;
1351
1352   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1353     return NULL;
1354
1355   if (b != NULL && b->value.logical != 0)
1356     back = 1;
1357   else
1358     back = 0;
1359
1360   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1361                                 &x->where);
1362
1363   len = x->value.character.length;
1364   lensub = y->value.character.length;
1365
1366   if (len < lensub)
1367     {
1368       mpz_set_si (result->value.integer, 0);
1369       return result;
1370     }
1371
1372   if (back == 0)
1373     {
1374
1375       if (lensub == 0)
1376         {
1377           mpz_set_si (result->value.integer, 1);
1378           return result;
1379         }
1380       else if (lensub == 1)
1381         {
1382           for (i = 0; i < len; i++)
1383             {
1384               for (j = 0; j < lensub; j++)
1385                 {
1386                   if (y->value.character.string[j] ==
1387                       x->value.character.string[i])
1388                     {
1389                       index = i + 1;
1390                       goto done;
1391                     }
1392                 }
1393             }
1394         }
1395       else
1396         {
1397           for (i = 0; i < len; i++)
1398             {
1399               for (j = 0; j < lensub; j++)
1400                 {
1401                   if (y->value.character.string[j] ==
1402                       x->value.character.string[i])
1403                     {
1404                       start = i;
1405                       count = 0;
1406
1407                       for (k = 0; k < lensub; k++)
1408                         {
1409                           if (y->value.character.string[k] ==
1410                               x->value.character.string[k + start])
1411                             count++;
1412                         }
1413
1414                       if (count == lensub)
1415                         {
1416                           index = start + 1;
1417                           goto done;
1418                         }
1419                     }
1420                 }
1421             }
1422         }
1423
1424     }
1425   else
1426     {
1427
1428       if (lensub == 0)
1429         {
1430           mpz_set_si (result->value.integer, len + 1);
1431           return result;
1432         }
1433       else if (lensub == 1)
1434         {
1435           for (i = 0; i < len; i++)
1436             {
1437               for (j = 0; j < lensub; j++)
1438                 {
1439                   if (y->value.character.string[j] ==
1440                       x->value.character.string[len - i])
1441                     {
1442                       index = len - i + 1;
1443                       goto done;
1444                     }
1445                 }
1446             }
1447         }
1448       else
1449         {
1450           for (i = 0; i < len; i++)
1451             {
1452               for (j = 0; j < lensub; j++)
1453                 {
1454                   if (y->value.character.string[j] ==
1455                       x->value.character.string[len - i])
1456                     {
1457                       start = len - i;
1458                       if (start <= len - lensub)
1459                         {
1460                           count = 0;
1461                           for (k = 0; k < lensub; k++)
1462                             if (y->value.character.string[k] ==
1463                                 x->value.character.string[k + start])
1464                               count++;
1465
1466                           if (count == lensub)
1467                             {
1468                               index = start + 1;
1469                               goto done;
1470                             }
1471                         }
1472                       else
1473                         {
1474                           continue;
1475                         }
1476                     }
1477                 }
1478             }
1479         }
1480     }
1481
1482 done:
1483   mpz_set_si (result->value.integer, index);
1484   return range_check (result, "INDEX");
1485 }
1486
1487
1488 gfc_expr *
1489 gfc_simplify_int (gfc_expr * e, gfc_expr * k)
1490 {
1491   gfc_expr *rpart, *rtrunc, *result;
1492   int kind;
1493
1494   kind = get_kind (BT_INTEGER, k, "INT", gfc_default_integer_kind);
1495   if (kind == -1)
1496     return &gfc_bad_expr;
1497
1498   if (e->expr_type != EXPR_CONSTANT)
1499     return NULL;
1500
1501   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
1502
1503   switch (e->ts.type)
1504     {
1505     case BT_INTEGER:
1506       mpz_set (result->value.integer, e->value.integer);
1507       break;
1508
1509     case BT_REAL:
1510       rtrunc = gfc_copy_expr (e);
1511       mpfr_trunc (rtrunc->value.real, e->value.real);
1512       gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real);
1513       gfc_free_expr (rtrunc);
1514       break;
1515
1516     case BT_COMPLEX:
1517       rpart = gfc_complex2real (e, kind);
1518       rtrunc = gfc_copy_expr (rpart);
1519       mpfr_trunc (rtrunc->value.real, rpart->value.real);
1520       gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real);
1521       gfc_free_expr (rpart);
1522       gfc_free_expr (rtrunc);
1523       break;
1524
1525     default:
1526       gfc_error ("Argument of INT at %L is not a valid type", &e->where);
1527       gfc_free_expr (result);
1528       return &gfc_bad_expr;
1529     }
1530
1531   return range_check (result, "INT");
1532 }
1533
1534
1535 gfc_expr *
1536 gfc_simplify_ifix (gfc_expr * e)
1537 {
1538   gfc_expr *rtrunc, *result;
1539
1540   if (e->expr_type != EXPR_CONSTANT)
1541     return NULL;
1542
1543   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1544                                 &e->where);
1545
1546   rtrunc = gfc_copy_expr (e);
1547
1548   mpfr_trunc (rtrunc->value.real, e->value.real);
1549   gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real);
1550
1551   gfc_free_expr (rtrunc);
1552   return range_check (result, "IFIX");
1553 }
1554
1555
1556 gfc_expr *
1557 gfc_simplify_idint (gfc_expr * e)
1558 {
1559   gfc_expr *rtrunc, *result;
1560
1561   if (e->expr_type != EXPR_CONSTANT)
1562     return NULL;
1563
1564   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1565                                 &e->where);
1566
1567   rtrunc = gfc_copy_expr (e);
1568
1569   mpfr_trunc (rtrunc->value.real, e->value.real);
1570   gfc_mpfr_to_mpz (result->value.integer, rtrunc->value.real);
1571
1572   gfc_free_expr (rtrunc);
1573   return range_check (result, "IDINT");
1574 }
1575
1576
1577 gfc_expr *
1578 gfc_simplify_ior (gfc_expr * x, gfc_expr * y)
1579 {
1580   gfc_expr *result;
1581
1582   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
1583     return NULL;
1584
1585   result = gfc_constant_result (BT_INTEGER, x->ts.kind, &x->where);
1586
1587   mpz_ior (result->value.integer, x->value.integer, y->value.integer);
1588   return range_check (result, "IOR");
1589 }
1590
1591
1592 gfc_expr *
1593 gfc_simplify_ishft (gfc_expr * e, gfc_expr * s)
1594 {
1595   gfc_expr *result;
1596   int shift, ashift, isize, k, *bits, i;
1597
1598   if (e->expr_type != EXPR_CONSTANT || s->expr_type != EXPR_CONSTANT)
1599     return NULL;
1600
1601   if (gfc_extract_int (s, &shift) != NULL)
1602     {
1603       gfc_error ("Invalid second argument of ISHFT at %L", &s->where);
1604       return &gfc_bad_expr;
1605     }
1606
1607   k = gfc_validate_kind (BT_INTEGER, e->ts.kind, false);
1608
1609   isize = gfc_integer_kinds[k].bit_size;
1610
1611   if (shift >= 0)
1612     ashift = shift;
1613   else
1614     ashift = -shift;
1615
1616   if (ashift > isize)
1617     {
1618       gfc_error
1619         ("Magnitude of second argument of ISHFT exceeds bit size at %L",
1620          &s->where);
1621       return &gfc_bad_expr;
1622     }
1623
1624   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
1625
1626   if (shift == 0)
1627     {
1628       mpz_set (result->value.integer, e->value.integer);
1629       return range_check (result, "ISHFT");
1630     }
1631
1632   bits = gfc_getmem (isize * sizeof (int));
1633
1634   for (i = 0; i < isize; i++)
1635     bits[i] = mpz_tstbit (e->value.integer, i);
1636
1637   if (shift > 0)
1638     {
1639       for (i = 0; i < shift; i++)
1640         mpz_clrbit (result->value.integer, i);
1641
1642       for (i = 0; i < isize - shift; i++)
1643         {
1644           if (bits[i] == 0)
1645             mpz_clrbit (result->value.integer, i + shift);
1646           else
1647             mpz_setbit (result->value.integer, i + shift);
1648         }
1649     }
1650   else
1651     {
1652       for (i = isize - 1; i >= isize - ashift; i--)
1653         mpz_clrbit (result->value.integer, i);
1654
1655       for (i = isize - 1; i >= ashift; i--)
1656         {
1657           if (bits[i] == 0)
1658             mpz_clrbit (result->value.integer, i - ashift);
1659           else
1660             mpz_setbit (result->value.integer, i - ashift);
1661         }
1662     }
1663
1664   twos_complement (result->value.integer, isize);
1665
1666   gfc_free (bits);
1667   return result;
1668 }
1669
1670
1671 gfc_expr *
1672 gfc_simplify_ishftc (gfc_expr * e, gfc_expr * s, gfc_expr * sz)
1673 {
1674   gfc_expr *result;
1675   int shift, ashift, isize, delta, k;
1676   int i, *bits;
1677
1678   if (e->expr_type != EXPR_CONSTANT || s->expr_type != EXPR_CONSTANT)
1679     return NULL;
1680
1681   if (gfc_extract_int (s, &shift) != NULL)
1682     {
1683       gfc_error ("Invalid second argument of ISHFTC at %L", &s->where);
1684       return &gfc_bad_expr;
1685     }
1686
1687   k = gfc_validate_kind (e->ts.type, e->ts.kind, false);
1688
1689   if (sz != NULL)
1690     {
1691       if (gfc_extract_int (sz, &isize) != NULL || isize < 0)
1692         {
1693           gfc_error ("Invalid third argument of ISHFTC at %L", &sz->where);
1694           return &gfc_bad_expr;
1695         }
1696     }
1697   else
1698     isize = gfc_integer_kinds[k].bit_size;
1699
1700   if (shift >= 0)
1701     ashift = shift;
1702   else
1703     ashift = -shift;
1704
1705   if (ashift > isize)
1706     {
1707       gfc_error
1708         ("Magnitude of second argument of ISHFTC exceeds third argument "
1709          "at %L", &s->where);
1710       return &gfc_bad_expr;
1711     }
1712
1713   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
1714
1715   if (shift == 0)
1716     {
1717       mpz_set (result->value.integer, e->value.integer);
1718       return result;
1719     }
1720
1721   bits = gfc_getmem (isize * sizeof (int));
1722
1723   for (i = 0; i < isize; i++)
1724     bits[i] = mpz_tstbit (e->value.integer, i);
1725
1726   delta = isize - ashift;
1727
1728   if (shift > 0)
1729     {
1730       for (i = 0; i < delta; i++)
1731         {
1732           if (bits[i] == 0)
1733             mpz_clrbit (result->value.integer, i + shift);
1734           else
1735             mpz_setbit (result->value.integer, i + shift);
1736         }
1737
1738       for (i = delta; i < isize; i++)
1739         {
1740           if (bits[i] == 0)
1741             mpz_clrbit (result->value.integer, i - delta);
1742           else
1743             mpz_setbit (result->value.integer, i - delta);
1744         }
1745     }
1746   else
1747     {
1748       for (i = 0; i < ashift; i++)
1749         {
1750           if (bits[i] == 0)
1751             mpz_clrbit (result->value.integer, i + delta);
1752           else
1753             mpz_setbit (result->value.integer, i + delta);
1754         }
1755
1756       for (i = ashift; i < isize; i++)
1757         {
1758           if (bits[i] == 0)
1759             mpz_clrbit (result->value.integer, i + shift);
1760           else
1761             mpz_setbit (result->value.integer, i + shift);
1762         }
1763     }
1764
1765   twos_complement (result->value.integer, isize);
1766
1767   gfc_free (bits);
1768   return result;
1769 }
1770
1771
1772 gfc_expr *
1773 gfc_simplify_kind (gfc_expr * e)
1774 {
1775
1776   if (e->ts.type == BT_DERIVED)
1777     {
1778       gfc_error ("Argument of KIND at %L is a DERIVED type", &e->where);
1779       return &gfc_bad_expr;
1780     }
1781
1782   return gfc_int_expr (e->ts.kind);
1783 }
1784
1785
1786 static gfc_expr *
1787 simplify_bound (gfc_expr * array, gfc_expr * dim, int upper)
1788 {
1789   gfc_ref *ref;
1790   gfc_array_spec *as;
1791   gfc_expr *e;
1792   int d;
1793
1794   if (array->expr_type != EXPR_VARIABLE)
1795     return NULL;
1796
1797   if (dim == NULL)
1798     /* TODO: Simplify constant multi-dimensional bounds.  */
1799     return NULL;
1800
1801   if (dim->expr_type != EXPR_CONSTANT)
1802     return NULL;
1803
1804   /* Follow any component references.  */
1805   as = array->symtree->n.sym->as;
1806   for (ref = array->ref; ref; ref = ref->next)
1807     {
1808       switch (ref->type)
1809         {
1810         case REF_ARRAY:
1811           switch (ref->u.ar.type)
1812             {
1813             case AR_ELEMENT:
1814               as = NULL;
1815               continue;
1816
1817             case AR_FULL:
1818               /* We're done because 'as' has already been set in the
1819                  previous iteration.  */
1820               goto done;
1821
1822             case AR_SECTION:
1823             case AR_UNKNOWN:
1824               return NULL;
1825             }
1826
1827           gcc_unreachable ();
1828
1829         case REF_COMPONENT:
1830           as = ref->u.c.component->as;
1831           continue;
1832
1833         case REF_SUBSTRING:
1834           continue;
1835         }
1836     }
1837
1838   gcc_unreachable ();
1839
1840  done:
1841   if (as->type == AS_DEFERRED || as->type == AS_ASSUMED_SHAPE)
1842     return NULL;
1843
1844   d = mpz_get_si (dim->value.integer);
1845
1846   if (d < 1 || d > as->rank
1847       || (d == as->rank && as->type == AS_ASSUMED_SIZE && upper))
1848     {
1849       gfc_error ("DIM argument at %L is out of bounds", &dim->where);
1850       return &gfc_bad_expr;
1851     }
1852
1853   e = upper ? as->upper[d-1] : as->lower[d-1];
1854
1855   if (e->expr_type != EXPR_CONSTANT)
1856     return NULL;
1857
1858   return gfc_copy_expr (e);
1859 }
1860
1861
1862 gfc_expr *
1863 gfc_simplify_lbound (gfc_expr * array, gfc_expr * dim)
1864 {
1865   return simplify_bound (array, dim, 0);
1866 }
1867
1868
1869 gfc_expr *
1870 gfc_simplify_len (gfc_expr * e)
1871 {
1872   gfc_expr *result;
1873
1874   if (e->expr_type != EXPR_CONSTANT)
1875     return NULL;
1876
1877   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1878                                 &e->where);
1879
1880   mpz_set_si (result->value.integer, e->value.character.length);
1881   return range_check (result, "LEN");
1882 }
1883
1884
1885 gfc_expr *
1886 gfc_simplify_len_trim (gfc_expr * e)
1887 {
1888   gfc_expr *result;
1889   int count, len, lentrim, i;
1890
1891   if (e->expr_type != EXPR_CONSTANT)
1892     return NULL;
1893
1894   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
1895                                 &e->where);
1896
1897   len = e->value.character.length;
1898
1899   for (count = 0, i = 1; i <= len; i++)
1900     if (e->value.character.string[len - i] == ' ')
1901       count++;
1902     else
1903       break;
1904
1905   lentrim = len - count;
1906
1907   mpz_set_si (result->value.integer, lentrim);
1908   return range_check (result, "LEN_TRIM");
1909 }
1910
1911
1912 gfc_expr *
1913 gfc_simplify_lge (gfc_expr * a, gfc_expr * b)
1914 {
1915
1916   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
1917     return NULL;
1918
1919   return gfc_logical_expr (gfc_compare_string (a, b, xascii_table) >= 0,
1920                            &a->where);
1921 }
1922
1923
1924 gfc_expr *
1925 gfc_simplify_lgt (gfc_expr * a, gfc_expr * b)
1926 {
1927
1928   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
1929     return NULL;
1930
1931   return gfc_logical_expr (gfc_compare_string (a, b, xascii_table) > 0,
1932                            &a->where);
1933 }
1934
1935
1936 gfc_expr *
1937 gfc_simplify_lle (gfc_expr * a, gfc_expr * b)
1938 {
1939
1940   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
1941     return NULL;
1942
1943   return gfc_logical_expr (gfc_compare_string (a, b, xascii_table) <= 0,
1944                            &a->where);
1945 }
1946
1947
1948 gfc_expr *
1949 gfc_simplify_llt (gfc_expr * a, gfc_expr * b)
1950 {
1951
1952   if (a->expr_type != EXPR_CONSTANT || b->expr_type != EXPR_CONSTANT)
1953     return NULL;
1954
1955   return gfc_logical_expr (gfc_compare_string (a, b, xascii_table) < 0,
1956                            &a->where);
1957 }
1958
1959
1960 gfc_expr *
1961 gfc_simplify_log (gfc_expr * x)
1962 {
1963   gfc_expr *result;
1964   mpfr_t xr, xi;
1965
1966   if (x->expr_type != EXPR_CONSTANT)
1967     return NULL;
1968
1969   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
1970
1971   gfc_set_model_kind (x->ts.kind);
1972
1973   switch (x->ts.type)
1974     {
1975     case BT_REAL:
1976       if (mpfr_sgn (x->value.real) <= 0)
1977         {
1978           gfc_error
1979             ("Argument of LOG at %L cannot be less than or equal to zero",
1980              &x->where);
1981           gfc_free_expr (result);
1982           return &gfc_bad_expr;
1983         }
1984
1985       mpfr_log(result->value.real, x->value.real, GFC_RND_MODE);
1986       break;
1987
1988     case BT_COMPLEX:
1989       if ((mpfr_sgn (x->value.complex.r) == 0)
1990           && (mpfr_sgn (x->value.complex.i) == 0))
1991         {
1992           gfc_error ("Complex argument of LOG at %L cannot be zero",
1993                      &x->where);
1994           gfc_free_expr (result);
1995           return &gfc_bad_expr;
1996         }
1997
1998       mpfr_init (xr);
1999       mpfr_init (xi);
2000
2001       arctangent2 (x->value.complex.i, x->value.complex.r,
2002                    result->value.complex.i);
2003
2004       mpfr_mul (xr, x->value.complex.r, x->value.complex.r, GFC_RND_MODE);
2005       mpfr_mul (xi, x->value.complex.i, x->value.complex.i, GFC_RND_MODE);
2006       mpfr_add (xr, xr, xi, GFC_RND_MODE);
2007       mpfr_sqrt (xr, xr, GFC_RND_MODE);
2008       mpfr_log (result->value.complex.r, xr, GFC_RND_MODE);
2009
2010       mpfr_clear (xr);
2011       mpfr_clear (xi);
2012
2013       break;
2014
2015     default:
2016       gfc_internal_error ("gfc_simplify_log: bad type");
2017     }
2018
2019   return range_check (result, "LOG");
2020 }
2021
2022
2023 gfc_expr *
2024 gfc_simplify_log10 (gfc_expr * x)
2025 {
2026   gfc_expr *result;
2027
2028   if (x->expr_type != EXPR_CONSTANT)
2029     return NULL;
2030
2031   gfc_set_model_kind (x->ts.kind);
2032
2033   if (mpfr_sgn (x->value.real) <= 0)
2034     {
2035       gfc_error
2036         ("Argument of LOG10 at %L cannot be less than or equal to zero",
2037          &x->where);
2038       return &gfc_bad_expr;
2039     }
2040
2041   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
2042
2043   mpfr_log10 (result->value.real, x->value.real, GFC_RND_MODE);
2044
2045   return range_check (result, "LOG10");
2046 }
2047
2048
2049 gfc_expr *
2050 gfc_simplify_logical (gfc_expr * e, gfc_expr * k)
2051 {
2052   gfc_expr *result;
2053   int kind;
2054
2055   kind = get_kind (BT_LOGICAL, k, "LOGICAL", gfc_default_logical_kind);
2056   if (kind < 0)
2057     return &gfc_bad_expr;
2058
2059   if (e->expr_type != EXPR_CONSTANT)
2060     return NULL;
2061
2062   result = gfc_constant_result (BT_LOGICAL, kind, &e->where);
2063
2064   result->value.logical = e->value.logical;
2065
2066   return result;
2067 }
2068
2069
2070 /* This function is special since MAX() can take any number of
2071    arguments.  The simplified expression is a rewritten version of the
2072    argument list containing at most one constant element.  Other
2073    constant elements are deleted.  Because the argument list has
2074    already been checked, this function always succeeds.  sign is 1 for
2075    MAX(), -1 for MIN().  */
2076
2077 static gfc_expr *
2078 simplify_min_max (gfc_expr * expr, int sign)
2079 {
2080   gfc_actual_arglist *arg, *last, *extremum;
2081   gfc_intrinsic_sym * specific;
2082
2083   last = NULL;
2084   extremum = NULL;
2085   specific = expr->value.function.isym;
2086
2087   arg = expr->value.function.actual;
2088
2089   for (; arg; last = arg, arg = arg->next)
2090     {
2091       if (arg->expr->expr_type != EXPR_CONSTANT)
2092         continue;
2093
2094       if (extremum == NULL)
2095         {
2096           extremum = arg;
2097           continue;
2098         }
2099
2100       switch (arg->expr->ts.type)
2101         {
2102         case BT_INTEGER:
2103           if (mpz_cmp (arg->expr->value.integer,
2104                        extremum->expr->value.integer) * sign > 0)
2105             mpz_set (extremum->expr->value.integer, arg->expr->value.integer);
2106
2107           break;
2108
2109         case BT_REAL:
2110           if (mpfr_cmp (arg->expr->value.real, extremum->expr->value.real) *
2111               sign > 0)
2112             mpfr_set (extremum->expr->value.real, arg->expr->value.real,
2113                       GFC_RND_MODE);
2114
2115           break;
2116
2117         default:
2118           gfc_internal_error ("gfc_simplify_max(): Bad type in arglist");
2119         }
2120
2121       /* Delete the extra constant argument.  */
2122       if (last == NULL)
2123         expr->value.function.actual = arg->next;
2124       else
2125         last->next = arg->next;
2126
2127       arg->next = NULL;
2128       gfc_free_actual_arglist (arg);
2129       arg = last;
2130     }
2131
2132   /* If there is one value left, replace the function call with the
2133      expression.  */
2134   if (expr->value.function.actual->next != NULL)
2135     return NULL;
2136
2137   /* Convert to the correct type and kind.  */
2138   if (expr->ts.type != BT_UNKNOWN)
2139     return gfc_convert_constant (expr->value.function.actual->expr,
2140         expr->ts.type, expr->ts.kind);
2141
2142   if (specific->ts.type != BT_UNKNOWN)
2143     return gfc_convert_constant (expr->value.function.actual->expr,
2144         specific->ts.type, specific->ts.kind);
2145
2146   return gfc_copy_expr (expr->value.function.actual->expr);
2147 }
2148
2149
2150 gfc_expr *
2151 gfc_simplify_min (gfc_expr * e)
2152 {
2153   return simplify_min_max (e, -1);
2154 }
2155
2156
2157 gfc_expr *
2158 gfc_simplify_max (gfc_expr * e)
2159 {
2160   return simplify_min_max (e, 1);
2161 }
2162
2163
2164 gfc_expr *
2165 gfc_simplify_maxexponent (gfc_expr * x)
2166 {
2167   gfc_expr *result;
2168   int i;
2169
2170   i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
2171
2172   result = gfc_int_expr (gfc_real_kinds[i].max_exponent);
2173   result->where = x->where;
2174
2175   return result;
2176 }
2177
2178
2179 gfc_expr *
2180 gfc_simplify_minexponent (gfc_expr * x)
2181 {
2182   gfc_expr *result;
2183   int i;
2184
2185   i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
2186
2187   result = gfc_int_expr (gfc_real_kinds[i].min_exponent);
2188   result->where = x->where;
2189
2190   return result;
2191 }
2192
2193
2194 gfc_expr *
2195 gfc_simplify_mod (gfc_expr * a, gfc_expr * p)
2196 {
2197   gfc_expr *result;
2198   mpfr_t quot, iquot, term;
2199
2200   if (a->expr_type != EXPR_CONSTANT || p->expr_type != EXPR_CONSTANT)
2201     return NULL;
2202
2203   result = gfc_constant_result (a->ts.type, a->ts.kind, &a->where);
2204
2205   switch (a->ts.type)
2206     {
2207     case BT_INTEGER:
2208       if (mpz_cmp_ui (p->value.integer, 0) == 0)
2209         {
2210           /* Result is processor-dependent.  */
2211           gfc_error ("Second argument MOD at %L is zero", &a->where);
2212           gfc_free_expr (result);
2213           return &gfc_bad_expr;
2214         }
2215       mpz_tdiv_r (result->value.integer, a->value.integer, p->value.integer);
2216       break;
2217
2218     case BT_REAL:
2219       if (mpfr_cmp_ui (p->value.real, 0) == 0)
2220         {
2221           /* Result is processor-dependent.  */
2222           gfc_error ("Second argument of MOD at %L is zero", &p->where);
2223           gfc_free_expr (result);
2224           return &gfc_bad_expr;
2225         }
2226
2227       gfc_set_model_kind (a->ts.kind);
2228       mpfr_init (quot);
2229       mpfr_init (iquot);
2230       mpfr_init (term);
2231
2232       mpfr_div (quot, a->value.real, p->value.real, GFC_RND_MODE);
2233       mpfr_trunc (iquot, quot);
2234       mpfr_mul (term, iquot, p->value.real, GFC_RND_MODE);
2235       mpfr_sub (result->value.real, a->value.real, term, GFC_RND_MODE);
2236
2237       mpfr_clear (quot);
2238       mpfr_clear (iquot);
2239       mpfr_clear (term);
2240       break;
2241
2242     default:
2243       gfc_internal_error ("gfc_simplify_mod(): Bad arguments");
2244     }
2245
2246   return range_check (result, "MOD");
2247 }
2248
2249
2250 gfc_expr *
2251 gfc_simplify_modulo (gfc_expr * a, gfc_expr * p)
2252 {
2253   gfc_expr *result;
2254   mpfr_t quot, iquot, term;
2255
2256   if (a->expr_type != EXPR_CONSTANT || p->expr_type != EXPR_CONSTANT)
2257     return NULL;
2258
2259   result = gfc_constant_result (a->ts.type, a->ts.kind, &a->where);
2260
2261   switch (a->ts.type)
2262     {
2263     case BT_INTEGER:
2264       if (mpz_cmp_ui (p->value.integer, 0) == 0)
2265         {
2266           /* Result is processor-dependent. This processor just opts
2267              to not handle it at all.  */
2268           gfc_error ("Second argument of MODULO at %L is zero", &a->where);
2269           gfc_free_expr (result);
2270           return &gfc_bad_expr;
2271         }
2272       mpz_fdiv_r (result->value.integer, a->value.integer, p->value.integer);
2273
2274       break;
2275
2276     case BT_REAL:
2277       if (mpfr_cmp_ui (p->value.real, 0) == 0)
2278         {
2279           /* Result is processor-dependent.  */
2280           gfc_error ("Second argument of MODULO at %L is zero", &p->where);
2281           gfc_free_expr (result);
2282           return &gfc_bad_expr;
2283         }
2284
2285       gfc_set_model_kind (a->ts.kind);
2286       mpfr_init (quot);
2287       mpfr_init (iquot);
2288       mpfr_init (term);
2289
2290       mpfr_div (quot, a->value.real, p->value.real, GFC_RND_MODE);
2291       mpfr_floor (iquot, quot);
2292       mpfr_mul (term, iquot, p->value.real, GFC_RND_MODE);
2293       mpfr_sub (result->value.real, a->value.real, term, GFC_RND_MODE);
2294
2295       mpfr_clear (quot);
2296       mpfr_clear (iquot);
2297       mpfr_clear (term);
2298       break;
2299
2300     default:
2301       gfc_internal_error ("gfc_simplify_modulo(): Bad arguments");
2302     }
2303
2304   return range_check (result, "MODULO");
2305 }
2306
2307
2308 /* Exists for the sole purpose of consistency with other intrinsics.  */
2309 gfc_expr *
2310 gfc_simplify_mvbits (gfc_expr * f  ATTRIBUTE_UNUSED,
2311                      gfc_expr * fp ATTRIBUTE_UNUSED,
2312                      gfc_expr * l  ATTRIBUTE_UNUSED,
2313                      gfc_expr * to ATTRIBUTE_UNUSED,
2314                      gfc_expr * tp ATTRIBUTE_UNUSED)
2315 {
2316   return NULL;
2317 }
2318
2319
2320 gfc_expr *
2321 gfc_simplify_nearest (gfc_expr * x, gfc_expr * s)
2322 {
2323   gfc_expr *result;
2324   mpfr_t tmp;
2325   int direction, sgn;
2326
2327   if (x->expr_type != EXPR_CONSTANT || s->expr_type != EXPR_CONSTANT)
2328     return NULL;
2329
2330   gfc_set_model_kind (x->ts.kind);
2331   result = gfc_copy_expr (x);
2332
2333   direction = mpfr_sgn (s->value.real);
2334
2335   if (direction == 0)
2336     {
2337       gfc_error ("Second argument of NEAREST at %L may not be zero",
2338                  &s->where);
2339       gfc_free (result);
2340       return &gfc_bad_expr;
2341     }
2342
2343   /* TODO: Use mpfr_nextabove and mpfr_nextbelow once we move to a
2344      newer version of mpfr.  */
2345
2346   sgn = mpfr_sgn (x->value.real);
2347
2348   if (sgn == 0)
2349     {
2350       int k = gfc_validate_kind (BT_REAL, x->ts.kind, 0);
2351
2352       if (direction > 0)
2353         mpfr_add (result->value.real,
2354                   x->value.real, gfc_real_kinds[k].subnormal, GFC_RND_MODE);
2355       else
2356         mpfr_sub (result->value.real,
2357                   x->value.real, gfc_real_kinds[k].subnormal, GFC_RND_MODE);
2358     }
2359   else
2360     {
2361       if (sgn < 0)
2362         {
2363           direction = -direction;
2364           mpfr_neg (result->value.real, result->value.real, GFC_RND_MODE);
2365         }
2366
2367       if (direction > 0)
2368         mpfr_add_one_ulp (result->value.real, GFC_RND_MODE);
2369       else
2370         {
2371           /* In this case the exponent can shrink, which makes us skip
2372              over one number because we subtract one ulp with the
2373              larger exponent.  Thus we need to compensate for this.  */
2374           mpfr_init_set (tmp, result->value.real, GFC_RND_MODE);
2375
2376           mpfr_sub_one_ulp (result->value.real, GFC_RND_MODE);
2377           mpfr_add_one_ulp (result->value.real, GFC_RND_MODE);
2378
2379           /* If we're back to where we started, the spacing is one
2380              ulp, and we get the correct result by subtracting.  */
2381           if (mpfr_cmp (tmp, result->value.real) == 0)
2382             mpfr_sub_one_ulp (result->value.real, GFC_RND_MODE);
2383
2384           mpfr_clear (tmp);
2385         }
2386
2387       if (sgn < 0)
2388         mpfr_neg (result->value.real, result->value.real, GFC_RND_MODE);
2389     }
2390
2391   return range_check (result, "NEAREST");
2392 }
2393
2394
2395 static gfc_expr *
2396 simplify_nint (const char *name, gfc_expr * e, gfc_expr * k)
2397 {
2398   gfc_expr *itrunc, *result;
2399   int kind;
2400
2401   kind = get_kind (BT_INTEGER, k, name, gfc_default_integer_kind);
2402   if (kind == -1)
2403     return &gfc_bad_expr;
2404
2405   if (e->expr_type != EXPR_CONSTANT)
2406     return NULL;
2407
2408   result = gfc_constant_result (BT_INTEGER, kind, &e->where);
2409
2410   itrunc = gfc_copy_expr (e);
2411
2412   mpfr_round(itrunc->value.real, e->value.real);
2413
2414   gfc_mpfr_to_mpz (result->value.integer, itrunc->value.real);
2415
2416   gfc_free_expr (itrunc);
2417
2418   return range_check (result, name);
2419 }
2420
2421
2422 gfc_expr *
2423 gfc_simplify_nint (gfc_expr * e, gfc_expr * k)
2424 {
2425   return simplify_nint ("NINT", e, k);
2426 }
2427
2428
2429 gfc_expr *
2430 gfc_simplify_idnint (gfc_expr * e)
2431 {
2432   return simplify_nint ("IDNINT", e, NULL);
2433 }
2434
2435
2436 gfc_expr *
2437 gfc_simplify_not (gfc_expr * e)
2438 {
2439   gfc_expr *result;
2440   int i;
2441
2442   if (e->expr_type != EXPR_CONSTANT)
2443     return NULL;
2444
2445   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
2446
2447   mpz_com (result->value.integer, e->value.integer);
2448
2449   /* Because of how GMP handles numbers, the result must be ANDed with
2450      the max_int mask.  For radices <> 2, this will require change.  */
2451
2452   i = gfc_validate_kind (BT_INTEGER, e->ts.kind, false);
2453
2454   mpz_and (result->value.integer, result->value.integer,
2455            gfc_integer_kinds[i].max_int);
2456
2457   return range_check (result, "NOT");
2458 }
2459
2460
2461 gfc_expr *
2462 gfc_simplify_null (gfc_expr * mold)
2463 {
2464   gfc_expr *result;
2465
2466   result = gfc_get_expr ();
2467   result->expr_type = EXPR_NULL;
2468
2469   if (mold == NULL)
2470     result->ts.type = BT_UNKNOWN;
2471   else
2472     {
2473       result->ts = mold->ts;
2474       result->where = mold->where;
2475     }
2476
2477   return result;
2478 }
2479
2480
2481 gfc_expr *
2482 gfc_simplify_precision (gfc_expr * e)
2483 {
2484   gfc_expr *result;
2485   int i;
2486
2487   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
2488
2489   result = gfc_int_expr (gfc_real_kinds[i].precision);
2490   result->where = e->where;
2491
2492   return result;
2493 }
2494
2495
2496 gfc_expr *
2497 gfc_simplify_radix (gfc_expr * e)
2498 {
2499   gfc_expr *result;
2500   int i;
2501
2502   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
2503   switch (e->ts.type)
2504     {
2505     case BT_INTEGER:
2506       i = gfc_integer_kinds[i].radix;
2507       break;
2508
2509     case BT_REAL:
2510       i = gfc_real_kinds[i].radix;
2511       break;
2512
2513     default:
2514       gcc_unreachable ();
2515     }
2516
2517   result = gfc_int_expr (i);
2518   result->where = e->where;
2519
2520   return result;
2521 }
2522
2523
2524 gfc_expr *
2525 gfc_simplify_range (gfc_expr * e)
2526 {
2527   gfc_expr *result;
2528   int i;
2529   long j;
2530
2531   i = gfc_validate_kind (e->ts.type, e->ts.kind, false);
2532
2533   switch (e->ts.type)
2534     {
2535     case BT_INTEGER:
2536       j = gfc_integer_kinds[i].range;
2537       break;
2538
2539     case BT_REAL:
2540     case BT_COMPLEX:
2541       j = gfc_real_kinds[i].range;
2542       break;
2543
2544     default:
2545       gcc_unreachable ();
2546     }
2547
2548   result = gfc_int_expr (j);
2549   result->where = e->where;
2550
2551   return result;
2552 }
2553
2554
2555 gfc_expr *
2556 gfc_simplify_real (gfc_expr * e, gfc_expr * k)
2557 {
2558   gfc_expr *result;
2559   int kind;
2560
2561   if (e->ts.type == BT_COMPLEX)
2562     kind = get_kind (BT_REAL, k, "REAL", e->ts.kind);
2563   else
2564     kind = get_kind (BT_REAL, k, "REAL", gfc_default_real_kind);
2565
2566   if (kind == -1)
2567     return &gfc_bad_expr;
2568
2569   if (e->expr_type != EXPR_CONSTANT)
2570     return NULL;
2571
2572   switch (e->ts.type)
2573     {
2574     case BT_INTEGER:
2575       result = gfc_int2real (e, kind);
2576       break;
2577
2578     case BT_REAL:
2579       result = gfc_real2real (e, kind);
2580       break;
2581
2582     case BT_COMPLEX:
2583       result = gfc_complex2real (e, kind);
2584       break;
2585
2586     default:
2587       gfc_internal_error ("bad type in REAL");
2588       /* Not reached */
2589     }
2590
2591   return range_check (result, "REAL");
2592 }
2593
2594 gfc_expr *
2595 gfc_simplify_repeat (gfc_expr * e, gfc_expr * n)
2596 {
2597   gfc_expr *result;
2598   int i, j, len, ncopies, nlen;
2599
2600   if (e->expr_type != EXPR_CONSTANT || n->expr_type != EXPR_CONSTANT)
2601     return NULL;
2602
2603   if (n != NULL && (gfc_extract_int (n, &ncopies) != NULL || ncopies < 0))
2604     {
2605       gfc_error ("Invalid second argument of REPEAT at %L", &n->where);
2606       return &gfc_bad_expr;
2607     }
2608
2609   len = e->value.character.length;
2610   nlen = ncopies * len;
2611
2612   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
2613
2614   if (ncopies == 0)
2615     {
2616       result->value.character.string = gfc_getmem (1);
2617       result->value.character.length = 0;
2618       result->value.character.string[0] = '\0';
2619       return result;
2620     }
2621
2622   result->value.character.length = nlen;
2623   result->value.character.string = gfc_getmem (nlen + 1);
2624
2625   for (i = 0; i < ncopies; i++)
2626     for (j = 0; j < len; j++)
2627       result->value.character.string[j + i * len] =
2628         e->value.character.string[j];
2629
2630   result->value.character.string[nlen] = '\0';  /* For debugger */
2631   return result;
2632 }
2633
2634
2635 /* This one is a bear, but mainly has to do with shuffling elements.  */
2636
2637 gfc_expr *
2638 gfc_simplify_reshape (gfc_expr * source, gfc_expr * shape_exp,
2639                       gfc_expr * pad, gfc_expr * order_exp)
2640 {
2641
2642   int order[GFC_MAX_DIMENSIONS], shape[GFC_MAX_DIMENSIONS];
2643   int i, rank, npad, x[GFC_MAX_DIMENSIONS];
2644   gfc_constructor *head, *tail;
2645   mpz_t index, size;
2646   unsigned long j;
2647   size_t nsource;
2648   gfc_expr *e;
2649
2650   /* Unpack the shape array.  */
2651   if (source->expr_type != EXPR_ARRAY || !gfc_is_constant_expr (source))
2652     return NULL;
2653
2654   if (shape_exp->expr_type != EXPR_ARRAY || !gfc_is_constant_expr (shape_exp))
2655     return NULL;
2656
2657   if (pad != NULL
2658       && (pad->expr_type != EXPR_ARRAY
2659           || !gfc_is_constant_expr (pad)))
2660     return NULL;
2661
2662   if (order_exp != NULL
2663       && (order_exp->expr_type != EXPR_ARRAY
2664           || !gfc_is_constant_expr (order_exp)))
2665     return NULL;
2666
2667   mpz_init (index);
2668   rank = 0;
2669   head = tail = NULL;
2670
2671   for (;;)
2672     {
2673       e = gfc_get_array_element (shape_exp, rank);
2674       if (e == NULL)
2675         break;
2676
2677       if (gfc_extract_int (e, &shape[rank]) != NULL)
2678         {
2679           gfc_error ("Integer too large in shape specification at %L",
2680                      &e->where);
2681           gfc_free_expr (e);
2682           goto bad_reshape;
2683         }
2684
2685       gfc_free_expr (e);
2686
2687       if (rank >= GFC_MAX_DIMENSIONS)
2688         {
2689           gfc_error ("Too many dimensions in shape specification for RESHAPE "
2690                      "at %L", &e->where);
2691
2692           goto bad_reshape;
2693         }
2694
2695       if (shape[rank] < 0)
2696         {
2697           gfc_error ("Shape specification at %L cannot be negative",
2698                      &e->where);
2699           goto bad_reshape;
2700         }
2701
2702       rank++;
2703     }
2704
2705   if (rank == 0)
2706     {
2707       gfc_error ("Shape specification at %L cannot be the null array",
2708                  &shape_exp->where);
2709       goto bad_reshape;
2710     }
2711
2712   /* Now unpack the order array if present.  */
2713   if (order_exp == NULL)
2714     {
2715       for (i = 0; i < rank; i++)
2716         order[i] = i;
2717
2718     }
2719   else
2720     {
2721
2722       for (i = 0; i < rank; i++)
2723         x[i] = 0;
2724
2725       for (i = 0; i < rank; i++)
2726         {
2727           e = gfc_get_array_element (order_exp, i);
2728           if (e == NULL)
2729             {
2730               gfc_error
2731                 ("ORDER parameter of RESHAPE at %L is not the same size "
2732                  "as SHAPE parameter", &order_exp->where);
2733               goto bad_reshape;
2734             }
2735
2736           if (gfc_extract_int (e, &order[i]) != NULL)
2737             {
2738               gfc_error ("Error in ORDER parameter of RESHAPE at %L",
2739                          &e->where);
2740               gfc_free_expr (e);
2741               goto bad_reshape;
2742             }
2743
2744           gfc_free_expr (e);
2745
2746           if (order[i] < 1 || order[i] > rank)
2747             {
2748               gfc_error ("ORDER parameter of RESHAPE at %L is out of range",
2749                          &e->where);
2750               goto bad_reshape;
2751             }
2752
2753           order[i]--;
2754
2755           if (x[order[i]])
2756             {
2757               gfc_error ("Invalid permutation in ORDER parameter at %L",
2758                          &e->where);
2759               goto bad_reshape;
2760             }
2761
2762           x[order[i]] = 1;
2763         }
2764     }
2765
2766   /* Count the elements in the source and padding arrays.  */
2767
2768   npad = 0;
2769   if (pad != NULL)
2770     {
2771       gfc_array_size (pad, &size);
2772       npad = mpz_get_ui (size);
2773       mpz_clear (size);
2774     }
2775
2776   gfc_array_size (source, &size);
2777   nsource = mpz_get_ui (size);
2778   mpz_clear (size);
2779
2780   /* If it weren't for that pesky permutation we could just loop
2781      through the source and round out any shortage with pad elements.
2782      But no, someone just had to have the compiler do something the
2783      user should be doing.  */
2784
2785   for (i = 0; i < rank; i++)
2786     x[i] = 0;
2787
2788   for (;;)
2789     {
2790       /* Figure out which element to extract.  */
2791       mpz_set_ui (index, 0);
2792
2793       for (i = rank - 1; i >= 0; i--)
2794         {
2795           mpz_add_ui (index, index, x[order[i]]);
2796           if (i != 0)
2797             mpz_mul_ui (index, index, shape[order[i - 1]]);
2798         }
2799
2800       if (mpz_cmp_ui (index, INT_MAX) > 0)
2801         gfc_internal_error ("Reshaped array too large at %L", &e->where);
2802
2803       j = mpz_get_ui (index);
2804
2805       if (j < nsource)
2806         e = gfc_get_array_element (source, j);
2807       else
2808         {
2809           j = j - nsource;
2810
2811           if (npad == 0)
2812             {
2813               gfc_error
2814                 ("PAD parameter required for short SOURCE parameter at %L",
2815                  &source->where);
2816               goto bad_reshape;
2817             }
2818
2819           j = j % npad;
2820           e = gfc_get_array_element (pad, j);
2821         }
2822
2823       if (head == NULL)
2824         head = tail = gfc_get_constructor ();
2825       else
2826         {
2827           tail->next = gfc_get_constructor ();
2828           tail = tail->next;
2829         }
2830
2831       if (e == NULL)
2832         goto bad_reshape;
2833
2834       tail->where = e->where;
2835       tail->expr = e;
2836
2837       /* Calculate the next element.  */
2838       i = 0;
2839
2840 inc:
2841       if (++x[i] < shape[i])
2842         continue;
2843       x[i++] = 0;
2844       if (i < rank)
2845         goto inc;
2846
2847       break;
2848     }
2849
2850   mpz_clear (index);
2851
2852   e = gfc_get_expr ();
2853   e->where = source->where;
2854   e->expr_type = EXPR_ARRAY;
2855   e->value.constructor = head;
2856   e->shape = gfc_get_shape (rank);
2857
2858   for (i = 0; i < rank; i++)
2859     mpz_init_set_ui (e->shape[i], shape[i]);
2860
2861   e->ts = head->expr->ts;
2862   e->rank = rank;
2863
2864   return e;
2865
2866 bad_reshape:
2867   gfc_free_constructor (head);
2868   mpz_clear (index);
2869   return &gfc_bad_expr;
2870 }
2871
2872
2873 gfc_expr *
2874 gfc_simplify_rrspacing (gfc_expr * x)
2875 {
2876   gfc_expr *result;
2877   mpfr_t absv, log2, exp, frac, pow2;
2878   int i, p;
2879
2880   if (x->expr_type != EXPR_CONSTANT)
2881     return NULL;
2882
2883   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
2884
2885   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
2886
2887   p = gfc_real_kinds[i].digits;
2888
2889   gfc_set_model_kind (x->ts.kind);
2890
2891   if (mpfr_sgn (x->value.real) == 0)
2892     {
2893       mpfr_ui_div (result->value.real, 1, gfc_real_kinds[i].tiny, GFC_RND_MODE);
2894       return result;
2895     }
2896
2897   mpfr_init (log2);
2898   mpfr_init (absv);
2899   mpfr_init (frac);
2900   mpfr_init (pow2);
2901
2902   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
2903   mpfr_log2 (log2, absv, GFC_RND_MODE);
2904
2905   mpfr_trunc (log2, log2);
2906   mpfr_add_ui (exp, log2, 1, GFC_RND_MODE);
2907
2908   mpfr_ui_pow (pow2, 2, exp, GFC_RND_MODE);
2909   mpfr_div (frac, absv, pow2, GFC_RND_MODE);
2910
2911   mpfr_mul_2exp (result->value.real, frac, (unsigned long)p, GFC_RND_MODE);
2912
2913   mpfr_clear (log2);
2914   mpfr_clear (absv);
2915   mpfr_clear (frac);
2916   mpfr_clear (pow2);
2917
2918   return range_check (result, "RRSPACING");
2919 }
2920
2921
2922 gfc_expr *
2923 gfc_simplify_scale (gfc_expr * x, gfc_expr * i)
2924 {
2925   int k, neg_flag, power, exp_range;
2926   mpfr_t scale, radix;
2927   gfc_expr *result;
2928
2929   if (x->expr_type != EXPR_CONSTANT || i->expr_type != EXPR_CONSTANT)
2930     return NULL;
2931
2932   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
2933
2934   if (mpfr_sgn (x->value.real) == 0)
2935     {
2936       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
2937       return result;
2938     }
2939
2940   k = gfc_validate_kind (BT_REAL, x->ts.kind, false);
2941
2942   exp_range = gfc_real_kinds[k].max_exponent - gfc_real_kinds[k].min_exponent;
2943
2944   /* This check filters out values of i that would overflow an int.  */
2945   if (mpz_cmp_si (i->value.integer, exp_range + 2) > 0
2946       || mpz_cmp_si (i->value.integer, -exp_range - 2) < 0)
2947     {
2948       gfc_error ("Result of SCALE overflows its kind at %L", &result->where);
2949       return &gfc_bad_expr;
2950     }
2951
2952   /* Compute scale = radix ** power.  */
2953   power = mpz_get_si (i->value.integer);
2954
2955   if (power >= 0)
2956     neg_flag = 0;
2957   else
2958     {
2959       neg_flag = 1;
2960       power = -power;
2961     }
2962
2963   gfc_set_model_kind (x->ts.kind);
2964   mpfr_init (scale);
2965   mpfr_init (radix);
2966   mpfr_set_ui (radix, gfc_real_kinds[k].radix, GFC_RND_MODE);
2967   mpfr_pow_ui (scale, radix, power, GFC_RND_MODE);
2968
2969   if (neg_flag)
2970     mpfr_div (result->value.real, x->value.real, scale, GFC_RND_MODE);
2971   else
2972     mpfr_mul (result->value.real, x->value.real, scale, GFC_RND_MODE);
2973
2974   mpfr_clear (scale);
2975   mpfr_clear (radix);
2976
2977   return range_check (result, "SCALE");
2978 }
2979
2980
2981 gfc_expr *
2982 gfc_simplify_scan (gfc_expr * e, gfc_expr * c, gfc_expr * b)
2983 {
2984   gfc_expr *result;
2985   int back;
2986   size_t i;
2987   size_t indx, len, lenc;
2988
2989   if (e->expr_type != EXPR_CONSTANT || c->expr_type != EXPR_CONSTANT)
2990     return NULL;
2991
2992   if (b != NULL && b->value.logical != 0)
2993     back = 1;
2994   else
2995     back = 0;
2996
2997   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
2998                                 &e->where);
2999
3000   len = e->value.character.length;
3001   lenc = c->value.character.length;
3002
3003   if (len == 0 || lenc == 0)
3004     {
3005       indx = 0;
3006     }
3007   else
3008     {
3009       if (back == 0)
3010         {
3011           indx =
3012             strcspn (e->value.character.string, c->value.character.string) + 1;
3013           if (indx > len)
3014             indx = 0;
3015         }
3016       else
3017         {
3018           i = 0;
3019           for (indx = len; indx > 0; indx--)
3020             {
3021               for (i = 0; i < lenc; i++)
3022                 {
3023                   if (c->value.character.string[i]
3024                         == e->value.character.string[indx - 1])
3025                     break;
3026                 }
3027               if (i < lenc)
3028                 break;
3029             }
3030         }
3031     }
3032   mpz_set_ui (result->value.integer, indx);
3033   return range_check (result, "SCAN");
3034 }
3035
3036
3037 gfc_expr *
3038 gfc_simplify_selected_int_kind (gfc_expr * e)
3039 {
3040   int i, kind, range;
3041   gfc_expr *result;
3042
3043   if (e->expr_type != EXPR_CONSTANT || gfc_extract_int (e, &range) != NULL)
3044     return NULL;
3045
3046   kind = INT_MAX;
3047
3048   for (i = 0; gfc_integer_kinds[i].kind != 0; i++)
3049     if (gfc_integer_kinds[i].range >= range
3050         && gfc_integer_kinds[i].kind < kind)
3051       kind = gfc_integer_kinds[i].kind;
3052
3053   if (kind == INT_MAX)
3054     kind = -1;
3055
3056   result = gfc_int_expr (kind);
3057   result->where = e->where;
3058
3059   return result;
3060 }
3061
3062
3063 gfc_expr *
3064 gfc_simplify_selected_real_kind (gfc_expr * p, gfc_expr * q)
3065 {
3066   int range, precision, i, kind, found_precision, found_range;
3067   gfc_expr *result;
3068
3069   if (p == NULL)
3070     precision = 0;
3071   else
3072     {
3073       if (p->expr_type != EXPR_CONSTANT
3074           || gfc_extract_int (p, &precision) != NULL)
3075         return NULL;
3076     }
3077
3078   if (q == NULL)
3079     range = 0;
3080   else
3081     {
3082       if (q->expr_type != EXPR_CONSTANT
3083           || gfc_extract_int (q, &range) != NULL)
3084         return NULL;
3085     }
3086
3087   kind = INT_MAX;
3088   found_precision = 0;
3089   found_range = 0;
3090
3091   for (i = 0; gfc_real_kinds[i].kind != 0; i++)
3092     {
3093       if (gfc_real_kinds[i].precision >= precision)
3094         found_precision = 1;
3095
3096       if (gfc_real_kinds[i].range >= range)
3097         found_range = 1;
3098
3099       if (gfc_real_kinds[i].precision >= precision
3100           && gfc_real_kinds[i].range >= range && gfc_real_kinds[i].kind < kind)
3101         kind = gfc_real_kinds[i].kind;
3102     }
3103
3104   if (kind == INT_MAX)
3105     {
3106       kind = 0;
3107
3108       if (!found_precision)
3109         kind = -1;
3110       if (!found_range)
3111         kind -= 2;
3112     }
3113
3114   result = gfc_int_expr (kind);
3115   result->where = (p != NULL) ? p->where : q->where;
3116
3117   return result;
3118 }
3119
3120
3121 gfc_expr *
3122 gfc_simplify_set_exponent (gfc_expr * x, gfc_expr * i)
3123 {
3124   gfc_expr *result;
3125   mpfr_t exp, absv, log2, pow2, frac;
3126   unsigned long exp2;
3127
3128   if (x->expr_type != EXPR_CONSTANT || i->expr_type != EXPR_CONSTANT)
3129     return NULL;
3130
3131   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
3132
3133   gfc_set_model_kind (x->ts.kind);
3134
3135   if (mpfr_sgn (x->value.real) == 0)
3136     {
3137       mpfr_set_ui (result->value.real, 0, GFC_RND_MODE);
3138       return result;
3139     }
3140
3141   mpfr_init (absv);
3142   mpfr_init (log2);
3143   mpfr_init (exp);
3144   mpfr_init (pow2);
3145   mpfr_init (frac);
3146
3147   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
3148   mpfr_log2 (log2, absv, GFC_RND_MODE);
3149
3150   mpfr_trunc (log2, log2);
3151   mpfr_add_ui (exp, log2, 1, GFC_RND_MODE);
3152
3153   /* Old exponent value, and fraction.  */
3154   mpfr_ui_pow (pow2, 2, exp, GFC_RND_MODE);
3155
3156   mpfr_div (frac, absv, pow2, GFC_RND_MODE);
3157
3158   /* New exponent.  */
3159   exp2 = (unsigned long) mpz_get_d (i->value.integer);
3160   mpfr_mul_2exp (result->value.real, frac, exp2, GFC_RND_MODE);
3161
3162   mpfr_clear (absv);
3163   mpfr_clear (log2);
3164   mpfr_clear (pow2);
3165   mpfr_clear (frac);
3166
3167   return range_check (result, "SET_EXPONENT");
3168 }
3169
3170
3171 gfc_expr *
3172 gfc_simplify_shape (gfc_expr * source)
3173 {
3174   mpz_t shape[GFC_MAX_DIMENSIONS];
3175   gfc_expr *result, *e, *f;
3176   gfc_array_ref *ar;
3177   int n;
3178   try t;
3179
3180   if (source->rank == 0 || source->expr_type != EXPR_VARIABLE)
3181     return NULL;
3182
3183   result = gfc_start_constructor (BT_INTEGER, gfc_default_integer_kind,
3184                                   &source->where);
3185
3186   ar = gfc_find_array_ref (source);
3187
3188   t = gfc_array_ref_shape (ar, shape);
3189
3190   for (n = 0; n < source->rank; n++)
3191     {
3192       e = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
3193                                &source->where);
3194
3195       if (t == SUCCESS)
3196         {
3197           mpz_set (e->value.integer, shape[n]);
3198           mpz_clear (shape[n]);
3199         }
3200       else
3201         {
3202           mpz_set_ui (e->value.integer, n + 1);
3203
3204           f = gfc_simplify_size (source, e);
3205           gfc_free_expr (e);
3206           if (f == NULL)
3207             {
3208               gfc_free_expr (result);
3209               return NULL;
3210             }
3211           else
3212             {
3213               e = f;
3214             }
3215         }
3216
3217       gfc_append_constructor (result, e);
3218     }
3219
3220   return result;
3221 }
3222
3223
3224 gfc_expr *
3225 gfc_simplify_size (gfc_expr * array, gfc_expr * dim)
3226 {
3227   mpz_t size;
3228   gfc_expr *result;
3229   int d;
3230
3231   if (dim == NULL)
3232     {
3233       if (gfc_array_size (array, &size) == FAILURE)
3234         return NULL;
3235     }
3236   else
3237     {
3238       if (dim->expr_type != EXPR_CONSTANT)
3239         return NULL;
3240
3241       d = mpz_get_ui (dim->value.integer) - 1;
3242       if (gfc_array_dimen_size (array, d, &size) == FAILURE)
3243         return NULL;
3244     }
3245
3246   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
3247                                 &array->where);
3248
3249   mpz_set (result->value.integer, size);
3250
3251   return result;
3252 }
3253
3254
3255 gfc_expr *
3256 gfc_simplify_sign (gfc_expr * x, gfc_expr * y)
3257 {
3258   gfc_expr *result;
3259
3260   if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
3261     return NULL;
3262
3263   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3264
3265   switch (x->ts.type)
3266     {
3267     case BT_INTEGER:
3268       mpz_abs (result->value.integer, x->value.integer);
3269       if (mpz_sgn (y->value.integer) < 0)
3270         mpz_neg (result->value.integer, result->value.integer);
3271
3272       break;
3273
3274     case BT_REAL:
3275       /* TODO: Handle -0.0 and +0.0 correctly on machines that support
3276          it.  */
3277       mpfr_abs (result->value.real, x->value.real, GFC_RND_MODE);
3278       if (mpfr_sgn (y->value.real) < 0)
3279         mpfr_neg (result->value.real, result->value.real, GFC_RND_MODE);
3280
3281       break;
3282
3283     default:
3284       gfc_internal_error ("Bad type in gfc_simplify_sign");
3285     }
3286
3287   return result;
3288 }
3289
3290
3291 gfc_expr *
3292 gfc_simplify_sin (gfc_expr * x)
3293 {
3294   gfc_expr *result;
3295   mpfr_t xp, xq;
3296
3297   if (x->expr_type != EXPR_CONSTANT)
3298     return NULL;
3299
3300   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3301
3302   switch (x->ts.type)
3303     {
3304     case BT_REAL:
3305       mpfr_sin (result->value.real, x->value.real, GFC_RND_MODE);
3306       break;
3307
3308     case BT_COMPLEX:
3309       gfc_set_model (x->value.real);
3310       mpfr_init (xp);
3311       mpfr_init (xq);
3312
3313       mpfr_sin  (xp, x->value.complex.r, GFC_RND_MODE);
3314       mpfr_cosh (xq, x->value.complex.i, GFC_RND_MODE);
3315       mpfr_mul (result->value.complex.r, xp, xq, GFC_RND_MODE);
3316
3317       mpfr_cos  (xp, x->value.complex.r, GFC_RND_MODE);
3318       mpfr_sinh (xq, x->value.complex.i, GFC_RND_MODE);
3319       mpfr_mul (result->value.complex.i, xp, xq, GFC_RND_MODE);
3320
3321       mpfr_clear (xp);
3322       mpfr_clear (xq);
3323       break;
3324
3325     default:
3326       gfc_internal_error ("in gfc_simplify_sin(): Bad type");
3327     }
3328
3329   return range_check (result, "SIN");
3330 }
3331
3332
3333 gfc_expr *
3334 gfc_simplify_sinh (gfc_expr * x)
3335 {
3336   gfc_expr *result;
3337
3338   if (x->expr_type != EXPR_CONSTANT)
3339     return NULL;
3340
3341   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3342
3343   mpfr_sinh(result->value.real, x->value.real, GFC_RND_MODE);
3344
3345   return range_check (result, "SINH");
3346 }
3347
3348
3349 /* The argument is always a double precision real that is converted to
3350    single precision.  TODO: Rounding!  */
3351
3352 gfc_expr *
3353 gfc_simplify_sngl (gfc_expr * a)
3354 {
3355   gfc_expr *result;
3356
3357   if (a->expr_type != EXPR_CONSTANT)
3358     return NULL;
3359
3360   result = gfc_real2real (a, gfc_default_real_kind);
3361   return range_check (result, "SNGL");
3362 }
3363
3364
3365 gfc_expr *
3366 gfc_simplify_spacing (gfc_expr * x)
3367 {
3368   gfc_expr *result;
3369   mpfr_t absv, log2;
3370   long diff;
3371   int i, p;
3372
3373   if (x->expr_type != EXPR_CONSTANT)
3374     return NULL;
3375
3376   i = gfc_validate_kind (x->ts.type, x->ts.kind, false);
3377
3378   p = gfc_real_kinds[i].digits;
3379
3380   result = gfc_constant_result (BT_REAL, x->ts.kind, &x->where);
3381
3382   gfc_set_model_kind (x->ts.kind);
3383
3384   if (mpfr_sgn (x->value.real) == 0)
3385     {
3386       mpfr_set (result->value.real, gfc_real_kinds[i].tiny, GFC_RND_MODE);
3387       return result;
3388     }
3389
3390   mpfr_init (log2);
3391   mpfr_init (absv);
3392
3393   mpfr_abs (absv, x->value.real, GFC_RND_MODE);
3394   mpfr_log2 (log2, absv, GFC_RND_MODE);
3395   mpfr_trunc (log2, log2);
3396
3397   mpfr_add_ui (log2, log2, 1, GFC_RND_MODE);
3398
3399   /* FIXME: We should be using mpfr_get_si here, but this function is
3400      not available with the version of mpfr distributed with gmp (as of
3401      2004-09-17). Replace once mpfr has been imported into the gcc cvs
3402      tree.  */
3403   diff = (long)mpfr_get_d (log2, GFC_RND_MODE) - (long)p;
3404   mpfr_set_ui (result->value.real, 1, GFC_RND_MODE);
3405   mpfr_mul_2si (result->value.real, result->value.real, diff, GFC_RND_MODE);
3406
3407   mpfr_clear (log2);
3408   mpfr_clear (absv);
3409
3410   if (mpfr_cmp (result->value.real, gfc_real_kinds[i].tiny) < 0)
3411     mpfr_set (result->value.real, gfc_real_kinds[i].tiny, GFC_RND_MODE);
3412
3413   return range_check (result, "SPACING");
3414 }
3415
3416
3417 gfc_expr *
3418 gfc_simplify_sqrt (gfc_expr * e)
3419 {
3420   gfc_expr *result;
3421   mpfr_t ac, ad, s, t, w;
3422
3423   if (e->expr_type != EXPR_CONSTANT)
3424     return NULL;
3425
3426   result = gfc_constant_result (e->ts.type, e->ts.kind, &e->where);
3427
3428   switch (e->ts.type)
3429     {
3430     case BT_REAL:
3431       if (mpfr_cmp_si (e->value.real, 0) < 0)
3432         goto negative_arg;
3433       mpfr_sqrt (result->value.real, e->value.real, GFC_RND_MODE);
3434
3435       break;
3436
3437     case BT_COMPLEX:
3438       /* Formula taken from Numerical Recipes to avoid over- and
3439          underflow.  */
3440
3441       gfc_set_model (e->value.real);
3442       mpfr_init (ac);
3443       mpfr_init (ad);
3444       mpfr_init (s);
3445       mpfr_init (t);
3446       mpfr_init (w);
3447
3448       if (mpfr_cmp_ui (e->value.complex.r, 0) == 0
3449           && mpfr_cmp_ui (e->value.complex.i, 0) == 0)
3450         {
3451
3452           mpfr_set_ui (result->value.complex.r, 0, GFC_RND_MODE);
3453           mpfr_set_ui (result->value.complex.i, 0, GFC_RND_MODE);
3454           break;
3455         }
3456
3457       mpfr_abs (ac, e->value.complex.r, GFC_RND_MODE);
3458       mpfr_abs (ad, e->value.complex.i, GFC_RND_MODE);
3459
3460       if (mpfr_cmp (ac, ad) >= 0)
3461         {
3462           mpfr_div (t, e->value.complex.i, e->value.complex.r, GFC_RND_MODE);
3463           mpfr_mul (t, t, t, GFC_RND_MODE);
3464           mpfr_add_ui (t, t, 1, GFC_RND_MODE);
3465           mpfr_sqrt (t, t, GFC_RND_MODE);
3466           mpfr_add_ui (t, t, 1, GFC_RND_MODE);
3467           mpfr_div_ui (t, t, 2, GFC_RND_MODE);
3468           mpfr_sqrt (t, t, GFC_RND_MODE);
3469           mpfr_sqrt (s, ac, GFC_RND_MODE);
3470           mpfr_mul (w, s, t, GFC_RND_MODE);
3471         }
3472       else
3473         {
3474           mpfr_div (s, e->value.complex.r, e->value.complex.i, GFC_RND_MODE);
3475           mpfr_mul (t, s, s, GFC_RND_MODE);
3476           mpfr_add_ui (t, t, 1, GFC_RND_MODE);
3477           mpfr_sqrt (t, t, GFC_RND_MODE);
3478           mpfr_abs (s, s, GFC_RND_MODE);
3479           mpfr_add (t, t, s, GFC_RND_MODE);
3480           mpfr_div_ui (t, t, 2, GFC_RND_MODE);
3481           mpfr_sqrt (t, t, GFC_RND_MODE);
3482           mpfr_sqrt (s, ad, GFC_RND_MODE);
3483           mpfr_mul (w, s, t, GFC_RND_MODE);
3484         }
3485
3486       if (mpfr_cmp_ui (w, 0) != 0 && mpfr_cmp_ui (e->value.complex.r, 0) >= 0)
3487         {
3488           mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
3489           mpfr_div (result->value.complex.i, e->value.complex.i, t, GFC_RND_MODE);
3490           mpfr_set (result->value.complex.r, w, GFC_RND_MODE);
3491         }
3492       else if (mpfr_cmp_ui (w, 0) != 0
3493                && mpfr_cmp_ui (e->value.complex.r, 0) < 0
3494                && mpfr_cmp_ui (e->value.complex.i, 0) >= 0)
3495         {
3496           mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
3497           mpfr_div (result->value.complex.r, e->value.complex.i, t, GFC_RND_MODE);
3498           mpfr_set (result->value.complex.i, w, GFC_RND_MODE);
3499         }
3500       else if (mpfr_cmp_ui (w, 0) != 0
3501                && mpfr_cmp_ui (e->value.complex.r, 0) < 0
3502                && mpfr_cmp_ui (e->value.complex.i, 0) < 0)
3503         {
3504           mpfr_mul_ui (t, w, 2, GFC_RND_MODE);
3505           mpfr_div (result->value.complex.r, ad, t, GFC_RND_MODE);
3506           mpfr_neg (w, w, GFC_RND_MODE);
3507           mpfr_set (result->value.complex.i, w, GFC_RND_MODE);
3508         }
3509       else
3510         gfc_internal_error ("invalid complex argument of SQRT at %L",
3511                             &e->where);
3512
3513       mpfr_clear (s);
3514       mpfr_clear (t);
3515       mpfr_clear (ac);
3516       mpfr_clear (ad);
3517       mpfr_clear (w);
3518
3519       break;
3520
3521     default:
3522       gfc_internal_error ("invalid argument of SQRT at %L", &e->where);
3523     }
3524
3525   return range_check (result, "SQRT");
3526
3527 negative_arg:
3528   gfc_free_expr (result);
3529   gfc_error ("Argument of SQRT at %L has a negative value", &e->where);
3530   return &gfc_bad_expr;
3531 }
3532
3533
3534 gfc_expr *
3535 gfc_simplify_tan (gfc_expr * x)
3536 {
3537   int i;
3538   gfc_expr *result;
3539
3540   if (x->expr_type != EXPR_CONSTANT)
3541     return NULL;
3542
3543   i = gfc_validate_kind (BT_REAL, x->ts.kind, false);
3544
3545   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3546
3547   mpfr_tan (result->value.real, x->value.real, GFC_RND_MODE);
3548
3549   return range_check (result, "TAN");
3550 }
3551
3552
3553 gfc_expr *
3554 gfc_simplify_tanh (gfc_expr * x)
3555 {
3556   gfc_expr *result;
3557
3558   if (x->expr_type != EXPR_CONSTANT)
3559     return NULL;
3560
3561   result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
3562
3563   mpfr_tanh (result->value.real, x->value.real, GFC_RND_MODE);
3564
3565   return range_check (result, "TANH");
3566
3567 }
3568
3569
3570 gfc_expr *
3571 gfc_simplify_tiny (gfc_expr * e)
3572 {
3573   gfc_expr *result;
3574   int i;
3575
3576   i = gfc_validate_kind (BT_REAL, e->ts.kind, false);
3577
3578   result = gfc_constant_result (BT_REAL, e->ts.kind, &e->where);
3579   mpfr_set (result->value.real, gfc_real_kinds[i].tiny, GFC_RND_MODE);
3580
3581   return result;
3582 }
3583
3584
3585 gfc_expr *
3586 gfc_simplify_trim (gfc_expr * e)
3587 {
3588   gfc_expr *result;
3589   int count, i, len, lentrim;
3590
3591   if (e->expr_type != EXPR_CONSTANT)
3592     return NULL;
3593
3594   len = e->value.character.length;
3595
3596   result = gfc_constant_result (BT_CHARACTER, e->ts.kind, &e->where);
3597
3598   for (count = 0, i = 1; i <= len; ++i)
3599     {
3600       if (e->value.character.string[len - i] == ' ')
3601         count++;
3602       else
3603         break;
3604     }
3605
3606   lentrim = len - count;
3607
3608   result->value.character.length = lentrim;
3609   result->value.character.string = gfc_getmem (lentrim + 1);
3610
3611   for (i = 0; i < lentrim; i++)
3612     result->value.character.string[i] = e->value.character.string[i];
3613
3614   result->value.character.string[lentrim] = '\0';       /* For debugger */
3615
3616   return result;
3617 }
3618
3619
3620 gfc_expr *
3621 gfc_simplify_ubound (gfc_expr * array, gfc_expr * dim)
3622 {
3623   return simplify_bound (array, dim, 1);
3624 }
3625
3626
3627 gfc_expr *
3628 gfc_simplify_verify (gfc_expr * s, gfc_expr * set, gfc_expr * b)
3629 {
3630   gfc_expr *result;
3631   int back;
3632   size_t index, len, lenset;
3633   size_t i;
3634
3635   if (s->expr_type != EXPR_CONSTANT || set->expr_type != EXPR_CONSTANT)
3636     return NULL;
3637
3638   if (b != NULL && b->value.logical != 0)
3639     back = 1;
3640   else
3641     back = 0;
3642
3643   result = gfc_constant_result (BT_INTEGER, gfc_default_integer_kind,
3644                                 &s->where);
3645
3646   len = s->value.character.length;
3647   lenset = set->value.character.length;
3648
3649   if (len == 0)
3650     {
3651       mpz_set_ui (result->value.integer, 0);
3652       return result;
3653     }
3654
3655   if (back == 0)
3656     {
3657       if (lenset == 0)
3658         {
3659           mpz_set_ui (result->value.integer, len);
3660           return result;
3661         }
3662
3663       index =
3664         strspn (s->value.character.string, set->value.character.string) + 1;
3665       if (index > len)
3666         index = 0;
3667
3668     }
3669   else
3670     {
3671       if (lenset == 0)
3672         {
3673           mpz_set_ui (result->value.integer, 1);
3674           return result;
3675         }
3676       for (index = len; index > 0; index --)
3677         {
3678           for (i = 0; i < lenset; i++)
3679             {
3680               if (s->value.character.string[index - 1]
3681                     == set->value.character.string[i])
3682                 break;
3683             }
3684           if (i == lenset)
3685             break;
3686         }
3687     }
3688
3689   mpz_set_ui (result->value.integer, index);
3690   return result;
3691 }
3692
3693 /****************** Constant simplification *****************/
3694
3695 /* Master function to convert one constant to another.  While this is
3696    used as a simplification function, it requires the destination type
3697    and kind information which is supplied by a special case in
3698    do_simplify().  */
3699
3700 gfc_expr *
3701 gfc_convert_constant (gfc_expr * e, bt type, int kind)
3702 {
3703   gfc_expr *g, *result, *(*f) (gfc_expr *, int);
3704   gfc_constructor *head, *c, *tail = NULL;
3705
3706   switch (e->ts.type)
3707     {
3708     case BT_INTEGER:
3709       switch (type)
3710         {
3711         case BT_INTEGER:
3712           f = gfc_int2int;
3713           break;
3714         case BT_REAL:
3715           f = gfc_int2real;
3716           break;
3717         case BT_COMPLEX:
3718           f = gfc_int2complex;
3719           break;
3720         case BT_LOGICAL:
3721           f = gfc_int2log;
3722           break;
3723         default:
3724           goto oops;
3725         }
3726       break;
3727
3728     case BT_REAL:
3729       switch (type)
3730         {
3731         case BT_INTEGER:
3732           f = gfc_real2int;
3733           break;
3734         case BT_REAL:
3735           f = gfc_real2real;
3736           break;
3737         case BT_COMPLEX:
3738           f = gfc_real2complex;
3739           break;
3740         default:
3741           goto oops;
3742         }
3743       break;
3744
3745     case BT_COMPLEX:
3746       switch (type)
3747         {
3748         case BT_INTEGER:
3749           f = gfc_complex2int;
3750           break;
3751         case BT_REAL:
3752           f = gfc_complex2real;
3753           break;
3754         case BT_COMPLEX:
3755           f = gfc_complex2complex;
3756           break;
3757
3758         default:
3759           goto oops;
3760         }
3761       break;
3762
3763     case BT_LOGICAL:
3764       switch (type)
3765         {
3766         case BT_INTEGER:
3767           f = gfc_log2int;
3768           break;
3769         case BT_LOGICAL:
3770           f = gfc_log2log;
3771           break;
3772         default:
3773           goto oops;
3774         }
3775       break;
3776
3777     default:
3778     oops:
3779       gfc_internal_error ("gfc_convert_constant(): Unexpected type");
3780     }
3781
3782   result = NULL;
3783
3784   switch (e->expr_type)
3785     {
3786     case EXPR_CONSTANT:
3787       result = f (e, kind);
3788       if (result == NULL)
3789         return &gfc_bad_expr;
3790       break;
3791
3792     case EXPR_ARRAY:
3793       if (!gfc_is_constant_expr (e))
3794         break;
3795
3796       head = NULL;
3797
3798       for (c = e->value.constructor; c; c = c->next)
3799         {
3800           if (head == NULL)
3801             head = tail = gfc_get_constructor ();
3802           else
3803             {
3804               tail->next = gfc_get_constructor ();
3805               tail = tail->next;
3806             }
3807
3808           tail->where = c->where;
3809
3810           if (c->iterator == NULL)
3811             tail->expr = f (c->expr, kind);
3812           else
3813             {
3814               g = gfc_convert_constant (c->expr, type, kind);
3815               if (g == &gfc_bad_expr)
3816                 return g;
3817               tail->expr = g;
3818             }
3819
3820           if (tail->expr == NULL)
3821             {
3822               gfc_free_constructor (head);
3823               return NULL;
3824             }
3825         }
3826
3827       result = gfc_get_expr ();
3828       result->ts.type = type;
3829       result->ts.kind = kind;
3830       result->expr_type = EXPR_ARRAY;
3831       result->value.constructor = head;
3832       result->shape = gfc_copy_shape (e->shape, e->rank);
3833       result->where = e->where;
3834       result->rank = e->rank;
3835       break;
3836
3837     default:
3838       break;
3839     }
3840
3841   return result;
3842 }
3843
3844
3845 /****************** Helper functions ***********************/
3846
3847 /* Given a collating table, create the inverse table.  */
3848
3849 static void
3850 invert_table (const int *table, int *xtable)
3851 {
3852   int i;
3853
3854   for (i = 0; i < 256; i++)
3855     xtable[i] = 0;
3856
3857   for (i = 0; i < 256; i++)
3858     xtable[table[i]] = i;
3859 }
3860
3861
3862 void
3863 gfc_simplify_init_1 (void)
3864 {
3865
3866   invert_table (ascii_table, xascii_table);
3867 }