2 Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008
3 Free Software Foundation, Inc.
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 3, 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 COPYING3. If not see
20 <http://www.gnu.org/licenses/>. */
22 /* Since target arithmetic must be done on the host, there has to
23 be some way of evaluating arithmetic expressions as the host
24 would evaluate them. We use the GNU MP library and the MPFR
25 library to do arithmetic, and this file provides the interface. */
32 #include "target-memory.h"
34 /* MPFR does not have a direct replacement for mpz_set_f() from GMP.
35 It's easily implemented with a few calls though. */
38 gfc_mpfr_to_mpz (mpz_t z
, mpfr_t x
)
42 e
= mpfr_get_z_exp (z
, x
);
43 /* MPFR 2.0.1 (included with GMP 4.1) has a bug whereby mpfr_get_z_exp
44 may set the sign of z incorrectly. Work around that here. */
45 if (mpfr_sgn (x
) != mpz_sgn (z
))
49 mpz_mul_2exp (z
, z
, e
);
51 mpz_tdiv_q_2exp (z
, z
, -e
);
55 /* Set the model number precision by the requested KIND. */
58 gfc_set_model_kind (int kind
)
60 int index
= gfc_validate_kind (BT_REAL
, kind
, false);
63 base2prec
= gfc_real_kinds
[index
].digits
;
64 if (gfc_real_kinds
[index
].radix
!= 2)
65 base2prec
*= gfc_real_kinds
[index
].radix
/ 2;
66 mpfr_set_default_prec (base2prec
);
70 /* Set the model number precision from mpfr_t x. */
73 gfc_set_model (mpfr_t x
)
75 mpfr_set_default_prec (mpfr_get_prec (x
));
79 /* Given an arithmetic error code, return a pointer to a string that
80 explains the error. */
83 gfc_arith_error (arith code
)
90 p
= _("Arithmetic OK at %L");
93 p
= _("Arithmetic overflow at %L");
96 p
= _("Arithmetic underflow at %L");
99 p
= _("Arithmetic NaN at %L");
102 p
= _("Division by zero at %L");
104 case ARITH_INCOMMENSURATE
:
105 p
= _("Array operands are incommensurate at %L");
107 case ARITH_ASYMMETRIC
:
109 _("Integer outside symmetric range implied by Standard Fortran at %L");
112 gfc_internal_error ("gfc_arith_error(): Bad error code");
119 /* Get things ready to do math. */
122 gfc_arith_init_1 (void)
124 gfc_integer_info
*int_info
;
125 gfc_real_info
*real_info
;
130 mpfr_set_default_prec (128);
134 /* Convert the minimum and maximum values for each kind into their
135 GNU MP representation. */
136 for (int_info
= gfc_integer_kinds
; int_info
->kind
!= 0; int_info
++)
139 mpz_set_ui (r
, int_info
->radix
);
140 mpz_pow_ui (r
, r
, int_info
->digits
);
142 mpz_init (int_info
->huge
);
143 mpz_sub_ui (int_info
->huge
, r
, 1);
145 /* These are the numbers that are actually representable by the
146 target. For bases other than two, this needs to be changed. */
147 if (int_info
->radix
!= 2)
148 gfc_internal_error ("Fix min_int calculation");
150 /* See PRs 13490 and 17912, related to integer ranges.
151 The pedantic_min_int exists for range checking when a program
152 is compiled with -pedantic, and reflects the belief that
153 Standard Fortran requires integers to be symmetrical, i.e.
154 every negative integer must have a representable positive
155 absolute value, and vice versa. */
157 mpz_init (int_info
->pedantic_min_int
);
158 mpz_neg (int_info
->pedantic_min_int
, int_info
->huge
);
160 mpz_init (int_info
->min_int
);
161 mpz_sub_ui (int_info
->min_int
, int_info
->pedantic_min_int
, 1);
164 mpfr_set_z (a
, int_info
->huge
, GFC_RND_MODE
);
165 mpfr_log10 (a
, a
, GFC_RND_MODE
);
167 gfc_mpfr_to_mpz (r
, a
);
168 int_info
->range
= mpz_get_si (r
);
173 for (real_info
= gfc_real_kinds
; real_info
->kind
!= 0; real_info
++)
175 gfc_set_model_kind (real_info
->kind
);
181 /* huge(x) = (1 - b**(-p)) * b**(emax-1) * b */
182 /* a = 1 - b**(-p) */
183 mpfr_set_ui (a
, 1, GFC_RND_MODE
);
184 mpfr_set_ui (b
, real_info
->radix
, GFC_RND_MODE
);
185 mpfr_pow_si (b
, b
, -real_info
->digits
, GFC_RND_MODE
);
186 mpfr_sub (a
, a
, b
, GFC_RND_MODE
);
188 /* c = b**(emax-1) */
189 mpfr_set_ui (b
, real_info
->radix
, GFC_RND_MODE
);
190 mpfr_pow_ui (c
, b
, real_info
->max_exponent
- 1, GFC_RND_MODE
);
192 /* a = a * c = (1 - b**(-p)) * b**(emax-1) */
193 mpfr_mul (a
, a
, c
, GFC_RND_MODE
);
195 /* a = (1 - b**(-p)) * b**(emax-1) * b */
196 mpfr_mul_ui (a
, a
, real_info
->radix
, GFC_RND_MODE
);
198 mpfr_init (real_info
->huge
);
199 mpfr_set (real_info
->huge
, a
, GFC_RND_MODE
);
201 /* tiny(x) = b**(emin-1) */
202 mpfr_set_ui (b
, real_info
->radix
, GFC_RND_MODE
);
203 mpfr_pow_si (b
, b
, real_info
->min_exponent
- 1, GFC_RND_MODE
);
205 mpfr_init (real_info
->tiny
);
206 mpfr_set (real_info
->tiny
, b
, GFC_RND_MODE
);
208 /* subnormal (x) = b**(emin - digit) */
209 mpfr_set_ui (b
, real_info
->radix
, GFC_RND_MODE
);
210 mpfr_pow_si (b
, b
, real_info
->min_exponent
- real_info
->digits
,
213 mpfr_init (real_info
->subnormal
);
214 mpfr_set (real_info
->subnormal
, b
, GFC_RND_MODE
);
216 /* epsilon(x) = b**(1-p) */
217 mpfr_set_ui (b
, real_info
->radix
, GFC_RND_MODE
);
218 mpfr_pow_si (b
, b
, 1 - real_info
->digits
, GFC_RND_MODE
);
220 mpfr_init (real_info
->epsilon
);
221 mpfr_set (real_info
->epsilon
, b
, GFC_RND_MODE
);
223 /* range(x) = int(min(log10(huge(x)), -log10(tiny)) */
224 mpfr_log10 (a
, real_info
->huge
, GFC_RND_MODE
);
225 mpfr_log10 (b
, real_info
->tiny
, GFC_RND_MODE
);
226 mpfr_neg (b
, b
, GFC_RND_MODE
);
229 mpfr_min (a
, a
, b
, GFC_RND_MODE
);
232 gfc_mpfr_to_mpz (r
, a
);
233 real_info
->range
= mpz_get_si (r
);
235 /* precision(x) = int((p - 1) * log10(b)) + k */
236 mpfr_set_ui (a
, real_info
->radix
, GFC_RND_MODE
);
237 mpfr_log10 (a
, a
, GFC_RND_MODE
);
239 mpfr_mul_ui (a
, a
, real_info
->digits
- 1, GFC_RND_MODE
);
241 gfc_mpfr_to_mpz (r
, a
);
242 real_info
->precision
= mpz_get_si (r
);
244 /* If the radix is an integral power of 10, add one to the precision. */
245 for (i
= 10; i
<= real_info
->radix
; i
*= 10)
246 if (i
== real_info
->radix
)
247 real_info
->precision
++;
258 /* Clean up, get rid of numeric constants. */
261 gfc_arith_done_1 (void)
263 gfc_integer_info
*ip
;
266 for (ip
= gfc_integer_kinds
; ip
->kind
; ip
++)
268 mpz_clear (ip
->min_int
);
269 mpz_clear (ip
->pedantic_min_int
);
270 mpz_clear (ip
->huge
);
273 for (rp
= gfc_real_kinds
; rp
->kind
; rp
++)
275 mpfr_clear (rp
->epsilon
);
276 mpfr_clear (rp
->huge
);
277 mpfr_clear (rp
->tiny
);
278 mpfr_clear (rp
->subnormal
);
283 /* Given an integer and a kind, make sure that the integer lies within
284 the range of the kind. Returns ARITH_OK, ARITH_ASYMMETRIC or
288 gfc_check_integer_range (mpz_t p
, int kind
)
293 i
= gfc_validate_kind (BT_INTEGER
, kind
, false);
298 if (mpz_cmp (p
, gfc_integer_kinds
[i
].pedantic_min_int
) < 0)
299 result
= ARITH_ASYMMETRIC
;
303 if (gfc_option
.flag_range_check
== 0)
306 if (mpz_cmp (p
, gfc_integer_kinds
[i
].min_int
) < 0
307 || mpz_cmp (p
, gfc_integer_kinds
[i
].huge
) > 0)
308 result
= ARITH_OVERFLOW
;
314 /* Given a real and a kind, make sure that the real lies within the
315 range of the kind. Returns ARITH_OK, ARITH_OVERFLOW or
319 gfc_check_real_range (mpfr_t p
, int kind
)
325 i
= gfc_validate_kind (BT_REAL
, kind
, false);
329 mpfr_abs (q
, p
, GFC_RND_MODE
);
333 if (gfc_option
.flag_range_check
== 0)
336 retval
= ARITH_OVERFLOW
;
338 else if (mpfr_nan_p (p
))
340 if (gfc_option
.flag_range_check
== 0)
345 else if (mpfr_sgn (q
) == 0)
347 else if (mpfr_cmp (q
, gfc_real_kinds
[i
].huge
) > 0)
349 if (gfc_option
.flag_range_check
== 0)
351 mpfr_set_inf (p
, mpfr_sgn (p
));
355 retval
= ARITH_OVERFLOW
;
357 else if (mpfr_cmp (q
, gfc_real_kinds
[i
].subnormal
) < 0)
359 if (gfc_option
.flag_range_check
== 0)
361 if (mpfr_sgn (p
) < 0)
363 mpfr_set_ui (p
, 0, GFC_RND_MODE
);
364 mpfr_set_si (q
, -1, GFC_RND_MODE
);
365 mpfr_copysign (p
, p
, q
, GFC_RND_MODE
);
368 mpfr_set_ui (p
, 0, GFC_RND_MODE
);
372 retval
= ARITH_UNDERFLOW
;
374 else if (mpfr_cmp (q
, gfc_real_kinds
[i
].tiny
) < 0)
379 /* Save current values of emin and emax. */
380 emin
= mpfr_get_emin ();
381 emax
= mpfr_get_emax ();
383 /* Set emin and emax for the current model number. */
384 en
= gfc_real_kinds
[i
].min_exponent
- gfc_real_kinds
[i
].digits
+ 1;
385 mpfr_set_emin ((mp_exp_t
) en
);
386 mpfr_set_emax ((mp_exp_t
) gfc_real_kinds
[i
].max_exponent
);
387 mpfr_subnormalize (q
, 0, GFC_RND_MODE
);
389 /* Reset emin and emax. */
390 mpfr_set_emin (emin
);
391 mpfr_set_emax (emax
);
393 /* Copy sign if needed. */
394 if (mpfr_sgn (p
) < 0)
395 mpfr_neg (p
, q
, GMP_RNDN
);
397 mpfr_set (p
, q
, GMP_RNDN
);
410 /* Function to return a constant expression node of a given type and kind. */
413 gfc_constant_result (bt type
, int kind
, locus
*where
)
418 gfc_internal_error ("gfc_constant_result(): locus 'where' cannot be NULL");
420 result
= gfc_get_expr ();
422 result
->expr_type
= EXPR_CONSTANT
;
423 result
->ts
.type
= type
;
424 result
->ts
.kind
= kind
;
425 result
->where
= *where
;
430 mpz_init (result
->value
.integer
);
434 gfc_set_model_kind (kind
);
435 mpfr_init (result
->value
.real
);
439 gfc_set_model_kind (kind
);
440 mpfr_init (result
->value
.complex.r
);
441 mpfr_init (result
->value
.complex.i
);
452 /* Low-level arithmetic functions. All of these subroutines assume
453 that all operands are of the same type and return an operand of the
454 same type. The other thing about these subroutines is that they
455 can fail in various ways -- overflow, underflow, division by zero,
456 zero raised to the zero, etc. */
459 gfc_arith_not (gfc_expr
*op1
, gfc_expr
**resultp
)
463 result
= gfc_constant_result (BT_LOGICAL
, op1
->ts
.kind
, &op1
->where
);
464 result
->value
.logical
= !op1
->value
.logical
;
472 gfc_arith_and (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
476 result
= gfc_constant_result (BT_LOGICAL
, gfc_kind_max (op1
, op2
),
478 result
->value
.logical
= op1
->value
.logical
&& op2
->value
.logical
;
486 gfc_arith_or (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
490 result
= gfc_constant_result (BT_LOGICAL
, gfc_kind_max (op1
, op2
),
492 result
->value
.logical
= op1
->value
.logical
|| op2
->value
.logical
;
500 gfc_arith_eqv (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
504 result
= gfc_constant_result (BT_LOGICAL
, gfc_kind_max (op1
, op2
),
506 result
->value
.logical
= op1
->value
.logical
== op2
->value
.logical
;
514 gfc_arith_neqv (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
518 result
= gfc_constant_result (BT_LOGICAL
, gfc_kind_max (op1
, op2
),
520 result
->value
.logical
= op1
->value
.logical
!= op2
->value
.logical
;
527 /* Make sure a constant numeric expression is within the range for
528 its type and kind. Note that there's also a gfc_check_range(),
529 but that one deals with the intrinsic RANGE function. */
532 gfc_range_check (gfc_expr
*e
)
540 rc
= gfc_check_integer_range (e
->value
.integer
, e
->ts
.kind
);
544 rc
= gfc_check_real_range (e
->value
.real
, e
->ts
.kind
);
545 if (rc
== ARITH_UNDERFLOW
)
546 mpfr_set_ui (e
->value
.real
, 0, GFC_RND_MODE
);
547 if (rc
== ARITH_OVERFLOW
)
548 mpfr_set_inf (e
->value
.real
, mpfr_sgn (e
->value
.real
));
550 mpfr_set_nan (e
->value
.real
);
554 rc
= gfc_check_real_range (e
->value
.complex.r
, e
->ts
.kind
);
555 if (rc
== ARITH_UNDERFLOW
)
556 mpfr_set_ui (e
->value
.complex.r
, 0, GFC_RND_MODE
);
557 if (rc
== ARITH_OVERFLOW
)
558 mpfr_set_inf (e
->value
.complex.r
, mpfr_sgn (e
->value
.complex.r
));
560 mpfr_set_nan (e
->value
.complex.r
);
562 rc2
= gfc_check_real_range (e
->value
.complex.i
, e
->ts
.kind
);
563 if (rc
== ARITH_UNDERFLOW
)
564 mpfr_set_ui (e
->value
.complex.i
, 0, GFC_RND_MODE
);
565 if (rc
== ARITH_OVERFLOW
)
566 mpfr_set_inf (e
->value
.complex.i
, mpfr_sgn (e
->value
.complex.i
));
568 mpfr_set_nan (e
->value
.complex.i
);
575 gfc_internal_error ("gfc_range_check(): Bad type");
582 /* Several of the following routines use the same set of statements to
583 check the validity of the result. Encapsulate the checking here. */
586 check_result (arith rc
, gfc_expr
*x
, gfc_expr
*r
, gfc_expr
**rp
)
590 if (val
== ARITH_UNDERFLOW
)
592 if (gfc_option
.warn_underflow
)
593 gfc_warning (gfc_arith_error (val
), &x
->where
);
597 if (val
== ARITH_ASYMMETRIC
)
599 gfc_warning (gfc_arith_error (val
), &x
->where
);
612 /* It may seem silly to have a subroutine that actually computes the
613 unary plus of a constant, but it prevents us from making exceptions
614 in the code elsewhere. Used for unary plus and parenthesized
618 gfc_arith_identity (gfc_expr
*op1
, gfc_expr
**resultp
)
620 *resultp
= gfc_copy_expr (op1
);
626 gfc_arith_uminus (gfc_expr
*op1
, gfc_expr
**resultp
)
631 result
= gfc_constant_result (op1
->ts
.type
, op1
->ts
.kind
, &op1
->where
);
633 switch (op1
->ts
.type
)
636 mpz_neg (result
->value
.integer
, op1
->value
.integer
);
640 mpfr_neg (result
->value
.real
, op1
->value
.real
, GFC_RND_MODE
);
644 mpfr_neg (result
->value
.complex.r
, op1
->value
.complex.r
, GFC_RND_MODE
);
645 mpfr_neg (result
->value
.complex.i
, op1
->value
.complex.i
, GFC_RND_MODE
);
649 gfc_internal_error ("gfc_arith_uminus(): Bad basic type");
652 rc
= gfc_range_check (result
);
654 return check_result (rc
, op1
, result
, resultp
);
659 gfc_arith_plus (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
664 result
= gfc_constant_result (op1
->ts
.type
, op1
->ts
.kind
, &op1
->where
);
666 switch (op1
->ts
.type
)
669 mpz_add (result
->value
.integer
, op1
->value
.integer
, op2
->value
.integer
);
673 mpfr_add (result
->value
.real
, op1
->value
.real
, op2
->value
.real
,
678 mpfr_add (result
->value
.complex.r
, op1
->value
.complex.r
,
679 op2
->value
.complex.r
, GFC_RND_MODE
);
681 mpfr_add (result
->value
.complex.i
, op1
->value
.complex.i
,
682 op2
->value
.complex.i
, GFC_RND_MODE
);
686 gfc_internal_error ("gfc_arith_plus(): Bad basic type");
689 rc
= gfc_range_check (result
);
691 return check_result (rc
, op1
, result
, resultp
);
696 gfc_arith_minus (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
701 result
= gfc_constant_result (op1
->ts
.type
, op1
->ts
.kind
, &op1
->where
);
703 switch (op1
->ts
.type
)
706 mpz_sub (result
->value
.integer
, op1
->value
.integer
, op2
->value
.integer
);
710 mpfr_sub (result
->value
.real
, op1
->value
.real
, op2
->value
.real
,
715 mpfr_sub (result
->value
.complex.r
, op1
->value
.complex.r
,
716 op2
->value
.complex.r
, GFC_RND_MODE
);
718 mpfr_sub (result
->value
.complex.i
, op1
->value
.complex.i
,
719 op2
->value
.complex.i
, GFC_RND_MODE
);
723 gfc_internal_error ("gfc_arith_minus(): Bad basic type");
726 rc
= gfc_range_check (result
);
728 return check_result (rc
, op1
, result
, resultp
);
733 gfc_arith_times (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
739 result
= gfc_constant_result (op1
->ts
.type
, op1
->ts
.kind
, &op1
->where
);
741 switch (op1
->ts
.type
)
744 mpz_mul (result
->value
.integer
, op1
->value
.integer
, op2
->value
.integer
);
748 mpfr_mul (result
->value
.real
, op1
->value
.real
, op2
->value
.real
,
753 gfc_set_model (op1
->value
.complex.r
);
757 mpfr_mul (x
, op1
->value
.complex.r
, op2
->value
.complex.r
, GFC_RND_MODE
);
758 mpfr_mul (y
, op1
->value
.complex.i
, op2
->value
.complex.i
, GFC_RND_MODE
);
759 mpfr_sub (result
->value
.complex.r
, x
, y
, GFC_RND_MODE
);
761 mpfr_mul (x
, op1
->value
.complex.r
, op2
->value
.complex.i
, GFC_RND_MODE
);
762 mpfr_mul (y
, op1
->value
.complex.i
, op2
->value
.complex.r
, GFC_RND_MODE
);
763 mpfr_add (result
->value
.complex.i
, x
, y
, GFC_RND_MODE
);
770 gfc_internal_error ("gfc_arith_times(): Bad basic type");
773 rc
= gfc_range_check (result
);
775 return check_result (rc
, op1
, result
, resultp
);
780 gfc_arith_divide (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
788 result
= gfc_constant_result (op1
->ts
.type
, op1
->ts
.kind
, &op1
->where
);
790 switch (op1
->ts
.type
)
793 if (mpz_sgn (op2
->value
.integer
) == 0)
799 mpz_tdiv_q (result
->value
.integer
, op1
->value
.integer
,
804 if (mpfr_sgn (op2
->value
.real
) == 0 && gfc_option
.flag_range_check
== 1)
810 mpfr_div (result
->value
.real
, op1
->value
.real
, op2
->value
.real
,
815 if (mpfr_sgn (op2
->value
.complex.r
) == 0
816 && mpfr_sgn (op2
->value
.complex.i
) == 0
817 && gfc_option
.flag_range_check
== 1)
823 gfc_set_model (op1
->value
.complex.r
);
828 mpfr_mul (x
, op2
->value
.complex.r
, op2
->value
.complex.r
, GFC_RND_MODE
);
829 mpfr_mul (y
, op2
->value
.complex.i
, op2
->value
.complex.i
, GFC_RND_MODE
);
830 mpfr_add (div
, x
, y
, GFC_RND_MODE
);
832 mpfr_mul (x
, op1
->value
.complex.r
, op2
->value
.complex.r
, GFC_RND_MODE
);
833 mpfr_mul (y
, op1
->value
.complex.i
, op2
->value
.complex.i
, GFC_RND_MODE
);
834 mpfr_add (result
->value
.complex.r
, x
, y
, GFC_RND_MODE
);
835 mpfr_div (result
->value
.complex.r
, result
->value
.complex.r
, div
,
838 mpfr_mul (x
, op1
->value
.complex.i
, op2
->value
.complex.r
, GFC_RND_MODE
);
839 mpfr_mul (y
, op1
->value
.complex.r
, op2
->value
.complex.i
, GFC_RND_MODE
);
840 mpfr_sub (result
->value
.complex.i
, x
, y
, GFC_RND_MODE
);
841 mpfr_div (result
->value
.complex.i
, result
->value
.complex.i
, div
,
850 gfc_internal_error ("gfc_arith_divide(): Bad basic type");
854 rc
= gfc_range_check (result
);
856 return check_result (rc
, op1
, result
, resultp
);
860 /* Compute the reciprocal of a complex number (guaranteed nonzero). */
863 complex_reciprocal (gfc_expr
*op
)
865 mpfr_t mod
, a
, re
, im
;
867 gfc_set_model (op
->value
.complex.r
);
873 mpfr_mul (mod
, op
->value
.complex.r
, op
->value
.complex.r
, GFC_RND_MODE
);
874 mpfr_mul (a
, op
->value
.complex.i
, op
->value
.complex.i
, GFC_RND_MODE
);
875 mpfr_add (mod
, mod
, a
, GFC_RND_MODE
);
877 mpfr_div (re
, op
->value
.complex.r
, mod
, GFC_RND_MODE
);
879 mpfr_neg (im
, op
->value
.complex.i
, GFC_RND_MODE
);
880 mpfr_div (im
, im
, mod
, GFC_RND_MODE
);
882 mpfr_set (op
->value
.complex.r
, re
, GFC_RND_MODE
);
883 mpfr_set (op
->value
.complex.i
, im
, GFC_RND_MODE
);
892 /* Raise a complex number to positive power (power > 0).
893 This function will modify the content of power.
895 Use Binary Method, which is not an optimal but a simple and reasonable
896 arithmetic. See section 4.6.3, "Evaluation of Powers" of Donald E. Knuth,
897 "Seminumerical Algorithms", Vol. 2, "The Art of Computer Programming",
898 3rd Edition, 1998. */
901 complex_pow (gfc_expr
*result
, gfc_expr
*base
, mpz_t power
)
903 mpfr_t x_r
, x_i
, tmp
, re
, im
;
905 gfc_set_model (base
->value
.complex.r
);
913 mpfr_set_ui (result
->value
.complex.r
, 1, GFC_RND_MODE
);
914 mpfr_set_ui (result
->value
.complex.i
, 0, GFC_RND_MODE
);
917 mpfr_set (x_r
, base
->value
.complex.r
, GFC_RND_MODE
);
918 mpfr_set (x_i
, base
->value
.complex.i
, GFC_RND_MODE
);
920 /* Macro for complex multiplication. We have to take care that
921 res_r/res_i and a_r/a_i can (and will) be the same variable. */
922 #define CMULT(res_r,res_i,a_r,a_i,b_r,b_i) \
923 mpfr_mul (re, a_r, b_r, GFC_RND_MODE), \
924 mpfr_mul (tmp, a_i, b_i, GFC_RND_MODE), \
925 mpfr_sub (re, re, tmp, GFC_RND_MODE), \
927 mpfr_mul (im, a_r, b_i, GFC_RND_MODE), \
928 mpfr_mul (tmp, a_i, b_r, GFC_RND_MODE), \
929 mpfr_add (res_i, im, tmp, GFC_RND_MODE), \
930 mpfr_set (res_r, re, GFC_RND_MODE)
932 #define res_r result->value.complex.r
933 #define res_i result->value.complex.i
935 /* for (; power > 0; x *= x) */
936 for (; mpz_cmp_si (power
, 0) > 0; CMULT(x_r
,x_i
,x_r
,x_i
,x_r
,x_i
))
938 /* if (power & 1) res = res * x; */
939 if (mpz_congruent_ui_p (power
, 1, 2))
940 CMULT(res_r
,res_i
,res_r
,res_i
,x_r
,x_i
);
943 mpz_fdiv_q_ui (power
, power
, 2);
958 /* Raise a number to an integer power. */
961 gfc_arith_power (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
967 gcc_assert (op2
->expr_type
== EXPR_CONSTANT
&& op2
->ts
.type
== BT_INTEGER
);
970 result
= gfc_constant_result (op1
->ts
.type
, op1
->ts
.kind
, &op1
->where
);
971 power_sign
= mpz_sgn (op2
->value
.integer
);
975 /* Handle something to the zeroth power. Since we're dealing
976 with integral exponents, there is no ambiguity in the
977 limiting procedure used to determine the value of 0**0. */
978 switch (op1
->ts
.type
)
981 mpz_set_ui (result
->value
.integer
, 1);
985 mpfr_set_ui (result
->value
.real
, 1, GFC_RND_MODE
);
989 mpfr_set_ui (result
->value
.complex.r
, 1, GFC_RND_MODE
);
990 mpfr_set_ui (result
->value
.complex.i
, 0, GFC_RND_MODE
);
994 gfc_internal_error ("gfc_arith_power(): Bad base");
999 switch (op1
->ts
.type
)
1005 /* First, we simplify the cases of op1 == 1, 0 or -1. */
1006 if (mpz_cmp_si (op1
->value
.integer
, 1) == 0)
1009 mpz_set_si (result
->value
.integer
, 1);
1011 else if (mpz_cmp_si (op1
->value
.integer
, 0) == 0)
1013 /* 0**op2 == 0, if op2 > 0
1014 0**op2 overflow, if op2 < 0 ; in that case, we
1015 set the result to 0 and return ARITH_DIV0. */
1016 mpz_set_si (result
->value
.integer
, 0);
1017 if (mpz_cmp_si (op2
->value
.integer
, 0) < 0)
1020 else if (mpz_cmp_si (op1
->value
.integer
, -1) == 0)
1022 /* (-1)**op2 == (-1)**(mod(op2,2)) */
1023 unsigned int odd
= mpz_fdiv_ui (op2
->value
.integer
, 2);
1025 mpz_set_si (result
->value
.integer
, -1);
1027 mpz_set_si (result
->value
.integer
, 1);
1029 /* Then, we take care of op2 < 0. */
1030 else if (mpz_cmp_si (op2
->value
.integer
, 0) < 0)
1032 /* if op2 < 0, op1**op2 == 0 because abs(op1) > 1. */
1033 mpz_set_si (result
->value
.integer
, 0);
1035 else if (gfc_extract_int (op2
, &power
) != NULL
)
1037 /* If op2 doesn't fit in an int, the exponentiation will
1038 overflow, because op2 > 0 and abs(op1) > 1. */
1040 int i
= gfc_validate_kind (BT_INTEGER
, result
->ts
.kind
, false);
1042 if (gfc_option
.flag_range_check
)
1043 rc
= ARITH_OVERFLOW
;
1045 /* Still, we want to give the same value as the processor. */
1047 mpz_add_ui (max
, gfc_integer_kinds
[i
].huge
, 1);
1048 mpz_mul_ui (max
, max
, 2);
1049 mpz_powm (result
->value
.integer
, op1
->value
.integer
,
1050 op2
->value
.integer
, max
);
1054 mpz_pow_ui (result
->value
.integer
, op1
->value
.integer
, power
);
1059 mpfr_pow_z (result
->value
.real
, op1
->value
.real
, op2
->value
.integer
,
1067 /* Compute op1**abs(op2) */
1069 mpz_abs (apower
, op2
->value
.integer
);
1070 complex_pow (result
, op1
, apower
);
1073 /* If (op2 < 0), compute the inverse. */
1075 complex_reciprocal (result
);
1086 rc
= gfc_range_check (result
);
1088 return check_result (rc
, op1
, result
, resultp
);
1092 /* Concatenate two string constants. */
1095 gfc_arith_concat (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
1100 result
= gfc_constant_result (BT_CHARACTER
, gfc_default_character_kind
,
1103 len
= op1
->value
.character
.length
+ op2
->value
.character
.length
;
1105 result
->value
.character
.string
= gfc_get_wide_string (len
+ 1);
1106 result
->value
.character
.length
= len
;
1108 memcpy (result
->value
.character
.string
, op1
->value
.character
.string
,
1109 op1
->value
.character
.length
* sizeof (gfc_char_t
));
1111 memcpy (&result
->value
.character
.string
[op1
->value
.character
.length
],
1112 op2
->value
.character
.string
,
1113 op2
->value
.character
.length
* sizeof (gfc_char_t
));
1115 result
->value
.character
.string
[len
] = '\0';
1122 /* Comparison between real values; returns 0 if (op1 .op. op2) is true.
1123 This function mimics mpr_cmp but takes NaN into account. */
1126 compare_real (gfc_expr
*op1
, gfc_expr
*op2
, gfc_intrinsic_op op
)
1132 rc
= mpfr_equal_p (op1
->value
.real
, op2
->value
.real
) ? 0 : 1;
1135 rc
= mpfr_greater_p (op1
->value
.real
, op2
->value
.real
) ? 1 : -1;
1138 rc
= mpfr_greaterequal_p (op1
->value
.real
, op2
->value
.real
) ? 1 : -1;
1141 rc
= mpfr_less_p (op1
->value
.real
, op2
->value
.real
) ? -1 : 1;
1144 rc
= mpfr_lessequal_p (op1
->value
.real
, op2
->value
.real
) ? -1 : 1;
1147 gfc_internal_error ("compare_real(): Bad operator");
1153 /* Comparison operators. Assumes that the two expression nodes
1154 contain two constants of the same type. The op argument is
1155 needed to handle NaN correctly. */
1158 gfc_compare_expr (gfc_expr
*op1
, gfc_expr
*op2
, gfc_intrinsic_op op
)
1162 switch (op1
->ts
.type
)
1165 rc
= mpz_cmp (op1
->value
.integer
, op2
->value
.integer
);
1169 rc
= compare_real (op1
, op2
, op
);
1173 rc
= gfc_compare_string (op1
, op2
);
1177 rc
= ((!op1
->value
.logical
&& op2
->value
.logical
)
1178 || (op1
->value
.logical
&& !op2
->value
.logical
));
1182 gfc_internal_error ("gfc_compare_expr(): Bad basic type");
1189 /* Compare a pair of complex numbers. Naturally, this is only for
1190 equality and nonequality. */
1193 compare_complex (gfc_expr
*op1
, gfc_expr
*op2
)
1195 return (mpfr_equal_p (op1
->value
.complex.r
, op2
->value
.complex.r
)
1196 && mpfr_equal_p (op1
->value
.complex.i
, op2
->value
.complex.i
));
1200 /* Given two constant strings and the inverse collating sequence, compare the
1201 strings. We return -1 for a < b, 0 for a == b and 1 for a > b.
1202 We use the processor's default collating sequence. */
1205 gfc_compare_string (gfc_expr
*a
, gfc_expr
*b
)
1207 int len
, alen
, blen
, i
;
1210 alen
= a
->value
.character
.length
;
1211 blen
= b
->value
.character
.length
;
1213 len
= MAX(alen
, blen
);
1215 for (i
= 0; i
< len
; i
++)
1217 ac
= ((i
< alen
) ? a
->value
.character
.string
[i
] : ' ');
1218 bc
= ((i
< blen
) ? b
->value
.character
.string
[i
] : ' ');
1226 /* Strings are equal */
1232 gfc_compare_with_Cstring (gfc_expr
*a
, const char *b
, bool case_sensitive
)
1234 int len
, alen
, blen
, i
;
1237 alen
= a
->value
.character
.length
;
1240 len
= MAX(alen
, blen
);
1242 for (i
= 0; i
< len
; i
++)
1244 ac
= ((i
< alen
) ? a
->value
.character
.string
[i
] : ' ');
1245 bc
= ((i
< blen
) ? b
[i
] : ' ');
1247 if (!case_sensitive
)
1259 /* Strings are equal */
1264 /* Specific comparison subroutines. */
1267 gfc_arith_eq (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
1271 result
= gfc_constant_result (BT_LOGICAL
, gfc_default_logical_kind
,
1273 result
->value
.logical
= (op1
->ts
.type
== BT_COMPLEX
)
1274 ? compare_complex (op1
, op2
)
1275 : (gfc_compare_expr (op1
, op2
, INTRINSIC_EQ
) == 0);
1283 gfc_arith_ne (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
1287 result
= gfc_constant_result (BT_LOGICAL
, gfc_default_logical_kind
,
1289 result
->value
.logical
= (op1
->ts
.type
== BT_COMPLEX
)
1290 ? !compare_complex (op1
, op2
)
1291 : (gfc_compare_expr (op1
, op2
, INTRINSIC_EQ
) != 0);
1299 gfc_arith_gt (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
1303 result
= gfc_constant_result (BT_LOGICAL
, gfc_default_logical_kind
,
1305 result
->value
.logical
= (gfc_compare_expr (op1
, op2
, INTRINSIC_GT
) > 0);
1313 gfc_arith_ge (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
1317 result
= gfc_constant_result (BT_LOGICAL
, gfc_default_logical_kind
,
1319 result
->value
.logical
= (gfc_compare_expr (op1
, op2
, INTRINSIC_GE
) >= 0);
1327 gfc_arith_lt (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
1331 result
= gfc_constant_result (BT_LOGICAL
, gfc_default_logical_kind
,
1333 result
->value
.logical
= (gfc_compare_expr (op1
, op2
, INTRINSIC_LT
) < 0);
1341 gfc_arith_le (gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**resultp
)
1345 result
= gfc_constant_result (BT_LOGICAL
, gfc_default_logical_kind
,
1347 result
->value
.logical
= (gfc_compare_expr (op1
, op2
, INTRINSIC_LE
) <= 0);
1355 reduce_unary (arith (*eval
) (gfc_expr
*, gfc_expr
**), gfc_expr
*op
,
1358 gfc_constructor
*c
, *head
;
1362 if (op
->expr_type
== EXPR_CONSTANT
)
1363 return eval (op
, result
);
1366 head
= gfc_copy_constructor (op
->value
.constructor
);
1368 for (c
= head
; c
; c
= c
->next
)
1370 rc
= reduce_unary (eval
, c
->expr
, &r
);
1375 gfc_replace_expr (c
->expr
, r
);
1379 gfc_free_constructor (head
);
1382 r
= gfc_get_expr ();
1383 r
->expr_type
= EXPR_ARRAY
;
1384 r
->value
.constructor
= head
;
1385 r
->shape
= gfc_copy_shape (op
->shape
, op
->rank
);
1387 r
->ts
= head
->expr
->ts
;
1388 r
->where
= op
->where
;
1399 reduce_binary_ac (arith (*eval
) (gfc_expr
*, gfc_expr
*, gfc_expr
**),
1400 gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**result
)
1402 gfc_constructor
*c
, *head
;
1406 head
= gfc_copy_constructor (op1
->value
.constructor
);
1409 for (c
= head
; c
; c
= c
->next
)
1411 if (c
->expr
->expr_type
== EXPR_CONSTANT
)
1412 rc
= eval (c
->expr
, op2
, &r
);
1414 rc
= reduce_binary_ac (eval
, c
->expr
, op2
, &r
);
1419 gfc_replace_expr (c
->expr
, r
);
1423 gfc_free_constructor (head
);
1426 r
= gfc_get_expr ();
1427 r
->expr_type
= EXPR_ARRAY
;
1428 r
->value
.constructor
= head
;
1429 r
->shape
= gfc_copy_shape (op1
->shape
, op1
->rank
);
1431 r
->ts
= head
->expr
->ts
;
1432 r
->where
= op1
->where
;
1433 r
->rank
= op1
->rank
;
1443 reduce_binary_ca (arith (*eval
) (gfc_expr
*, gfc_expr
*, gfc_expr
**),
1444 gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**result
)
1446 gfc_constructor
*c
, *head
;
1450 head
= gfc_copy_constructor (op2
->value
.constructor
);
1453 for (c
= head
; c
; c
= c
->next
)
1455 if (c
->expr
->expr_type
== EXPR_CONSTANT
)
1456 rc
= eval (op1
, c
->expr
, &r
);
1458 rc
= reduce_binary_ca (eval
, op1
, c
->expr
, &r
);
1463 gfc_replace_expr (c
->expr
, r
);
1467 gfc_free_constructor (head
);
1470 r
= gfc_get_expr ();
1471 r
->expr_type
= EXPR_ARRAY
;
1472 r
->value
.constructor
= head
;
1473 r
->shape
= gfc_copy_shape (op2
->shape
, op2
->rank
);
1475 r
->ts
= head
->expr
->ts
;
1476 r
->where
= op2
->where
;
1477 r
->rank
= op2
->rank
;
1486 /* We need a forward declaration of reduce_binary. */
1487 static arith
reduce_binary (arith (*eval
) (gfc_expr
*, gfc_expr
*, gfc_expr
**),
1488 gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**result
);
1492 reduce_binary_aa (arith (*eval
) (gfc_expr
*, gfc_expr
*, gfc_expr
**),
1493 gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**result
)
1495 gfc_constructor
*c
, *d
, *head
;
1499 head
= gfc_copy_constructor (op1
->value
.constructor
);
1502 d
= op2
->value
.constructor
;
1504 if (gfc_check_conformance ("elemental binary operation", op1
, op2
)
1506 rc
= ARITH_INCOMMENSURATE
;
1509 for (c
= head
; c
; c
= c
->next
, d
= d
->next
)
1513 rc
= ARITH_INCOMMENSURATE
;
1517 rc
= reduce_binary (eval
, c
->expr
, d
->expr
, &r
);
1521 gfc_replace_expr (c
->expr
, r
);
1525 rc
= ARITH_INCOMMENSURATE
;
1529 gfc_free_constructor (head
);
1532 r
= gfc_get_expr ();
1533 r
->expr_type
= EXPR_ARRAY
;
1534 r
->value
.constructor
= head
;
1535 r
->shape
= gfc_copy_shape (op1
->shape
, op1
->rank
);
1537 r
->ts
= head
->expr
->ts
;
1538 r
->where
= op1
->where
;
1539 r
->rank
= op1
->rank
;
1549 reduce_binary (arith (*eval
) (gfc_expr
*, gfc_expr
*, gfc_expr
**),
1550 gfc_expr
*op1
, gfc_expr
*op2
, gfc_expr
**result
)
1552 if (op1
->expr_type
== EXPR_CONSTANT
&& op2
->expr_type
== EXPR_CONSTANT
)
1553 return eval (op1
, op2
, result
);
1555 if (op1
->expr_type
== EXPR_CONSTANT
&& op2
->expr_type
== EXPR_ARRAY
)
1556 return reduce_binary_ca (eval
, op1
, op2
, result
);
1558 if (op1
->expr_type
== EXPR_ARRAY
&& op2
->expr_type
== EXPR_CONSTANT
)
1559 return reduce_binary_ac (eval
, op1
, op2
, result
);
1561 return reduce_binary_aa (eval
, op1
, op2
, result
);
1567 arith (*f2
)(gfc_expr
*, gfc_expr
**);
1568 arith (*f3
)(gfc_expr
*, gfc_expr
*, gfc_expr
**);
1572 /* High level arithmetic subroutines. These subroutines go into
1573 eval_intrinsic(), which can do one of several things to its
1574 operands. If the operands are incompatible with the intrinsic
1575 operation, we return a node pointing to the operands and hope that
1576 an operator interface is found during resolution.
1578 If the operands are compatible and are constants, then we try doing
1579 the arithmetic. We also handle the cases where either or both
1580 operands are array constructors. */
1583 eval_intrinsic (gfc_intrinsic_op
operator,
1584 eval_f eval
, gfc_expr
*op1
, gfc_expr
*op2
)
1586 gfc_expr temp
, *result
;
1590 gfc_clear_ts (&temp
.ts
);
1596 if (op1
->ts
.type
!= BT_LOGICAL
)
1599 temp
.ts
.type
= BT_LOGICAL
;
1600 temp
.ts
.kind
= gfc_default_logical_kind
;
1604 /* Logical binary operators */
1607 case INTRINSIC_NEQV
:
1609 if (op1
->ts
.type
!= BT_LOGICAL
|| op2
->ts
.type
!= BT_LOGICAL
)
1612 temp
.ts
.type
= BT_LOGICAL
;
1613 temp
.ts
.kind
= gfc_default_logical_kind
;
1618 case INTRINSIC_UPLUS
:
1619 case INTRINSIC_UMINUS
:
1620 if (!gfc_numeric_ts (&op1
->ts
))
1627 case INTRINSIC_PARENTHESES
:
1632 /* Additional restrictions for ordering relations. */
1634 case INTRINSIC_GE_OS
:
1636 case INTRINSIC_LT_OS
:
1638 case INTRINSIC_LE_OS
:
1640 case INTRINSIC_GT_OS
:
1641 if (op1
->ts
.type
== BT_COMPLEX
|| op2
->ts
.type
== BT_COMPLEX
)
1643 temp
.ts
.type
= BT_LOGICAL
;
1644 temp
.ts
.kind
= gfc_default_logical_kind
;
1650 case INTRINSIC_EQ_OS
:
1652 case INTRINSIC_NE_OS
:
1653 if (op1
->ts
.type
== BT_CHARACTER
&& op2
->ts
.type
== BT_CHARACTER
)
1656 temp
.ts
.type
= BT_LOGICAL
;
1657 temp
.ts
.kind
= gfc_default_logical_kind
;
1662 /* Numeric binary */
1663 case INTRINSIC_PLUS
:
1664 case INTRINSIC_MINUS
:
1665 case INTRINSIC_TIMES
:
1666 case INTRINSIC_DIVIDE
:
1667 case INTRINSIC_POWER
:
1668 if (!gfc_numeric_ts (&op1
->ts
) || !gfc_numeric_ts (&op2
->ts
))
1671 /* Insert any necessary type conversions to make the operands
1674 temp
.expr_type
= EXPR_OP
;
1675 gfc_clear_ts (&temp
.ts
);
1676 temp
.value
.op
.operator = operator;
1678 temp
.value
.op
.op1
= op1
;
1679 temp
.value
.op
.op2
= op2
;
1681 gfc_type_convert_binary (&temp
);
1683 if (operator == INTRINSIC_EQ
|| operator == INTRINSIC_NE
1684 || operator == INTRINSIC_GE
|| operator == INTRINSIC_GT
1685 || operator == INTRINSIC_LE
|| operator == INTRINSIC_LT
1686 || operator == INTRINSIC_EQ_OS
|| operator == INTRINSIC_NE_OS
1687 || operator == INTRINSIC_GE_OS
|| operator == INTRINSIC_GT_OS
1688 || operator == INTRINSIC_LE_OS
|| operator == INTRINSIC_LT_OS
)
1690 temp
.ts
.type
= BT_LOGICAL
;
1691 temp
.ts
.kind
= gfc_default_logical_kind
;
1697 /* Character binary */
1698 case INTRINSIC_CONCAT
:
1699 if (op1
->ts
.type
!= BT_CHARACTER
|| op2
->ts
.type
!= BT_CHARACTER
)
1702 temp
.ts
.type
= BT_CHARACTER
;
1703 temp
.ts
.kind
= gfc_default_character_kind
;
1707 case INTRINSIC_USER
:
1711 gfc_internal_error ("eval_intrinsic(): Bad operator");
1714 /* Try to combine the operators. */
1715 if (operator == INTRINSIC_POWER
&& op2
->ts
.type
!= BT_INTEGER
)
1718 if (op1
->expr_type
!= EXPR_CONSTANT
1719 && (op1
->expr_type
!= EXPR_ARRAY
1720 || !gfc_is_constant_expr (op1
) || !gfc_expanded_ac (op1
)))
1724 && op2
->expr_type
!= EXPR_CONSTANT
1725 && (op2
->expr_type
!= EXPR_ARRAY
1726 || !gfc_is_constant_expr (op2
) || !gfc_expanded_ac (op2
)))
1730 rc
= reduce_unary (eval
.f2
, op1
, &result
);
1732 rc
= reduce_binary (eval
.f3
, op1
, op2
, &result
);
1735 { /* Something went wrong. */
1736 gfc_error (gfc_arith_error (rc
), &op1
->where
);
1740 gfc_free_expr (op1
);
1741 gfc_free_expr (op2
);
1745 /* Create a run-time expression. */
1746 result
= gfc_get_expr ();
1747 result
->ts
= temp
.ts
;
1749 result
->expr_type
= EXPR_OP
;
1750 result
->value
.op
.operator = operator;
1752 result
->value
.op
.op1
= op1
;
1753 result
->value
.op
.op2
= op2
;
1755 result
->where
= op1
->where
;
1761 /* Modify type of expression for zero size array. */
1764 eval_type_intrinsic0 (gfc_intrinsic_op
operator, gfc_expr
*op
)
1767 gfc_internal_error ("eval_type_intrinsic0(): op NULL");
1772 case INTRINSIC_GE_OS
:
1774 case INTRINSIC_LT_OS
:
1776 case INTRINSIC_LE_OS
:
1778 case INTRINSIC_GT_OS
:
1780 case INTRINSIC_EQ_OS
:
1782 case INTRINSIC_NE_OS
:
1783 op
->ts
.type
= BT_LOGICAL
;
1784 op
->ts
.kind
= gfc_default_logical_kind
;
1795 /* Return nonzero if the expression is a zero size array. */
1798 gfc_zero_size_array (gfc_expr
*e
)
1800 if (e
->expr_type
!= EXPR_ARRAY
)
1803 return e
->value
.constructor
== NULL
;
1807 /* Reduce a binary expression where at least one of the operands
1808 involves a zero-length array. Returns NULL if neither of the
1809 operands is a zero-length array. */
1812 reduce_binary0 (gfc_expr
*op1
, gfc_expr
*op2
)
1814 if (gfc_zero_size_array (op1
))
1816 gfc_free_expr (op2
);
1820 if (gfc_zero_size_array (op2
))
1822 gfc_free_expr (op1
);
1831 eval_intrinsic_f2 (gfc_intrinsic_op
operator,
1832 arith (*eval
) (gfc_expr
*, gfc_expr
**),
1833 gfc_expr
*op1
, gfc_expr
*op2
)
1840 if (gfc_zero_size_array (op1
))
1841 return eval_type_intrinsic0 (operator, op1
);
1845 result
= reduce_binary0 (op1
, op2
);
1847 return eval_type_intrinsic0 (operator, result
);
1851 return eval_intrinsic (operator, f
, op1
, op2
);
1856 eval_intrinsic_f3 (gfc_intrinsic_op
operator,
1857 arith (*eval
) (gfc_expr
*, gfc_expr
*, gfc_expr
**),
1858 gfc_expr
*op1
, gfc_expr
*op2
)
1863 result
= reduce_binary0 (op1
, op2
);
1865 return eval_type_intrinsic0(operator, result
);
1868 return eval_intrinsic (operator, f
, op1
, op2
);
1873 gfc_parentheses (gfc_expr
*op
)
1875 if (gfc_is_constant_expr (op
))
1878 return eval_intrinsic_f2 (INTRINSIC_PARENTHESES
, gfc_arith_identity
,
1883 gfc_uplus (gfc_expr
*op
)
1885 return eval_intrinsic_f2 (INTRINSIC_UPLUS
, gfc_arith_identity
, op
, NULL
);
1890 gfc_uminus (gfc_expr
*op
)
1892 return eval_intrinsic_f2 (INTRINSIC_UMINUS
, gfc_arith_uminus
, op
, NULL
);
1897 gfc_add (gfc_expr
*op1
, gfc_expr
*op2
)
1899 return eval_intrinsic_f3 (INTRINSIC_PLUS
, gfc_arith_plus
, op1
, op2
);
1904 gfc_subtract (gfc_expr
*op1
, gfc_expr
*op2
)
1906 return eval_intrinsic_f3 (INTRINSIC_MINUS
, gfc_arith_minus
, op1
, op2
);
1911 gfc_multiply (gfc_expr
*op1
, gfc_expr
*op2
)
1913 return eval_intrinsic_f3 (INTRINSIC_TIMES
, gfc_arith_times
, op1
, op2
);
1918 gfc_divide (gfc_expr
*op1
, gfc_expr
*op2
)
1920 return eval_intrinsic_f3 (INTRINSIC_DIVIDE
, gfc_arith_divide
, op1
, op2
);
1925 gfc_power (gfc_expr
*op1
, gfc_expr
*op2
)
1927 return eval_intrinsic_f3 (INTRINSIC_POWER
, gfc_arith_power
, op1
, op2
);
1932 gfc_concat (gfc_expr
*op1
, gfc_expr
*op2
)
1934 return eval_intrinsic_f3 (INTRINSIC_CONCAT
, gfc_arith_concat
, op1
, op2
);
1939 gfc_and (gfc_expr
*op1
, gfc_expr
*op2
)
1941 return eval_intrinsic_f3 (INTRINSIC_AND
, gfc_arith_and
, op1
, op2
);
1946 gfc_or (gfc_expr
*op1
, gfc_expr
*op2
)
1948 return eval_intrinsic_f3 (INTRINSIC_OR
, gfc_arith_or
, op1
, op2
);
1953 gfc_not (gfc_expr
*op1
)
1955 return eval_intrinsic_f2 (INTRINSIC_NOT
, gfc_arith_not
, op1
, NULL
);
1960 gfc_eqv (gfc_expr
*op1
, gfc_expr
*op2
)
1962 return eval_intrinsic_f3 (INTRINSIC_EQV
, gfc_arith_eqv
, op1
, op2
);
1967 gfc_neqv (gfc_expr
*op1
, gfc_expr
*op2
)
1969 return eval_intrinsic_f3 (INTRINSIC_NEQV
, gfc_arith_neqv
, op1
, op2
);
1974 gfc_eq (gfc_expr
*op1
, gfc_expr
*op2
, gfc_intrinsic_op op
)
1976 return eval_intrinsic_f3 (op
, gfc_arith_eq
, op1
, op2
);
1981 gfc_ne (gfc_expr
*op1
, gfc_expr
*op2
, gfc_intrinsic_op op
)
1983 return eval_intrinsic_f3 (op
, gfc_arith_ne
, op1
, op2
);
1988 gfc_gt (gfc_expr
*op1
, gfc_expr
*op2
, gfc_intrinsic_op op
)
1990 return eval_intrinsic_f3 (op
, gfc_arith_gt
, op1
, op2
);
1995 gfc_ge (gfc_expr
*op1
, gfc_expr
*op2
, gfc_intrinsic_op op
)
1997 return eval_intrinsic_f3 (op
, gfc_arith_ge
, op1
, op2
);
2002 gfc_lt (gfc_expr
*op1
, gfc_expr
*op2
, gfc_intrinsic_op op
)
2004 return eval_intrinsic_f3 (op
, gfc_arith_lt
, op1
, op2
);
2009 gfc_le (gfc_expr
*op1
, gfc_expr
*op2
, gfc_intrinsic_op op
)
2011 return eval_intrinsic_f3 (op
, gfc_arith_le
, op1
, op2
);
2015 /* Convert an integer string to an expression node. */
2018 gfc_convert_integer (const char *buffer
, int kind
, int radix
, locus
*where
)
2023 e
= gfc_constant_result (BT_INTEGER
, kind
, where
);
2024 /* A leading plus is allowed, but not by mpz_set_str. */
2025 if (buffer
[0] == '+')
2029 mpz_set_str (e
->value
.integer
, t
, radix
);
2035 /* Convert a real string to an expression node. */
2038 gfc_convert_real (const char *buffer
, int kind
, locus
*where
)
2042 e
= gfc_constant_result (BT_REAL
, kind
, where
);
2043 mpfr_set_str (e
->value
.real
, buffer
, 10, GFC_RND_MODE
);
2049 /* Convert a pair of real, constant expression nodes to a single
2050 complex expression node. */
2053 gfc_convert_complex (gfc_expr
*real
, gfc_expr
*imag
, int kind
)
2057 e
= gfc_constant_result (BT_COMPLEX
, kind
, &real
->where
);
2058 mpfr_set (e
->value
.complex.r
, real
->value
.real
, GFC_RND_MODE
);
2059 mpfr_set (e
->value
.complex.i
, imag
->value
.real
, GFC_RND_MODE
);
2065 /******* Simplification of intrinsic functions with constant arguments *****/
2068 /* Deal with an arithmetic error. */
2071 arith_error (arith rc
, gfc_typespec
*from
, gfc_typespec
*to
, locus
*where
)
2076 gfc_error ("Arithmetic OK converting %s to %s at %L",
2077 gfc_typename (from
), gfc_typename (to
), where
);
2079 case ARITH_OVERFLOW
:
2080 gfc_error ("Arithmetic overflow converting %s to %s at %L. This check "
2081 "can be disabled with the option -fno-range-check",
2082 gfc_typename (from
), gfc_typename (to
), where
);
2084 case ARITH_UNDERFLOW
:
2085 gfc_error ("Arithmetic underflow converting %s to %s at %L",
2086 gfc_typename (from
), gfc_typename (to
), where
);
2089 gfc_error ("Arithmetic NaN converting %s to %s at %L",
2090 gfc_typename (from
), gfc_typename (to
), where
);
2093 gfc_error ("Division by zero converting %s to %s at %L",
2094 gfc_typename (from
), gfc_typename (to
), where
);
2096 case ARITH_INCOMMENSURATE
:
2097 gfc_error ("Array operands are incommensurate converting %s to %s at %L",
2098 gfc_typename (from
), gfc_typename (to
), where
);
2100 case ARITH_ASYMMETRIC
:
2101 gfc_error ("Integer outside symmetric range implied by Standard Fortran"
2102 " converting %s to %s at %L",
2103 gfc_typename (from
), gfc_typename (to
), where
);
2106 gfc_internal_error ("gfc_arith_error(): Bad error code");
2109 /* TODO: Do something about the error, ie, throw exception, return
2114 /* Convert integers to integers. */
2117 gfc_int2int (gfc_expr
*src
, int kind
)
2122 result
= gfc_constant_result (BT_INTEGER
, kind
, &src
->where
);
2124 mpz_set (result
->value
.integer
, src
->value
.integer
);
2126 if ((rc
= gfc_check_integer_range (result
->value
.integer
, kind
)) != ARITH_OK
)
2128 if (rc
== ARITH_ASYMMETRIC
)
2130 gfc_warning (gfc_arith_error (rc
), &src
->where
);
2134 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2135 gfc_free_expr (result
);
2144 /* Convert integers to reals. */
2147 gfc_int2real (gfc_expr
*src
, int kind
)
2152 result
= gfc_constant_result (BT_REAL
, kind
, &src
->where
);
2154 mpfr_set_z (result
->value
.real
, src
->value
.integer
, GFC_RND_MODE
);
2156 if ((rc
= gfc_check_real_range (result
->value
.real
, kind
)) != ARITH_OK
)
2158 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2159 gfc_free_expr (result
);
2167 /* Convert default integer to default complex. */
2170 gfc_int2complex (gfc_expr
*src
, int kind
)
2175 result
= gfc_constant_result (BT_COMPLEX
, kind
, &src
->where
);
2177 mpfr_set_z (result
->value
.complex.r
, src
->value
.integer
, GFC_RND_MODE
);
2178 mpfr_set_ui (result
->value
.complex.i
, 0, GFC_RND_MODE
);
2180 if ((rc
= gfc_check_real_range (result
->value
.complex.r
, kind
)) != ARITH_OK
)
2182 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2183 gfc_free_expr (result
);
2191 /* Convert default real to default integer. */
2194 gfc_real2int (gfc_expr
*src
, int kind
)
2199 result
= gfc_constant_result (BT_INTEGER
, kind
, &src
->where
);
2201 gfc_mpfr_to_mpz (result
->value
.integer
, src
->value
.real
);
2203 if ((rc
= gfc_check_integer_range (result
->value
.integer
, kind
)) != ARITH_OK
)
2205 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2206 gfc_free_expr (result
);
2214 /* Convert real to real. */
2217 gfc_real2real (gfc_expr
*src
, int kind
)
2222 result
= gfc_constant_result (BT_REAL
, kind
, &src
->where
);
2224 mpfr_set (result
->value
.real
, src
->value
.real
, GFC_RND_MODE
);
2226 rc
= gfc_check_real_range (result
->value
.real
, kind
);
2228 if (rc
== ARITH_UNDERFLOW
)
2230 if (gfc_option
.warn_underflow
)
2231 gfc_warning (gfc_arith_error (rc
), &src
->where
);
2232 mpfr_set_ui (result
->value
.real
, 0, GFC_RND_MODE
);
2234 else if (rc
!= ARITH_OK
)
2236 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2237 gfc_free_expr (result
);
2245 /* Convert real to complex. */
2248 gfc_real2complex (gfc_expr
*src
, int kind
)
2253 result
= gfc_constant_result (BT_COMPLEX
, kind
, &src
->where
);
2255 mpfr_set (result
->value
.complex.r
, src
->value
.real
, GFC_RND_MODE
);
2256 mpfr_set_ui (result
->value
.complex.i
, 0, GFC_RND_MODE
);
2258 rc
= gfc_check_real_range (result
->value
.complex.r
, kind
);
2260 if (rc
== ARITH_UNDERFLOW
)
2262 if (gfc_option
.warn_underflow
)
2263 gfc_warning (gfc_arith_error (rc
), &src
->where
);
2264 mpfr_set_ui (result
->value
.complex.r
, 0, GFC_RND_MODE
);
2266 else if (rc
!= ARITH_OK
)
2268 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2269 gfc_free_expr (result
);
2277 /* Convert complex to integer. */
2280 gfc_complex2int (gfc_expr
*src
, int kind
)
2285 result
= gfc_constant_result (BT_INTEGER
, kind
, &src
->where
);
2287 gfc_mpfr_to_mpz (result
->value
.integer
, src
->value
.complex.r
);
2289 if ((rc
= gfc_check_integer_range (result
->value
.integer
, kind
)) != ARITH_OK
)
2291 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2292 gfc_free_expr (result
);
2300 /* Convert complex to real. */
2303 gfc_complex2real (gfc_expr
*src
, int kind
)
2308 result
= gfc_constant_result (BT_REAL
, kind
, &src
->where
);
2310 mpfr_set (result
->value
.real
, src
->value
.complex.r
, GFC_RND_MODE
);
2312 rc
= gfc_check_real_range (result
->value
.real
, kind
);
2314 if (rc
== ARITH_UNDERFLOW
)
2316 if (gfc_option
.warn_underflow
)
2317 gfc_warning (gfc_arith_error (rc
), &src
->where
);
2318 mpfr_set_ui (result
->value
.real
, 0, GFC_RND_MODE
);
2322 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2323 gfc_free_expr (result
);
2331 /* Convert complex to complex. */
2334 gfc_complex2complex (gfc_expr
*src
, int kind
)
2339 result
= gfc_constant_result (BT_COMPLEX
, kind
, &src
->where
);
2341 mpfr_set (result
->value
.complex.r
, src
->value
.complex.r
, GFC_RND_MODE
);
2342 mpfr_set (result
->value
.complex.i
, src
->value
.complex.i
, GFC_RND_MODE
);
2344 rc
= gfc_check_real_range (result
->value
.complex.r
, kind
);
2346 if (rc
== ARITH_UNDERFLOW
)
2348 if (gfc_option
.warn_underflow
)
2349 gfc_warning (gfc_arith_error (rc
), &src
->where
);
2350 mpfr_set_ui (result
->value
.complex.r
, 0, GFC_RND_MODE
);
2352 else if (rc
!= ARITH_OK
)
2354 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2355 gfc_free_expr (result
);
2359 rc
= gfc_check_real_range (result
->value
.complex.i
, kind
);
2361 if (rc
== ARITH_UNDERFLOW
)
2363 if (gfc_option
.warn_underflow
)
2364 gfc_warning (gfc_arith_error (rc
), &src
->where
);
2365 mpfr_set_ui (result
->value
.complex.i
, 0, GFC_RND_MODE
);
2367 else if (rc
!= ARITH_OK
)
2369 arith_error (rc
, &src
->ts
, &result
->ts
, &src
->where
);
2370 gfc_free_expr (result
);
2378 /* Logical kind conversion. */
2381 gfc_log2log (gfc_expr
*src
, int kind
)
2385 result
= gfc_constant_result (BT_LOGICAL
, kind
, &src
->where
);
2386 result
->value
.logical
= src
->value
.logical
;
2392 /* Convert logical to integer. */
2395 gfc_log2int (gfc_expr
*src
, int kind
)
2399 result
= gfc_constant_result (BT_INTEGER
, kind
, &src
->where
);
2400 mpz_set_si (result
->value
.integer
, src
->value
.logical
);
2406 /* Convert integer to logical. */
2409 gfc_int2log (gfc_expr
*src
, int kind
)
2413 result
= gfc_constant_result (BT_LOGICAL
, kind
, &src
->where
);
2414 result
->value
.logical
= (mpz_cmp_si (src
->value
.integer
, 0) != 0);
2420 /* Helper function to set the representation in a Hollerith conversion.
2421 This assumes that the ts.type and ts.kind of the result have already
2425 hollerith2representation (gfc_expr
*result
, gfc_expr
*src
)
2427 int src_len
, result_len
;
2429 src_len
= src
->representation
.length
;
2430 result_len
= gfc_target_expr_size (result
);
2432 if (src_len
> result_len
)
2434 gfc_warning ("The Hollerith constant at %L is too long to convert to %s",
2435 &src
->where
, gfc_typename(&result
->ts
));
2438 result
->representation
.string
= gfc_getmem (result_len
+ 1);
2439 memcpy (result
->representation
.string
, src
->representation
.string
,
2440 MIN (result_len
, src_len
));
2442 if (src_len
< result_len
)
2443 memset (&result
->representation
.string
[src_len
], ' ', result_len
- src_len
);
2445 result
->representation
.string
[result_len
] = '\0'; /* For debugger */
2446 result
->representation
.length
= result_len
;
2450 /* Convert Hollerith to integer. The constant will be padded or truncated. */
2453 gfc_hollerith2int (gfc_expr
*src
, int kind
)
2457 result
= gfc_get_expr ();
2458 result
->expr_type
= EXPR_CONSTANT
;
2459 result
->ts
.type
= BT_INTEGER
;
2460 result
->ts
.kind
= kind
;
2461 result
->where
= src
->where
;
2463 hollerith2representation (result
, src
);
2464 gfc_interpret_integer (kind
, (unsigned char *) result
->representation
.string
,
2465 result
->representation
.length
, result
->value
.integer
);
2471 /* Convert Hollerith to real. The constant will be padded or truncated. */
2474 gfc_hollerith2real (gfc_expr
*src
, int kind
)
2479 len
= src
->value
.character
.length
;
2481 result
= gfc_get_expr ();
2482 result
->expr_type
= EXPR_CONSTANT
;
2483 result
->ts
.type
= BT_REAL
;
2484 result
->ts
.kind
= kind
;
2485 result
->where
= src
->where
;
2487 hollerith2representation (result
, src
);
2488 gfc_interpret_float (kind
, (unsigned char *) result
->representation
.string
,
2489 result
->representation
.length
, result
->value
.real
);
2495 /* Convert Hollerith to complex. The constant will be padded or truncated. */
2498 gfc_hollerith2complex (gfc_expr
*src
, int kind
)
2503 len
= src
->value
.character
.length
;
2505 result
= gfc_get_expr ();
2506 result
->expr_type
= EXPR_CONSTANT
;
2507 result
->ts
.type
= BT_COMPLEX
;
2508 result
->ts
.kind
= kind
;
2509 result
->where
= src
->where
;
2511 hollerith2representation (result
, src
);
2512 gfc_interpret_complex (kind
, (unsigned char *) result
->representation
.string
,
2513 result
->representation
.length
, result
->value
.complex.r
,
2514 result
->value
.complex.i
);
2520 /* Convert Hollerith to character. */
2523 gfc_hollerith2character (gfc_expr
*src
, int kind
)
2527 result
= gfc_copy_expr (src
);
2528 result
->ts
.type
= BT_CHARACTER
;
2529 result
->ts
.kind
= kind
;
2531 result
->value
.character
.length
= result
->representation
.length
;
2532 result
->value
.character
.string
2533 = gfc_char_to_widechar (result
->representation
.string
);
2539 /* Convert Hollerith to logical. The constant will be padded or truncated. */
2542 gfc_hollerith2logical (gfc_expr
*src
, int kind
)
2547 len
= src
->value
.character
.length
;
2549 result
= gfc_get_expr ();
2550 result
->expr_type
= EXPR_CONSTANT
;
2551 result
->ts
.type
= BT_LOGICAL
;
2552 result
->ts
.kind
= kind
;
2553 result
->where
= src
->where
;
2555 hollerith2representation (result
, src
);
2556 gfc_interpret_logical (kind
, (unsigned char *) result
->representation
.string
,
2557 result
->representation
.length
, &result
->value
.logical
);
2563 /* Returns an initializer whose value is one higher than the value of the
2564 LAST_INITIALIZER argument. If the argument is NULL, the
2565 initializers value will be set to zero. The initializer's kind
2566 will be set to gfc_c_int_kind.
2568 If -fshort-enums is given, the appropriate kind will be selected
2569 later after all enumerators have been parsed. A warning is issued
2570 here if an initializer exceeds gfc_c_int_kind. */
2573 gfc_enum_initializer (gfc_expr
*last_initializer
, locus where
)
2577 result
= gfc_get_expr ();
2578 result
->expr_type
= EXPR_CONSTANT
;
2579 result
->ts
.type
= BT_INTEGER
;
2580 result
->ts
.kind
= gfc_c_int_kind
;
2581 result
->where
= where
;
2583 mpz_init (result
->value
.integer
);
2585 if (last_initializer
!= NULL
)
2587 mpz_add_ui (result
->value
.integer
, last_initializer
->value
.integer
, 1);
2588 result
->where
= last_initializer
->where
;
2590 if (gfc_check_integer_range (result
->value
.integer
,
2591 gfc_c_int_kind
) != ARITH_OK
)
2593 gfc_error ("Enumerator exceeds the C integer type at %C");
2599 /* Control comes here, if it's the very first enumerator and no
2600 initializer has been given. It will be initialized to zero. */
2601 mpz_set_si (result
->value
.integer
, 0);