2 Copyright (C) 2000, 2001, 2002, 2003, 2004 Free Software Foundation,
4 Contributed by Andy Vaught
6 This file is part of GCC.
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
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
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
23 /* Since target arithmetic must be done on the host, there has to
24 be some way of evaluating arithmetic expressions as the host
25 would evaluate them. We use the GNU MP library to do arithmetic,
26 and this file provides the interface. */
35 mpf_t pi, half_pi, two_pi, e;
37 /* The gfc_(integer|real)_kinds[] structures have everything the front
38 end needs to know about integers and real numbers on the target.
39 Other entries of the structure are calculated from these values.
40 The first entry is the default kind, the second entry of the real
41 structure is the default double kind. */
43 #define MPZ_NULL {{0,0,0}}
44 #define MPF_NULL {{0,0,0,0}}
46 #define DEF_GFC_INTEGER_KIND(KIND,RADIX,DIGITS,BIT_SIZE) \
47 {KIND, RADIX, DIGITS, BIT_SIZE, 0, MPZ_NULL, MPZ_NULL, MPZ_NULL}
49 #define DEF_GFC_LOGICAL_KIND(KIND,BIT_SIZE) \
52 #define DEF_GFC_REAL_KIND(KIND,RADIX,DIGITS,MIN_EXP, MAX_EXP) \
53 {KIND, RADIX, DIGITS, MIN_EXP, MAX_EXP, \
54 0, 0, MPF_NULL, MPF_NULL, MPF_NULL}
56 gfc_integer_info gfc_integer_kinds[] = {
57 DEF_GFC_INTEGER_KIND (4, 2, 31, 32),
58 DEF_GFC_INTEGER_KIND (8, 2, 63, 64),
59 DEF_GFC_INTEGER_KIND (2, 2, 15, 16),
60 DEF_GFC_INTEGER_KIND (1, 2, 7, 8),
61 DEF_GFC_INTEGER_KIND (0, 0, 0, 0)
64 gfc_logical_info gfc_logical_kinds[] = {
65 DEF_GFC_LOGICAL_KIND (4, 32),
66 DEF_GFC_LOGICAL_KIND (8, 64),
67 DEF_GFC_LOGICAL_KIND (2, 16),
68 DEF_GFC_LOGICAL_KIND (1, 8),
69 DEF_GFC_LOGICAL_KIND (0, 0)
72 gfc_real_info gfc_real_kinds[] = {
73 DEF_GFC_REAL_KIND (4, 2, 24, -125, 128),
74 DEF_GFC_REAL_KIND (8, 2, 53, -1021, 1024),
75 DEF_GFC_REAL_KIND (0, 0, 0, 0, 0)
79 /* The integer kind to use for array indices. This will be set to the
80 proper value based on target information from the backend. */
82 int gfc_index_integer_kind;
85 /* Compute the natural log of arg.
87 We first get the argument into the range 0.5 to 1.5 by successive
88 multiplications or divisions by e. Then we use the series:
90 ln(x) = (x-1) - (x-1)^/2 + (x-1)^3/3 - (x-1)^4/4 + ...
92 Because we are expanding in powers of (x-1), and 0.5 < x < 1.5, we
93 have -0.5 < (x-1) < 0.5. Ignoring the harmonic term, this means
94 that each term is at most 1/(2^i), meaning one bit is gained per
97 Not very efficient, but it doesn't have to be. */
100 natural_logarithm (mpf_t * arg, mpf_t * result)
105 mpf_init_set (x, *arg);
110 mpf_set_str (t, "0.5", 10);
111 while (mpf_cmp (x, t) < 0)
117 mpf_set_str (t, "1.5", 10);
118 while (mpf_cmp (x, t) > 0)
124 mpf_sub_ui (x, x, 1);
125 mpf_init_set_ui (log, 0);
126 mpf_init_set_ui (xp, 1);
128 for (i = 1; i < GFC_REAL_BITS; i++)
131 mpf_div_ui (t, xp, i);
134 mpf_sub (log, log, t);
136 mpf_add (log, log, t);
139 /* Add in the log (e^p) = p */
142 mpf_sub_ui (log, log, -p);
144 mpf_add_ui (log, log, p);
150 mpf_set (*result, log);
155 /* Calculate the common logarithm of arg. We use the natural
156 logaritm of arg and of 10:
158 log10(arg) = log(arg)/log(10) */
161 common_logarithm (mpf_t * arg, mpf_t * result)
165 natural_logarithm (arg, result);
167 mpf_init_set_ui (i10, 10);
169 natural_logarithm (&i10, &log10);
171 mpf_div (*result, *result, log10);
176 /* Calculate exp(arg).
178 We use a reduction of the form
182 Then we obtain exp(r) from the McLaurin series.
183 exp(x) is then recovered from the identity
185 exp(x) = 2^N*exp(r). */
188 exponential (mpf_t * arg, mpf_t * result)
190 mpf_t two, ln2, power, q, r, num, denom, term, x, xp;
196 mpf_init_set (x, *arg);
198 if (mpf_cmp_ui (x, 0) == 0)
200 mpf_set_ui (*result, 1);
202 else if (mpf_cmp_ui (x, 1) == 0)
204 mpf_set (*result, e);
208 mpf_init_set_ui (two, 2);
215 natural_logarithm (&two, &ln2);
218 mpf_floor (power, q);
219 mpf_mul (q, power, ln2);
222 mpf_init_set_ui (xp, 1);
223 mpf_init_set_ui (num, 1);
224 mpf_init_set_ui (denom, 1);
226 for (i = 1; i <= GFC_REAL_BITS + 10; i++)
228 mpf_mul (num, num, r);
229 mpf_mul_ui (denom, denom, i);
230 mpf_div (term, num, denom);
231 mpf_add (xp, xp, term);
234 /* Reconstruction step */
235 n = (long) mpf_get_d (power);
239 p = (unsigned int) n;
240 mpf_mul_2exp (*result, xp, p);
244 mp = (unsigned int) (-n);
245 mpf_div_2exp (*result, xp, mp);
263 /* Calculate sin(arg).
265 We use a reduction of the form
269 Then we obtain sin(r) from the McLaurin series. */
272 sine (mpf_t * arg, mpf_t * result)
274 mpf_t factor, q, r, num, denom, term, x, xp;
277 mpf_init_set (x, *arg);
279 /* Special case (we do not treat multiples of pi due to roundoff issues) */
280 if (mpf_cmp_ui (x, 0) == 0)
282 mpf_set_ui (*result, 0);
291 mpf_div (q, x, two_pi);
292 mpf_floor (factor, q);
293 mpf_mul (q, factor, two_pi);
296 mpf_init_set_ui (xp, 0);
297 mpf_init_set_ui (num, 1);
298 mpf_init_set_ui (denom, 1);
301 for (i = 1; i < GFC_REAL_BITS + 10; i++)
303 mpf_mul (num, num, r);
304 mpf_mul_ui (denom, denom, i);
309 mpf_div (term, num, denom);
311 mpf_add (xp, xp, term);
313 mpf_sub (xp, xp, term);
316 mpf_set (*result, xp);
331 /* Calculate cos(arg).
336 cosine (mpf_t * arg, mpf_t * result)
338 mpf_t factor, q, r, num, denom, term, x, xp;
341 mpf_init_set (x, *arg);
343 /* Special case (we do not treat multiples of pi due to roundoff issues) */
344 if (mpf_cmp_ui (x, 0) == 0)
346 mpf_set_ui (*result, 1);
355 mpf_div (q, x, two_pi);
356 mpf_floor (factor, q);
357 mpf_mul (q, factor, two_pi);
360 mpf_init_set_ui (xp, 1);
361 mpf_init_set_ui (num, 1);
362 mpf_init_set_ui (denom, 1);
365 for (i = 1; i < GFC_REAL_BITS + 10; i++)
367 mpf_mul (num, num, r);
368 mpf_mul_ui (denom, denom, i);
373 mpf_div (term, num, denom);
375 mpf_add (xp, xp, term);
377 mpf_sub (xp, xp, term);
379 mpf_set (*result, xp);
394 /* Calculate atan(arg).
396 Similar to sine but requires special handling for x near 1. */
399 arctangent (mpf_t * arg, mpf_t * result)
401 mpf_t absval, convgu, convgl, num, term, x, xp;
404 mpf_init_set (x, *arg);
407 if (mpf_cmp_ui (x, 0) == 0)
409 mpf_set_ui (*result, 0);
411 else if (mpf_cmp_ui (x, 1) == 0)
414 mpf_div_ui (num, half_pi, 2);
415 mpf_set (*result, num);
418 else if (mpf_cmp_si (x, -1) == 0)
421 mpf_div_ui (num, half_pi, 2);
422 mpf_neg (*result, num);
426 { /* General cases */
431 mpf_init_set_d (convgu, 1.5);
432 mpf_init_set_d (convgl, 0.5);
433 mpf_init_set_ui (num, 1);
436 if (mpf_cmp (absval, convgl) < 0)
438 mpf_init_set_ui (xp, 0);
440 for (i = 1; i < GFC_REAL_BITS + 10; i++)
442 mpf_mul (num, num, absval);
447 mpf_div_ui (term, num, i);
449 mpf_add (xp, xp, term);
451 mpf_sub (xp, xp, term);
454 else if (mpf_cmp (absval, convgu) >= 0)
456 mpf_init_set (xp, half_pi);
458 for (i = 1; i < GFC_REAL_BITS + 10; i++)
460 mpf_div (num, num, absval);
465 mpf_div_ui (term, num, i);
467 mpf_add (xp, xp, term);
469 mpf_sub (xp, xp, term);
474 mpf_init_set_ui (xp, 0);
476 mpf_sub_ui (num, absval, 1);
477 mpf_add_ui (term, absval, 1);
478 mpf_div (absval, num, term);
483 for (i = 1; i < GFC_REAL_BITS + 10; i++)
485 mpf_mul (num, num, absval);
489 mpf_div_ui (term, num, i);
491 mpf_add (xp, xp, term);
493 mpf_sub (xp, xp, term);
496 mpf_div_ui (term, half_pi, 2);
497 mpf_add (xp, term, xp);
500 /* This makes sure to preserve the identity arctan(-x) = -arctan(x)
501 and improves accuracy to boot. */
503 if (mpf_cmp_ui (x, 0) > 0)
504 mpf_set (*result, xp);
506 mpf_neg (*result, xp);
519 /* Calculate atan2 (y, x)
521 atan2(y, x) = atan(y/x) if x > 0,
522 sign(y)*(pi - atan(|y/x|)) if x < 0,
523 0 if x = 0 && y == 0,
524 sign(y)*pi/2 if x = 0 && y != 0.
528 arctangent2 (mpf_t * y, mpf_t * x, mpf_t * result)
534 switch (mpf_sgn (*x))
538 arctangent (&t, result);
544 mpf_sub (*result, pi, t);
545 if (mpf_sgn (*y) == -1)
546 mpf_neg (*result, *result);
549 if (mpf_sgn (*y) == 0)
550 mpf_set_ui (*result, 0);
553 mpf_set (*result, half_pi);
554 if (mpf_sgn (*y) == -1)
555 mpf_neg (*result, *result);
562 /* Calculate cosh(arg). */
565 hypercos (mpf_t * arg, mpf_t * result)
567 mpf_t neg, term1, term2, x, xp;
569 mpf_init_set (x, *arg);
578 exponential (&x, &term1);
579 exponential (&neg, &term2);
581 mpf_add (xp, term1, term2);
582 mpf_div_ui (*result, xp, 2);
592 /* Calculate sinh(arg). */
595 hypersine (mpf_t * arg, mpf_t * result)
597 mpf_t neg, term1, term2, x, xp;
599 mpf_init_set (x, *arg);
608 exponential (&x, &term1);
609 exponential (&neg, &term2);
611 mpf_sub (xp, term1, term2);
612 mpf_div_ui (*result, xp, 2);
622 /* Given an arithmetic error code, return a pointer to a string that
623 explains the error. */
626 gfc_arith_error (arith code)
636 p = "Arithmetic overflow";
638 case ARITH_UNDERFLOW:
639 p = "Arithmetic underflow";
642 p = "Division by zero";
645 p = "Indeterminate form 0 ** 0";
647 case ARITH_INCOMMENSURATE:
648 p = "Array operands are incommensurate";
651 gfc_internal_error ("gfc_arith_error(): Bad error code");
658 /* Get things ready to do math. */
661 gfc_arith_init_1 (void)
663 gfc_integer_info *int_info;
664 gfc_real_info *real_info;
669 /* Set the default precision for GMP computations. */
670 mpf_set_default_prec (GFC_REAL_BITS + 30);
672 /* Calculate e, needed by the natural_logarithm() subroutine. */
674 mpf_init_set_ui (e, 0);
675 mpf_init_set_ui (a, 1);
677 for (i = 1; i < 100; i++)
680 mpf_div_ui (a, a, i); /* 1/(i!) */
683 /* Calculate pi, 2pi, pi/2, and -pi/2, needed for trigonometric
686 We use the Bailey, Borwein and Plouffe formula:
688 pi = \sum{n=0}^\infty (1/16)^n [4/(8n+1) - 2/(8n+4) - 1/(8n+5) - 1/(8n+6)]
690 which gives about four bits per iteration. */
692 mpf_init_set_ui (pi, 0);
697 limit = (GFC_REAL_BITS / 4) + 10; /* (1/16)^n gives 4 bits per iteration */
699 for (n = 0; n < limit; n++)
702 mpf_div_ui (b, b, 8 * n + 1); /* 4/(8n+1) */
705 mpf_div_ui (a, a, 8 * n + 4); /* 2/(8n+4) */
709 mpf_div_ui (a, a, 8 * n + 5); /* 1/(8n+5) */
713 mpf_div_ui (a, a, 8 * n + 6); /* 1/(8n+6) */
717 mpf_pow_ui (a, a, n); /* 16^n */
724 mpf_mul_ui (two_pi, pi, 2);
725 mpf_div_ui (half_pi, pi, 2);
727 /* Convert the minimum/maximum values for each kind into their
728 GNU MP representation. */
731 for (int_info = gfc_integer_kinds; int_info->kind != 0; int_info++)
734 mpz_set_ui (r, int_info->radix);
735 mpz_pow_ui (r, r, int_info->digits);
737 mpz_init (int_info->huge);
738 mpz_sub_ui (int_info->huge, r, 1);
740 /* These are the numbers that are actually representable by the
741 target. For bases other than two, this needs to be changed. */
742 if (int_info->radix != 2)
743 gfc_internal_error ("Fix min_int, max_int calculation");
745 mpz_init (int_info->min_int);
746 mpz_neg (int_info->min_int, int_info->huge);
747 /* No -1 here, because the representation is symmetric. */
749 mpz_init (int_info->max_int);
750 mpz_add (int_info->max_int, int_info->huge, int_info->huge);
751 mpz_add_ui (int_info->max_int, int_info->max_int, 1);
754 mpf_set_z (a, int_info->huge);
755 common_logarithm (&a, &a);
758 int_info->range = mpz_get_si (r);
761 /* mpf_set_default_prec(GFC_REAL_BITS); */
762 for (real_info = gfc_real_kinds; real_info->kind != 0; real_info++)
765 mpf_set_ui (a, real_info->radix);
766 mpf_set_ui (b, real_info->radix);
768 mpf_pow_ui (a, a, real_info->max_exponent);
769 mpf_pow_ui (b, b, real_info->max_exponent - real_info->digits);
771 mpf_init (real_info->huge);
772 mpf_sub (real_info->huge, a, b);
775 mpf_set_ui (b, real_info->radix);
776 mpf_pow_ui (b, b, 1 - real_info->min_exponent);
778 mpf_init (real_info->tiny);
779 mpf_ui_div (real_info->tiny, 1, b);
782 mpf_set_ui (b, real_info->radix);
783 mpf_pow_ui (b, b, real_info->digits - 1);
785 mpf_init (real_info->epsilon);
786 mpf_ui_div (real_info->epsilon, 1, b);
789 common_logarithm (&real_info->huge, &a);
790 common_logarithm (&real_info->tiny, &b);
793 if (mpf_cmp (a, b) > 0)
794 mpf_set (a, b); /* a = min(a, b) */
798 real_info->range = mpz_get_si (r);
801 mpf_set_ui (a, real_info->radix);
802 common_logarithm (&a, &a);
804 mpf_mul_ui (a, a, real_info->digits - 1);
807 real_info->precision = mpz_get_si (r);
809 /* If the radix is an integral power of 10, add one to the
811 for (i = 10; i <= real_info->radix; i *= 10)
812 if (i == real_info->radix)
813 real_info->precision++;
822 /* Clean up, get rid of numeric constants. */
825 gfc_arith_done_1 (void)
827 gfc_integer_info *ip;
836 for (ip = gfc_integer_kinds; ip->kind; ip++)
838 mpz_clear (ip->min_int);
839 mpz_clear (ip->max_int);
840 mpz_clear (ip->huge);
843 for (rp = gfc_real_kinds; rp->kind; rp++)
845 mpf_clear (rp->epsilon);
846 mpf_clear (rp->huge);
847 mpf_clear (rp->tiny);
852 /* Return default kinds. */
855 gfc_default_integer_kind (void)
857 return gfc_integer_kinds[gfc_option.i8 ? 1 : 0].kind;
861 gfc_default_real_kind (void)
863 return gfc_real_kinds[gfc_option.r8 ? 1 : 0].kind;
867 gfc_default_double_kind (void)
869 return gfc_real_kinds[1].kind;
873 gfc_default_character_kind (void)
879 gfc_default_logical_kind (void)
881 return gfc_logical_kinds[gfc_option.i8 ? 1 : 0].kind;
885 gfc_default_complex_kind (void)
887 return gfc_default_real_kind ();
891 /* Make sure that a valid kind is present. Returns an index into the
892 gfc_integer_kinds array, -1 if the kind is not present. */
895 validate_integer (int kind)
901 if (gfc_integer_kinds[i].kind == 0)
906 if (gfc_integer_kinds[i].kind == kind)
915 validate_real (int kind)
921 if (gfc_real_kinds[i].kind == 0)
926 if (gfc_real_kinds[i].kind == kind)
935 validate_logical (int kind)
941 if (gfc_logical_kinds[i].kind == 0)
946 if (gfc_logical_kinds[i].kind == kind)
955 validate_character (int kind)
958 if (kind == gfc_default_character_kind ())
964 /* Validate a kind given a basic type. The return value is the same
965 for the child functions, with -1 indicating nonexistence of the
969 gfc_validate_kind (bt type, int kind)
975 case BT_REAL: /* Fall through */
977 rc = validate_real (kind);
980 rc = validate_integer (kind);
983 rc = validate_logical (kind);
986 rc = validate_character (kind);
990 gfc_internal_error ("gfc_validate_kind(): Got bad type");
997 /* Given an integer and a kind, make sure that the integer lies within
998 the range of the kind. Returns ARITH_OK or ARITH_OVERFLOW. */
1001 gfc_check_integer_range (mpz_t p, int kind)
1006 i = validate_integer (kind);
1008 gfc_internal_error ("gfc_check_integer_range(): Bad kind");
1012 if (mpz_cmp (p, gfc_integer_kinds[i].min_int) < 0
1013 || mpz_cmp (p, gfc_integer_kinds[i].max_int) > 0)
1014 result = ARITH_OVERFLOW;
1020 /* Given a real and a kind, make sure that the real lies within the
1021 range of the kind. Returns ARITH_OK, ARITH_OVERFLOW or
1025 gfc_check_real_range (mpf_t p, int kind)
1034 i = validate_real (kind);
1036 gfc_internal_error ("gfc_check_real_range(): Bad kind");
1039 if (mpf_sgn (q) == 0)
1042 if (mpf_cmp (q, gfc_real_kinds[i].huge) == 1)
1044 retval = ARITH_OVERFLOW;
1048 if (mpf_cmp (q, gfc_real_kinds[i].tiny) == -1)
1049 retval = ARITH_UNDERFLOW;
1058 /* Function to return a constant expression node of a given type and
1062 gfc_constant_result (bt type, int kind, locus * where)
1068 ("gfc_constant_result(): locus 'where' cannot be NULL");
1070 result = gfc_get_expr ();
1072 result->expr_type = EXPR_CONSTANT;
1073 result->ts.type = type;
1074 result->ts.kind = kind;
1075 result->where = *where;
1080 mpz_init (result->value.integer);
1084 mpf_init (result->value.real);
1088 mpf_init (result->value.complex.r);
1089 mpf_init (result->value.complex.i);
1100 /* Low-level arithmetic functions. All of these subroutines assume
1101 that all operands are of the same type and return an operand of the
1102 same type. The other thing about these subroutines is that they
1103 can fail in various ways -- overflow, underflow, division by zero,
1104 zero raised to the zero, etc. */
1107 gfc_arith_not (gfc_expr * op1, gfc_expr ** resultp)
1111 result = gfc_constant_result (BT_LOGICAL, op1->ts.kind, &op1->where);
1112 result->value.logical = !op1->value.logical;
1120 gfc_arith_and (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1124 result = gfc_constant_result (BT_LOGICAL, gfc_kind_max (op1, op2),
1126 result->value.logical = op1->value.logical && op2->value.logical;
1134 gfc_arith_or (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1138 result = gfc_constant_result (BT_LOGICAL, gfc_kind_max (op1, op2),
1140 result->value.logical = op1->value.logical || op2->value.logical;
1148 gfc_arith_eqv (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1152 result = gfc_constant_result (BT_LOGICAL, gfc_kind_max (op1, op2),
1154 result->value.logical = op1->value.logical == op2->value.logical;
1162 gfc_arith_neqv (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1166 result = gfc_constant_result (BT_LOGICAL, gfc_kind_max (op1, op2),
1168 result->value.logical = op1->value.logical != op2->value.logical;
1175 /* Make sure a constant numeric expression is within the range for
1176 it's type and kind. Note that there's also a gfc_check_range(),
1177 but that one deals with the intrinsic RANGE function. */
1180 gfc_range_check (gfc_expr * e)
1187 rc = gfc_check_integer_range (e->value.integer, e->ts.kind);
1191 rc = gfc_check_real_range (e->value.real, e->ts.kind);
1195 rc = gfc_check_real_range (e->value.complex.r, e->ts.kind);
1197 rc = gfc_check_real_range (e->value.complex.i, e->ts.kind);
1202 gfc_internal_error ("gfc_range_check(): Bad type");
1209 /* It may seem silly to have a subroutine that actually computes the
1210 unary plus of a constant, but it prevents us from making exceptions
1211 in the code elsewhere. */
1214 gfc_arith_uplus (gfc_expr * op1, gfc_expr ** resultp)
1217 *resultp = gfc_copy_expr (op1);
1223 gfc_arith_uminus (gfc_expr * op1, gfc_expr ** resultp)
1228 result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
1230 switch (op1->ts.type)
1233 mpz_neg (result->value.integer, op1->value.integer);
1237 mpf_neg (result->value.real, op1->value.real);
1241 mpf_neg (result->value.complex.r, op1->value.complex.r);
1242 mpf_neg (result->value.complex.i, op1->value.complex.i);
1246 gfc_internal_error ("gfc_arith_uminus(): Bad basic type");
1249 rc = gfc_range_check (result);
1252 gfc_free_expr (result);
1261 gfc_arith_plus (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1266 result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
1268 switch (op1->ts.type)
1271 mpz_add (result->value.integer, op1->value.integer, op2->value.integer);
1275 mpf_add (result->value.real, op1->value.real, op2->value.real);
1279 mpf_add (result->value.complex.r, op1->value.complex.r,
1280 op2->value.complex.r);
1282 mpf_add (result->value.complex.i, op1->value.complex.i,
1283 op2->value.complex.i);
1287 gfc_internal_error ("gfc_arith_plus(): Bad basic type");
1290 rc = gfc_range_check (result);
1293 gfc_free_expr (result);
1302 gfc_arith_minus (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1307 result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
1309 switch (op1->ts.type)
1312 mpz_sub (result->value.integer, op1->value.integer, op2->value.integer);
1316 mpf_sub (result->value.real, op1->value.real, op2->value.real);
1320 mpf_sub (result->value.complex.r, op1->value.complex.r,
1321 op2->value.complex.r);
1323 mpf_sub (result->value.complex.i, op1->value.complex.i,
1324 op2->value.complex.i);
1329 gfc_internal_error ("gfc_arith_minus(): Bad basic type");
1332 rc = gfc_range_check (result);
1335 gfc_free_expr (result);
1344 gfc_arith_times (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1350 result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
1352 switch (op1->ts.type)
1355 mpz_mul (result->value.integer, op1->value.integer, op2->value.integer);
1359 mpf_mul (result->value.real, op1->value.real, op2->value.real);
1366 mpf_mul (x, op1->value.complex.r, op2->value.complex.r);
1367 mpf_mul (y, op1->value.complex.i, op2->value.complex.i);
1368 mpf_sub (result->value.complex.r, x, y);
1370 mpf_mul (x, op1->value.complex.r, op2->value.complex.i);
1371 mpf_mul (y, op1->value.complex.i, op2->value.complex.r);
1372 mpf_add (result->value.complex.i, x, y);
1380 gfc_internal_error ("gfc_arith_times(): Bad basic type");
1383 rc = gfc_range_check (result);
1386 gfc_free_expr (result);
1395 gfc_arith_divide (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1403 result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
1405 switch (op1->ts.type)
1408 if (mpz_sgn (op2->value.integer) == 0)
1414 mpz_tdiv_q (result->value.integer, op1->value.integer,
1415 op2->value.integer);
1419 if (mpf_sgn (op2->value.real) == 0)
1425 mpf_div (result->value.real, op1->value.real, op2->value.real);
1429 if (mpf_sgn (op2->value.complex.r) == 0
1430 && mpf_sgn (op2->value.complex.i) == 0)
1440 mpf_mul (x, op2->value.complex.r, op2->value.complex.r);
1441 mpf_mul (y, op2->value.complex.i, op2->value.complex.i);
1442 mpf_add (div, x, y);
1444 mpf_mul (x, op1->value.complex.r, op2->value.complex.r);
1445 mpf_mul (y, op1->value.complex.i, op2->value.complex.i);
1446 mpf_add (result->value.complex.r, x, y);
1447 mpf_div (result->value.complex.r, result->value.complex.r, div);
1449 mpf_mul (x, op1->value.complex.i, op2->value.complex.r);
1450 mpf_mul (y, op1->value.complex.r, op2->value.complex.i);
1451 mpf_sub (result->value.complex.i, x, y);
1452 mpf_div (result->value.complex.i, result->value.complex.i, div);
1461 gfc_internal_error ("gfc_arith_divide(): Bad basic type");
1465 rc = gfc_range_check (result);
1468 gfc_free_expr (result);
1476 /* Compute the reciprocal of a complex number (guaranteed nonzero). */
1479 complex_reciprocal (gfc_expr * op)
1481 mpf_t mod, a, result_r, result_i;
1486 mpf_mul (mod, op->value.complex.r, op->value.complex.r);
1487 mpf_mul (a, op->value.complex.i, op->value.complex.i);
1488 mpf_add (mod, mod, a);
1490 mpf_init (result_r);
1491 mpf_div (result_r, op->value.complex.r, mod);
1493 mpf_init (result_i);
1494 mpf_neg (result_i, op->value.complex.i);
1495 mpf_div (result_i, result_i, mod);
1497 mpf_set (op->value.complex.r, result_r);
1498 mpf_set (op->value.complex.i, result_i);
1500 mpf_clear (result_r);
1501 mpf_clear (result_i);
1508 /* Raise a complex number to positive power. */
1511 complex_pow_ui (gfc_expr * base, int power, gfc_expr * result)
1513 mpf_t temp_r, temp_i, a;
1515 mpf_set_ui (result->value.complex.r, 1);
1516 mpf_set_ui (result->value.complex.i, 0);
1522 for (; power > 0; power--)
1524 mpf_mul (temp_r, base->value.complex.r, result->value.complex.r);
1525 mpf_mul (a, base->value.complex.i, result->value.complex.i);
1526 mpf_sub (temp_r, temp_r, a);
1528 mpf_mul (temp_i, base->value.complex.r, result->value.complex.i);
1529 mpf_mul (a, base->value.complex.i, result->value.complex.r);
1530 mpf_add (temp_i, temp_i, a);
1532 mpf_set (result->value.complex.r, temp_r);
1533 mpf_set (result->value.complex.i, temp_i);
1542 /* Raise a number to an integer power. */
1545 gfc_arith_power (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1555 if (gfc_extract_int (op2, &power) != NULL)
1556 gfc_internal_error ("gfc_arith_power(): Bad exponent");
1558 result = gfc_constant_result (op1->ts.type, op1->ts.kind, &op1->where);
1561 { /* Handle something to the zeroth power */
1562 switch (op1->ts.type)
1565 if (mpz_sgn (op1->value.integer) == 0)
1568 mpz_set_ui (result->value.integer, 1);
1573 if (mpf_sgn (op1->value.real) == 0)
1576 mpf_set_ui (result->value.real, 1);
1581 if (mpf_sgn (op1->value.complex.r) == 0
1582 && mpf_sgn (op1->value.complex.i) == 0)
1586 mpf_set_ui (result->value.complex.r, 1);
1587 mpf_set_ui (result->value.complex.r, 0);
1593 gfc_internal_error ("gfc_arith_power(): Bad base");
1603 switch (op1->ts.type)
1606 mpz_pow_ui (result->value.integer, op1->value.integer, apower);
1610 mpz_init_set_ui (unity_z, 1);
1611 mpz_tdiv_q (result->value.integer, unity_z,
1612 result->value.integer);
1613 mpz_clear (unity_z);
1619 mpf_pow_ui (result->value.real, op1->value.real, apower);
1623 mpf_init_set_ui (unity_f, 1);
1624 mpf_div (result->value.real, unity_f, result->value.real);
1625 mpf_clear (unity_f);
1631 complex_pow_ui (op1, apower, result);
1633 complex_reciprocal (result);
1643 rc = gfc_range_check (result);
1646 gfc_free_expr (result);
1654 /* Concatenate two string constants. */
1657 gfc_arith_concat (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1662 result = gfc_constant_result (BT_CHARACTER, gfc_default_character_kind (),
1665 len = op1->value.character.length + op2->value.character.length;
1667 result->value.character.string = gfc_getmem (len + 1);
1668 result->value.character.length = len;
1670 memcpy (result->value.character.string, op1->value.character.string,
1671 op1->value.character.length);
1673 memcpy (result->value.character.string + op1->value.character.length,
1674 op2->value.character.string, op2->value.character.length);
1676 result->value.character.string[len] = '\0';
1684 /* Comparison operators. Assumes that the two expression nodes
1685 contain two constants of the same type. */
1688 gfc_compare_expr (gfc_expr * op1, gfc_expr * op2)
1692 switch (op1->ts.type)
1695 rc = mpz_cmp (op1->value.integer, op2->value.integer);
1699 rc = mpf_cmp (op1->value.real, op2->value.real);
1703 rc = gfc_compare_string (op1, op2, NULL);
1707 rc = ((!op1->value.logical && op2->value.logical)
1708 || (op1->value.logical && !op2->value.logical));
1712 gfc_internal_error ("gfc_compare_expr(): Bad basic type");
1719 /* Compare a pair of complex numbers. Naturally, this is only for
1720 equality/nonequality. */
1723 compare_complex (gfc_expr * op1, gfc_expr * op2)
1726 return (mpf_cmp (op1->value.complex.r, op2->value.complex.r) == 0
1727 && mpf_cmp (op1->value.complex.i, op2->value.complex.i) == 0);
1731 /* Given two constant strings and the inverse collating sequence,
1732 compare the strings. We return -1 for a<b, 0 for a==b and 1 for
1733 a>b. If the xcoll_table is NULL, we use the processor's default
1734 collating sequence. */
1737 gfc_compare_string (gfc_expr * a, gfc_expr * b, const int *xcoll_table)
1739 int len, alen, blen, i, ac, bc;
1741 alen = a->value.character.length;
1742 blen = b->value.character.length;
1744 len = (alen > blen) ? alen : blen;
1746 for (i = 0; i < len; i++)
1748 ac = (i < alen) ? a->value.character.string[i] : ' ';
1749 bc = (i < blen) ? b->value.character.string[i] : ' ';
1751 if (xcoll_table != NULL)
1753 ac = xcoll_table[ac];
1754 bc = xcoll_table[bc];
1763 /* Strings are equal */
1769 /* Specific comparison subroutines. */
1772 gfc_arith_eq (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1776 result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind (),
1778 result->value.logical = (op1->ts.type == BT_COMPLEX) ?
1779 compare_complex (op1, op2) : (gfc_compare_expr (op1, op2) == 0);
1787 gfc_arith_ne (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1791 result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind (),
1793 result->value.logical = (op1->ts.type == BT_COMPLEX) ?
1794 !compare_complex (op1, op2) : (gfc_compare_expr (op1, op2) != 0);
1802 gfc_arith_gt (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1806 result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind (),
1808 result->value.logical = (gfc_compare_expr (op1, op2) > 0);
1816 gfc_arith_ge (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1820 result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind (),
1822 result->value.logical = (gfc_compare_expr (op1, op2) >= 0);
1830 gfc_arith_lt (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1834 result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind (),
1836 result->value.logical = (gfc_compare_expr (op1, op2) < 0);
1844 gfc_arith_le (gfc_expr * op1, gfc_expr * op2, gfc_expr ** resultp)
1848 result = gfc_constant_result (BT_LOGICAL, gfc_default_logical_kind (),
1850 result->value.logical = (gfc_compare_expr (op1, op2) <= 0);
1858 reduce_unary (arith (*eval) (gfc_expr *, gfc_expr **), gfc_expr * op,
1861 gfc_constructor *c, *head;
1865 if (op->expr_type == EXPR_CONSTANT)
1866 return eval (op, result);
1869 head = gfc_copy_constructor (op->value.constructor);
1871 for (c = head; c; c = c->next)
1873 rc = eval (c->expr, &r);
1877 gfc_replace_expr (c->expr, r);
1881 gfc_free_constructor (head);
1884 r = gfc_get_expr ();
1885 r->expr_type = EXPR_ARRAY;
1886 r->value.constructor = head;
1887 r->shape = gfc_copy_shape (op->shape, op->rank);
1889 r->ts = head->expr->ts;
1890 r->where = op->where;
1901 reduce_binary_ac (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1902 gfc_expr * op1, gfc_expr * op2,
1905 gfc_constructor *c, *head;
1909 head = gfc_copy_constructor (op1->value.constructor);
1912 for (c = head; c; c = c->next)
1914 rc = eval (c->expr, op2, &r);
1918 gfc_replace_expr (c->expr, r);
1922 gfc_free_constructor (head);
1925 r = gfc_get_expr ();
1926 r->expr_type = EXPR_ARRAY;
1927 r->value.constructor = head;
1928 r->shape = gfc_copy_shape (op1->shape, op1->rank);
1930 r->ts = head->expr->ts;
1931 r->where = op1->where;
1932 r->rank = op1->rank;
1942 reduce_binary_ca (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1943 gfc_expr * op1, gfc_expr * op2,
1946 gfc_constructor *c, *head;
1950 head = gfc_copy_constructor (op2->value.constructor);
1953 for (c = head; c; c = c->next)
1955 rc = eval (op1, c->expr, &r);
1959 gfc_replace_expr (c->expr, r);
1963 gfc_free_constructor (head);
1966 r = gfc_get_expr ();
1967 r->expr_type = EXPR_ARRAY;
1968 r->value.constructor = head;
1969 r->shape = gfc_copy_shape (op2->shape, op2->rank);
1971 r->ts = head->expr->ts;
1972 r->where = op2->where;
1973 r->rank = op2->rank;
1983 reduce_binary_aa (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
1984 gfc_expr * op1, gfc_expr * op2,
1987 gfc_constructor *c, *d, *head;
1991 head = gfc_copy_constructor (op1->value.constructor);
1994 d = op2->value.constructor;
1996 if (gfc_check_conformance ("Elemental binary operation", op1, op2)
1998 rc = ARITH_INCOMMENSURATE;
2002 for (c = head; c; c = c->next, d = d->next)
2006 rc = ARITH_INCOMMENSURATE;
2010 rc = eval (c->expr, d->expr, &r);
2014 gfc_replace_expr (c->expr, r);
2018 rc = ARITH_INCOMMENSURATE;
2022 gfc_free_constructor (head);
2025 r = gfc_get_expr ();
2026 r->expr_type = EXPR_ARRAY;
2027 r->value.constructor = head;
2028 r->shape = gfc_copy_shape (op1->shape, op1->rank);
2030 r->ts = head->expr->ts;
2031 r->where = op1->where;
2032 r->rank = op1->rank;
2042 reduce_binary (arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
2043 gfc_expr * op1, gfc_expr * op2,
2047 if (op1->expr_type == EXPR_CONSTANT && op2->expr_type == EXPR_CONSTANT)
2048 return eval (op1, op2, result);
2050 if (op1->expr_type == EXPR_CONSTANT && op2->expr_type == EXPR_ARRAY)
2051 return reduce_binary_ca (eval, op1, op2, result);
2053 if (op1->expr_type == EXPR_ARRAY && op2->expr_type == EXPR_CONSTANT)
2054 return reduce_binary_ac (eval, op1, op2, result);
2056 return reduce_binary_aa (eval, op1, op2, result);
2062 arith (*f2)(gfc_expr *, gfc_expr **);
2063 arith (*f3)(gfc_expr *, gfc_expr *, gfc_expr **);
2067 /* High level arithmetic subroutines. These subroutines go into
2068 eval_intrinsic(), which can do one of several things to its
2069 operands. If the operands are incompatible with the intrinsic
2070 operation, we return a node pointing to the operands and hope that
2071 an operator interface is found during resolution.
2073 If the operands are compatible and are constants, then we try doing
2074 the arithmetic. We also handle the cases where either or both
2075 operands are array constructors. */
2078 eval_intrinsic (gfc_intrinsic_op operator,
2079 eval_f eval, gfc_expr * op1, gfc_expr * op2)
2081 gfc_expr temp, *result;
2085 gfc_clear_ts (&temp.ts);
2089 case INTRINSIC_NOT: /* Logical unary */
2090 if (op1->ts.type != BT_LOGICAL)
2093 temp.ts.type = BT_LOGICAL;
2094 temp.ts.kind = gfc_default_logical_kind ();
2099 /* Logical binary operators */
2102 case INTRINSIC_NEQV:
2104 if (op1->ts.type != BT_LOGICAL || op2->ts.type != BT_LOGICAL)
2107 temp.ts.type = BT_LOGICAL;
2108 temp.ts.kind = gfc_default_logical_kind ();
2113 case INTRINSIC_UPLUS:
2114 case INTRINSIC_UMINUS: /* Numeric unary */
2115 if (!gfc_numeric_ts (&op1->ts))
2124 case INTRINSIC_LT: /* Additional restrictions */
2125 case INTRINSIC_LE: /* for ordering relations. */
2127 if (op1->ts.type == BT_COMPLEX || op2->ts.type == BT_COMPLEX)
2129 temp.ts.type = BT_LOGICAL;
2130 temp.ts.kind = gfc_default_logical_kind();
2134 /* else fall through */
2138 if (op1->ts.type == BT_CHARACTER && op2->ts.type == BT_CHARACTER)
2141 temp.ts.type = BT_LOGICAL;
2142 temp.ts.kind = gfc_default_logical_kind();
2146 /* else fall through */
2148 case INTRINSIC_PLUS:
2149 case INTRINSIC_MINUS:
2150 case INTRINSIC_TIMES:
2151 case INTRINSIC_DIVIDE:
2152 case INTRINSIC_POWER: /* Numeric binary */
2153 if (!gfc_numeric_ts (&op1->ts) || !gfc_numeric_ts (&op2->ts))
2156 /* Insert any necessary type conversions to make the operands compatible. */
2158 temp.expr_type = EXPR_OP;
2159 gfc_clear_ts (&temp.ts);
2160 temp.operator = operator;
2165 gfc_type_convert_binary (&temp);
2167 if (operator == INTRINSIC_EQ || operator == INTRINSIC_NE
2168 || operator == INTRINSIC_GE || operator == INTRINSIC_GT
2169 || operator == INTRINSIC_LE || operator == INTRINSIC_LT)
2171 temp.ts.type = BT_LOGICAL;
2172 temp.ts.kind = gfc_default_logical_kind ();
2178 case INTRINSIC_CONCAT: /* Character binary */
2179 if (op1->ts.type != BT_CHARACTER || op2->ts.type != BT_CHARACTER)
2182 temp.ts.type = BT_CHARACTER;
2183 temp.ts.kind = gfc_default_character_kind ();
2188 case INTRINSIC_USER:
2192 gfc_internal_error ("eval_intrinsic(): Bad operator");
2195 /* Try to combine the operators. */
2196 if (operator == INTRINSIC_POWER && op2->ts.type != BT_INTEGER)
2199 if (op1->expr_type != EXPR_CONSTANT
2200 && (op1->expr_type != EXPR_ARRAY
2201 || !gfc_is_constant_expr (op1)
2202 || !gfc_expanded_ac (op1)))
2206 && op2->expr_type != EXPR_CONSTANT
2207 && (op2->expr_type != EXPR_ARRAY
2208 || !gfc_is_constant_expr (op2)
2209 || !gfc_expanded_ac (op2)))
2213 rc = reduce_unary (eval.f2, op1, &result);
2215 rc = reduce_binary (eval.f3, op1, op2, &result);
2218 { /* Something went wrong */
2219 gfc_error ("%s at %L", gfc_arith_error (rc), &op1->where);
2223 gfc_free_expr (op1);
2224 gfc_free_expr (op2);
2228 /* Create a run-time expression */
2229 result = gfc_get_expr ();
2230 result->ts = temp.ts;
2232 result->expr_type = EXPR_OP;
2233 result->operator = operator;
2238 result->where = op1->where;
2244 /* Modify type of expression for zero size array. */
2246 eval_type_intrinsic0 (gfc_intrinsic_op operator, gfc_expr *op)
2249 gfc_internal_error("eval_type_intrinsic0(): op NULL");
2259 op->ts.type = BT_LOGICAL;
2260 op->ts.kind = gfc_default_logical_kind();
2271 /* Return nonzero if the expression is a zero size array. */
2274 gfc_zero_size_array (gfc_expr * e)
2277 if (e->expr_type != EXPR_ARRAY)
2280 return e->value.constructor == NULL;
2284 /* Reduce a binary expression where at least one of the operands
2285 involves a zero-length array. Returns NULL if neither of the
2286 operands is a zero-length array. */
2289 reduce_binary0 (gfc_expr * op1, gfc_expr * op2)
2292 if (gfc_zero_size_array (op1))
2294 gfc_free_expr (op2);
2298 if (gfc_zero_size_array (op2))
2300 gfc_free_expr (op1);
2309 eval_intrinsic_f2 (gfc_intrinsic_op operator,
2310 arith (*eval) (gfc_expr *, gfc_expr **),
2311 gfc_expr * op1, gfc_expr * op2)
2318 if (gfc_zero_size_array (op1))
2319 return eval_type_intrinsic0(operator, op1);
2323 result = reduce_binary0 (op1, op2);
2325 return eval_type_intrinsic0(operator, result);
2329 return eval_intrinsic (operator, f, op1, op2);
2334 eval_intrinsic_f3 (gfc_intrinsic_op operator,
2335 arith (*eval) (gfc_expr *, gfc_expr *, gfc_expr **),
2336 gfc_expr * op1, gfc_expr * op2)
2341 result = reduce_binary0 (op1, op2);
2343 return eval_type_intrinsic0(operator, result);
2346 return eval_intrinsic (operator, f, op1, op2);
2352 gfc_uplus (gfc_expr * op)
2354 return eval_intrinsic_f2 (INTRINSIC_UPLUS, gfc_arith_uplus, op, NULL);
2358 gfc_uminus (gfc_expr * op)
2360 return eval_intrinsic_f2 (INTRINSIC_UMINUS, gfc_arith_uminus, op, NULL);
2364 gfc_add (gfc_expr * op1, gfc_expr * op2)
2366 return eval_intrinsic_f3 (INTRINSIC_PLUS, gfc_arith_plus, op1, op2);
2370 gfc_subtract (gfc_expr * op1, gfc_expr * op2)
2372 return eval_intrinsic_f3 (INTRINSIC_MINUS, gfc_arith_minus, op1, op2);
2376 gfc_multiply (gfc_expr * op1, gfc_expr * op2)
2378 return eval_intrinsic_f3 (INTRINSIC_TIMES, gfc_arith_times, op1, op2);
2382 gfc_divide (gfc_expr * op1, gfc_expr * op2)
2384 return eval_intrinsic_f3 (INTRINSIC_DIVIDE, gfc_arith_divide, op1, op2);
2388 gfc_power (gfc_expr * op1, gfc_expr * op2)
2390 return eval_intrinsic_f3 (INTRINSIC_POWER, gfc_arith_power, op1, op2);
2394 gfc_concat (gfc_expr * op1, gfc_expr * op2)
2396 return eval_intrinsic_f3 (INTRINSIC_CONCAT, gfc_arith_concat, op1, op2);
2400 gfc_and (gfc_expr * op1, gfc_expr * op2)
2402 return eval_intrinsic_f3 (INTRINSIC_AND, gfc_arith_and, op1, op2);
2406 gfc_or (gfc_expr * op1, gfc_expr * op2)
2408 return eval_intrinsic_f3 (INTRINSIC_OR, gfc_arith_or, op1, op2);
2412 gfc_not (gfc_expr * op1)
2414 return eval_intrinsic_f2 (INTRINSIC_NOT, gfc_arith_not, op1, NULL);
2418 gfc_eqv (gfc_expr * op1, gfc_expr * op2)
2420 return eval_intrinsic_f3 (INTRINSIC_EQV, gfc_arith_eqv, op1, op2);
2424 gfc_neqv (gfc_expr * op1, gfc_expr * op2)
2426 return eval_intrinsic_f3 (INTRINSIC_NEQV, gfc_arith_neqv, op1, op2);
2430 gfc_eq (gfc_expr * op1, gfc_expr * op2)
2432 return eval_intrinsic_f3 (INTRINSIC_EQ, gfc_arith_eq, op1, op2);
2436 gfc_ne (gfc_expr * op1, gfc_expr * op2)
2438 return eval_intrinsic_f3 (INTRINSIC_NE, gfc_arith_ne, op1, op2);
2442 gfc_gt (gfc_expr * op1, gfc_expr * op2)
2444 return eval_intrinsic_f3 (INTRINSIC_GT, gfc_arith_gt, op1, op2);
2448 gfc_ge (gfc_expr * op1, gfc_expr * op2)
2450 return eval_intrinsic_f3 (INTRINSIC_GE, gfc_arith_ge, op1, op2);
2454 gfc_lt (gfc_expr * op1, gfc_expr * op2)
2456 return eval_intrinsic_f3 (INTRINSIC_LT, gfc_arith_lt, op1, op2);
2460 gfc_le (gfc_expr * op1, gfc_expr * op2)
2462 return eval_intrinsic_f3 (INTRINSIC_LE, gfc_arith_le, op1, op2);
2466 /* Convert an integer string to an expression node. */
2469 gfc_convert_integer (const char *buffer, int kind, int radix, locus * where)
2474 e = gfc_constant_result (BT_INTEGER, kind, where);
2475 /* a leading plus is allowed, but not by mpz_set_str */
2476 if (buffer[0] == '+')
2480 mpz_set_str (e->value.integer, t, radix);
2486 /* Convert a real string to an expression node. */
2489 gfc_convert_real (const char *buffer, int kind, locus * where)
2494 e = gfc_constant_result (BT_REAL, kind, where);
2495 /* a leading plus is allowed, but not by mpf_set_str */
2496 if (buffer[0] == '+')
2500 mpf_set_str (e->value.real, t, 10);
2506 /* Convert a pair of real, constant expression nodes to a single
2507 complex expression node. */
2510 gfc_convert_complex (gfc_expr * real, gfc_expr * imag, int kind)
2514 e = gfc_constant_result (BT_COMPLEX, kind, &real->where);
2515 mpf_set (e->value.complex.r, real->value.real);
2516 mpf_set (e->value.complex.i, imag->value.real);
2522 /******* Simplification of intrinsic functions with constant arguments *****/
2525 /* Deal with an arithmetic error. */
2528 arith_error (arith rc, gfc_typespec * from, gfc_typespec * to, locus * where)
2531 gfc_error ("%s converting %s to %s at %L", gfc_arith_error (rc),
2532 gfc_typename (from), gfc_typename (to), where);
2534 /* TODO: Do something about the error, ie underflow rounds to 0,
2535 throw exception, return NaN, etc. */
2538 /* Convert integers to integers. */
2541 gfc_int2int (gfc_expr * src, int kind)
2546 result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2548 mpz_set (result->value.integer, src->value.integer);
2550 if ((rc = gfc_check_integer_range (result->value.integer, kind))
2553 arith_error (rc, &src->ts, &result->ts, &src->where);
2554 gfc_free_expr (result);
2562 /* Convert integers to reals. */
2565 gfc_int2real (gfc_expr * src, int kind)
2570 result = gfc_constant_result (BT_REAL, kind, &src->where);
2572 mpf_set_z (result->value.real, src->value.integer);
2574 if ((rc = gfc_check_real_range (result->value.real, kind)) != ARITH_OK)
2576 arith_error (rc, &src->ts, &result->ts, &src->where);
2577 gfc_free_expr (result);
2585 /* Convert default integer to default complex. */
2588 gfc_int2complex (gfc_expr * src, int kind)
2593 result = gfc_constant_result (BT_COMPLEX, kind, &src->where);
2595 mpf_set_z (result->value.complex.r, src->value.integer);
2596 mpf_set_ui (result->value.complex.i, 0);
2598 if ((rc = gfc_check_real_range (result->value.complex.i, kind)) != ARITH_OK)
2600 arith_error (rc, &src->ts, &result->ts, &src->where);
2601 gfc_free_expr (result);
2609 /* Convert default real to default integer. */
2612 gfc_real2int (gfc_expr * src, int kind)
2617 result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2619 mpz_set_f (result->value.integer, src->value.real);
2621 if ((rc = gfc_check_integer_range (result->value.integer, kind))
2624 arith_error (rc, &src->ts, &result->ts, &src->where);
2625 gfc_free_expr (result);
2633 /* Convert real to real. */
2636 gfc_real2real (gfc_expr * src, int kind)
2641 result = gfc_constant_result (BT_REAL, kind, &src->where);
2643 mpf_set (result->value.real, src->value.real);
2645 if ((rc = gfc_check_real_range (result->value.real, kind)) != ARITH_OK)
2647 arith_error (rc, &src->ts, &result->ts, &src->where);
2648 gfc_free_expr (result);
2656 /* Convert real to complex. */
2659 gfc_real2complex (gfc_expr * src, int kind)
2664 result = gfc_constant_result (BT_COMPLEX, kind, &src->where);
2666 mpf_set (result->value.complex.r, src->value.real);
2667 mpf_set_ui (result->value.complex.i, 0);
2669 if ((rc = gfc_check_real_range (result->value.complex.i, kind)) != ARITH_OK)
2671 arith_error (rc, &src->ts, &result->ts, &src->where);
2672 gfc_free_expr (result);
2680 /* Convert complex to integer. */
2683 gfc_complex2int (gfc_expr * src, int kind)
2688 result = gfc_constant_result (BT_INTEGER, kind, &src->where);
2690 mpz_set_f (result->value.integer, src->value.complex.r);
2692 if ((rc = gfc_check_integer_range (result->value.integer, kind))
2695 arith_error (rc, &src->ts, &result->ts, &src->where);
2696 gfc_free_expr (result);
2704 /* Convert complex to real. */
2707 gfc_complex2real (gfc_expr * src, int kind)
2712 result = gfc_constant_result (BT_REAL, kind, &src->where);
2714 mpf_set (result->value.real, src->value.complex.r);
2716 if ((rc = gfc_check_real_range (result->value.real, kind)) != ARITH_OK)
2718 arith_error (rc, &src->ts, &result->ts, &src->where);
2719 gfc_free_expr (result);
2727 /* Convert complex to complex. */
2730 gfc_complex2complex (gfc_expr * src, int kind)
2735 result = gfc_constant_result (BT_COMPLEX, kind, &src->where);
2737 mpf_set (result->value.complex.r, src->value.complex.r);
2738 mpf_set (result->value.complex.i, src->value.complex.i);
2740 if ((rc = gfc_check_real_range (result->value.complex.r, kind)) != ARITH_OK
2742 gfc_check_real_range (result->value.complex.i, kind)) != ARITH_OK)
2744 arith_error (rc, &src->ts, &result->ts, &src->where);
2745 gfc_free_expr (result);
2753 /* Logical kind conversion. */
2756 gfc_log2log (gfc_expr * src, int kind)
2760 result = gfc_constant_result (BT_LOGICAL, kind, &src->where);
2761 result->value.logical = src->value.logical;