This is the mail archive of the
fortran@gcc.gnu.org
mailing list for the GNU Fortran project.
[gfortran] Patch for 14921, 14540
- From: Feng Wang <wf_cs at yahoo dot com>
- To: fortran <fortran at gcc dot gnu dot org>, patch <gcc-patches at gcc dot gnu dot org>
- Date: Tue, 13 Apr 2004 21:15:37 +0800 (CST)
- Subject: [gfortran] Patch for 14921, 14540
Hi, all
This patch fix intrinsic function atan2 and clog. The original implementation
arithmetic was error. No additional error on testsuits.
Feng Wang
fortran/ChangeLog entry:
2004-04-13 Feng Wang <fengwang@nudt.edu.cn>
PR fortran/14921
PR fortran/14540
* arith.c (arctangent2): Add arctangent2 implementation.
* arith.h (arctangent2): Add function prototype.
* simplify.c (gfc_simplify_atan2): Use it.
(gfc_simplify_log): Use it.
testsuits/ChangeLog.tree-ssa entry:
2004-04-13 Feng Wang <fengwang@nudt.edu.cn>
PR fortran/14921
PR fortran/14540
* gfortran.fortran-torture/execute/math.f90: Add atan2 and clog
simplify test.
_________________________________________________________
Do You Yahoo!?
惠普TT游戏剧,玩游戏,中大奖!
http://cn.rd.yahoo.com/mail_cn/tag/*http://hp.allyes.com/laserjet/gamestory/index.html?jumpid=ex_hphqapcn_MongooseLJ1010/201073CN407016/Yahoo
diff -c3p fortran/arith.c gcc/gcc/fortran/arith.c
*** fortran/arith.c Mon Apr 12 10:22:45 2004
--- gcc/gcc/fortran/arith.c Tue Apr 13 11:01:38 2004
*************** cosine (mpf_t * arg, mpf_t * result)
*** 390,396 ****
}
-
/* Calculate atan(arg).
Similar to sine but requires special handling for x near 1. */
--- 390,395 ----
*************** arctangent (mpf_t * arg, mpf_t * result)
*** 515,520 ****
--- 514,561 ----
mpf_clear (x);
}
+ /* Calculate atan2 (y, x)
+
+ atan2(y, x) = atan(y/x) if x > 0,
+ sign(y)*(pi - atan(|y/x|)) if x < 0,
+ 0 if x = 0 && y == 0,
+ sign(y)*pi/2 if x = 0 && y != 0.
+ */
+
+ void
+ arctangent2 (mpf_t * y, mpf_t * x, mpf_t * result)
+ {
+ mpf_t t;
+
+ mpf_init (t);
+
+ switch (mpf_sgn (*x))
+ {
+ case 1:
+ mpf_div (t, *y, *x);
+ arctangent (&t, result);
+ break;
+ case -1:
+ mpf_div (t, *y, *x);
+ mpf_abs (t, t);
+ arctangent (&t, &t);
+ mpf_sub (*result, pi, t);
+ if (mpf_sgn (*y) == -1)
+ mpf_neg (*result, *result);
+ break;
+ case 0:
+ if (mpf_sgn (*y) == 0)
+ mpf_set_ui (*result, 0);
+ else
+ {
+ mpf_set (*result, half_pi);
+ if (mpf_sgn (*y) == -1)
+ mpf_neg (*result, *result);
+ }
+ break;
+ }
+ mpf_clear (t);
+ }
/* Calculate cosh(arg). */
diff -c3p fortran/arith.h gcc/gcc/fortran/arith.h
*** fortran/arith.h Mon Apr 12 10:22:45 2004
--- gcc/gcc/fortran/arith.h Tue Apr 13 11:01:38 2004
*************** void exponential (mpf_t *, mpf_t *);
*** 35,40 ****
--- 35,41 ----
void sine (mpf_t *, mpf_t *);
void cosine (mpf_t *, mpf_t *);
void arctangent (mpf_t *, mpf_t *);
+ void arctangent2 (mpf_t *, mpf_t *, mpf_t *);
void hypercos (mpf_t *, mpf_t *);
void hypersine (mpf_t *, mpf_t *);
diff -c3p fortran/simplify.c gcc/gcc/fortran/simplify.c
*** fortran/simplify.c Mon Apr 12 10:22:45 2004
--- gcc/gcc/fortran/simplify.c Tue Apr 13 11:01:38 2004
*************** gfc_expr *
*** 562,622 ****
gfc_simplify_atan2 (gfc_expr * y, gfc_expr * x)
{
gfc_expr *result;
- mpf_t term;
if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
return NULL;
result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
- mpf_init (term);
! if (mpf_cmp_ui (y->value.real, 0) == 0)
{
! if (mpf_cmp_ui (x->value.real, 0) == 0)
! {
! mpf_clear (term);
! gfc_error
! ("If first argument of ATAN2 %L is zero, the second argument "
! "must not be zero", &x->where);
! gfc_free_expr (result);
! return &gfc_bad_expr;
! }
! else if (mpf_cmp_si (x->value.real, 0) < 0)
! {
! mpf_set (result->value.real, pi);
! mpf_clear (term);
! return result;
! }
! else if (mpf_cmp_si (x->value.real, -1) == 0)
! {
! mpf_set_ui (result->value.real, 0);
! mpf_clear (term);
! return range_check (result, "ATAN2");
! }
! }
!
! if (mpf_cmp_ui (x->value.real, 0) == 0)
! {
! if (mpf_cmp_si (y->value.real, 0) < 0)
! {
! mpf_neg (term, half_pi);
! mpf_set (result->value.real, term);
! mpf_clear (term);
! return range_check (result, "ATAN2");
! }
! else if (mpf_cmp_si (y->value.real, 0) > 0)
! {
! mpf_set (result->value.real, half_pi);
! mpf_clear (term);
! return range_check (result, "ATAN2");
! }
! }
!
! mpf_div (term, y->value.real, x->value.real);
! arctangent (&term, &result->value.real);
! mpf_clear (term);
return range_check (result, "ATAN2");
--- 562,584 ----
gfc_simplify_atan2 (gfc_expr * y, gfc_expr * x)
{
gfc_expr *result;
if (x->expr_type != EXPR_CONSTANT || y->expr_type != EXPR_CONSTANT)
return NULL;
result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
! if (mpf_sgn (y->value.real) == 0 && mpf_sgn (x->value.real) == 0)
{
! gfc_error
! ("If first argument of ATAN2 %L is zero, the second argument "
! "must not be zero", &x->where);
! gfc_free_expr (result);
! return &gfc_bad_expr;
! }
! arctangent2 (&y->value.real, &x->value.real, &result->value.real);
return range_check (result, "ATAN2");
*************** gfc_simplify_log (gfc_expr * x)
*** 2073,2079 ****
mpf_init (xi);
mpf_div (xr, x->value.complex.i, x->value.complex.r);
! arctangent (&xr, &result->value.complex.i);
mpf_mul (xr, x->value.complex.r, x->value.complex.r);
mpf_mul (xi, x->value.complex.i, x->value.complex.i);
--- 2035,2042 ----
mpf_init (xi);
mpf_div (xr, x->value.complex.i, x->value.complex.r);
! arctangent2 (&x->value.complex.i, &x->value.complex.r,
! &result->value.complex.i);
mpf_mul (xr, x->value.complex.r, x->value.complex.r);
mpf_mul (xi, x->value.complex.i, x->value.complex.i);
*** /gcc/gcc/testsuite/gfortran.fortran-torture/execute/math.f90 Sat Jul 26 08:27:49 2003
--- math.f90 Tue Apr 13 11:42:22 2004
*************** subroutine dotest (n, val4, val8, known)
*** 9,19 ****
if (abs (real (val8, kind=4) - known) .gt. 0.001) call abort
end subroutine
program testmath
implicit none
real(kind=4) r, two4, half4
real(kind=8) q, two8, half8
! external dotest
two4 = 2.0
two8 = 2.0_8
--- 9,30 ----
if (abs (real (val8, kind=4) - known) .gt. 0.001) call abort
end subroutine
+ subroutine dotestc (n, val4, val8, known)
+ implicit none
+ complex(kind=4) val4, known
+ complex(kind=8) val8
+ integer n
+ if (abs (val4 - known) .gt. 0.001) call abort
+ if (abs (cmplx (val8, kind=4) - known) .gt. 0.001) call abort
+ end subroutine
+
program testmath
implicit none
real(kind=4) r, two4, half4
real(kind=8) q, two8, half8
! complex(kind=4) cr
! complex(kind=8) cq
! external dotest, dotest2
two4 = 2.0
two8 = 2.0_8
*************** program testmath
*** 61,65 ****
--- 72,100 ----
r = sqrt (two4)
q = sqrt (two8)
call dotest (14, r, q, 1.4142)
+
+ r = atan2 (0.0, 1.0)
+ q = atan2 (0.0_8, 1.0_8)
+ call dotest (15, r, q, 0.0)
+ r = atan2 (-1.0, 1.0)
+ q = atan2 (-1.0_8, 1.0_8)
+ call dotest (16, r, q, -0.7854)
+ r = atan2 (0.0, -1.0)
+ q = atan2 (0.0_8, -1.0_8)
+ call dotest (17, r, q, 3.1416)
+ r = atan2 (-1.0, -1.0)
+ q = atan2 (-1.0_8, -1.0_8)
+ call dotest (18, r, q, -2.3562)
+ r = atan2 (1.0, 0.0)
+ q = atan2 (1.0_8, 0.0_8)
+ call dotest (19, r, q, 1.5708)
+ r = atan2 (-1.0, 0.0)
+ q = atan2 (-1.0_8, 0.0_8)
+ call dotest (20, r, q, -1.5708)
+
+ cr = log ((-1.0, -1.0))
+ cq = log ((-1.0_8, -1.0_8))
+ call dotestc (21, cr, cq, (0.3466, -2.3562))
+
end program