Re: [Rd] branch cuts in atan(), asin(), acos() (PR#7871)

2005-05-17 Thread r . hankin

On May 17, 2005, at 07:01 am, Prof Brian Ripley wrote:

> On Mon, 16 May 2005 [EMAIL PROTECTED] wrote:
>
>> The complex versions of atain(), asin(), acos() are not documented
>> except in the source
>> code.
>
> Thank you, but your patch does not address that so they would become 
> undocumented (they no longer follow the comment in the source code, 
> which BTW is regarded as documentation).
>
> Could you please supply documentation to match this change, including 
> the exact definition you have implemented.
>
> Also: have you looked at atan2?
>


OK, this follows PR#7871, discussing branch cuts of inverse trig 
functions.
In this email, I provide a context diff between complex.c (R-2.1.0) and 
a version that
moves the branch cuts to the standard places.  The changes are 
documented in the code.
I also give updated Trig.Rd and Hyperbolic.Rd files which document the 
behaviour
of the new functions.   R-2.1.0 passes "make check"  with the new files 
(after editing
Rdutils.Rd, which caused a problem).

I have not attempted to alter z_atan2(): there appear to be unrelated 
issues here.
For example, atan2(0,0) gives 0 as documented, but atan2(0i,0i) gives 
NA.
Also, atan2(0,1i) returns an error (I would expect a zero here too).

I feel it is better to consider z_atan2() as a separate issue.  I will 
file a separate bugreport
for this.



diff -c  


*** /Users/rksh/downloads/R-2.1.0/src/main/complex.cMon Apr 18 
22:30:00 2005
--- complex.c   Tue May 17 13:35:27 2005
***
*** 355,360 
--- 355,365 

/* Complex Arcsin and Arccos Functions */
/* Equation (4.4.37) Abramowitz and Stegun */
+   /* with additional terms to force the branch
+   /* to agree with figure 4.4, p79.  Continuity */
+   /* on the branch cuts (real axis; y==0, |x| > 1) is */
+   /* standard: z_asin() is continuous from below if x >= 1 */
+   /* and continuous from above if x <= -1. */

   static void z_asin(Rcomplex *r, Rcomplex *z)
   {
***
*** 367,372 
--- 372,382 
   bet = t1 - t2;
   r->r = asin(bet);
   r->i = log(alpha + sqrt(alpha*alpha - 1));
+
+ if(y<0 || (y==0 && x > 1)){
+   r->i *= -1;
+   }
+
   }

   static void z_acos(Rcomplex *r, Rcomplex *z)
***
*** 379,384 
--- 389,399 

/* Complex Arctangent Function */
/* Equation (4.4.39) Abramowitz and Stegun */
+   /* with additional terms to force the branch cuts */
+   /* to agree with figure 4.4, p79.  Continuity */
+   /* on the branch cuts (pure imaginary axis; x==0, |y|>1) */
+   /* is standard: z_asin() is continuous from the right */
+   /*  if y >= 1, and continuous from the left if y <= -1. */

   static void z_atan(Rcomplex *r, Rcomplex *z)
   {
***
*** 388,393 
--- 403,416 
   r->r = 0.5 * atan(2 * x / ( 1 - x * x - y * y));
   r->i = 0.25 * log((x * x + (y + 1) * (y + 1)) /
  (x * x + (y - 1) * (y - 1)));
+
+ if(x*x + y*y > 1){
+   r->r += M_PI_2;
+   if(x < 0 || (x==0 && y<0)){
+   r->r -= M_PI;
+   }
+ }
+
   }

   static void z_atan2(Rcomplex *r, Rcomplex *csn, Rcomplex *ccs)





FILE Hyperbolic.Rd STARTS
\name{Hyperbolic}
\title{Hyperbolic Functions}
\usage{
cosh(x)
sinh(x)
tanh(x)
acosh(x)
asinh(x)
atanh(x)
}
\alias{cosh}
\alias{sinh}
\alias{tanh}
\alias{acosh}
\alias{asinh}
\alias{atanh}
\description{
   These functions give the obvious hyperbolic functions.  They
   respectively compute the hyperbolic cosine, sine, tangent, and their
   inverses, arc-cosine, arc-sine, arc-tangent (or \dQuote{\emph{area 
cosine}},
   etc).
}
\arguments{
   \item{x}{a numeric or complex vector}
}
\details{
   These are generic functions: methods can be defined for them
   individually or via the \code{\link{Math}} group generic.

   Branch cuts are consistent with the inverse trigonometric functions
   \code{asin()} et seq, and agree with those defined in Abramowitz and
   Stegun, figure 4.7, page 86.

}
\seealso{
   The trigonometric functions, \code{\link{cos}}, \code{\link{sin}},
   \code{\link{tan}}, and their inverses
   \code{\link{acos}}, \code{\link{asin}}, \code{\link{atan}}.

   The logistic distribution function \code{\link{plogis}} is a shifted
   version of \code{tanh()}.
}
\keyword{math}
FILE Hyperbolic.Rd ENDS



File Trig.Rd STARTS
\name{Trig}
\alias{Trig}
\alias{cos}
\alias{sin}
\alias{tan}
\alias{acos}
\alias{asin}
\alias{atan}
\alias{atan2}
\title{Trigonometric Functions}
\description{
   These functions give the obvious trigonometric functions.  They
   respectively compute the cosine, sine, tangent, arc-cosine, arc-sine,
   arc-tangent, and the two-argument arc-tangent.
}
\usage{
cos(x)
sin(x)
tan(x)
acos(x)
asin(x)
atan(x)
atan2(y, x)
}
\arguments{
   \item{x, y}{numeric or complex vector}
}
\details{
   The arc-tangent of two arguments \code{atan2(y,x)} returns the angle
   between the x-axis and the vector fro

[Rd] branch cuts in atan(), asin(), acos() (PR#7871)

2005-05-16 Thread r . hankin
The complex versions of atain(), asin(), acos() are not documented 
except in the source
code.

The branch cuts that are implemented do not tally with those of 
Abramowitz and Stegun,
figure 4.4, p79.  The current branch cuts appear to be inherited from 
sqrt() and log(),
in a non-obvious manner.


Also, the current behaviour of atan() leads to the following:


R> tan(atan(2))
[1] 2
R> tan(atan(2+0i))
[1] -0.5+0i

(I would expect  a value close to 2 in both these)

R> atan(2)
[1] 1.107149
R> atan(2+0i)
[1] -0.4636476+0i
 >

(I would expect these two expressions to evaluate to numbers with small 
absolute
difference)


 > atan(1.0001+0i)
[1] -0.7853482+0i
 > atan(0.+0i)
[1] 0.7853482+0i
 >

(I would expect these two expressions to evaluate to numbers with small 
absolute
difference)


The following patch to src/main/complex.c
appears to implement the desired behaviour:<


octopus:~/downloads/R-2.1.0/src/main% diff -c complex.c 
/Users/rksh/unix/rksh/complex.c
*** complex.c   Mon Apr 18 22:30:00 2005
--- /Users/rksh/unix/rksh/complex.c Mon May 16 14:33:15 2005
***
*** 367,372 
--- 367,376 
   bet = t1 - t2;
   r->r = asin(bet);
   r->i = log(alpha + sqrt(alpha*alpha - 1));
+
+ if(y<0){
+   r->i *= -1;
+   }
   }

   static void z_acos(Rcomplex *r, Rcomplex *z)
***
*** 388,393 
--- 392,404 
   r->r = 0.5 * atan(2 * x / ( 1 - x * x - y * y));
   r->i = 0.25 * log((x * x + (y + 1) * (y + 1)) /
   (x * x + (y - 1) * (y - 1)));
+
+ if((x*x+y*y > 1) && x>0){
+   r->r += M_PI_2;
+   }
+ if((x*x+y*y > 1) && x<0){
+   r->r -= M_PI_2;
+   }
   }

   static void z_atan2(Rcomplex *r, Rcomplex *csn, Rcomplex *ccs)
octopus:~/downloads/R-2.1.0/src/main%








--
Robin Hankin
Uncertainty Analyst
National Oceanography Centre, Southampton
European Way, Southampton SO14 3ZH, UK
  tel  023-8059-7743

__
R-devel@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel