This is the mail archive of the gcc-patches@gcc.gnu.org mailing list for the GCC project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

[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
  

Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]