R-beta: (0+0i)^2

[Martin Maechler]

>>>>> "Osman" == U-E59264-Osman F Buyukisik <absd00t@c1186.ae.ge.com> writes:

Osman> I am not a "c" programmer but this is very simple and seems to be
Osman> working (just a check for 0): in src/main/complex.c 
Osman> ...
Osman> static void complex_pow(complex *r, complex *a, complex *b)
Osman> {
Osman>     double logr, logi, x, y, mag;
Osman>     if ((mag = hypot(a->r, a->i)) == 0.0) {
Osman>         r->r = r->i = 0.0;
Osman>     }
Osman>     else {
Osman>     logr = log(mag );
Osman>     logi = atan2(a->i, a->r);
Osman>     x = exp( logr * b->r - logi * b->i );
Osman>     y = logr * b->i + logi * b->r;
Osman>     r->r = x * cos(y);
Osman>     r->i = x * sin(y);
Osman>     }
Osman> }

It's not *that* simple:  
We want
	0^ 0   => 1	(since we have that for reals)
	0^ -1  => Inf	(a^-1 = 1/a)
	0^ 1i  => NaN

I've committed a longer patch for 0.63,
which also leads to   1i^(-5:5)  be 100% accurate
but the following fixes the bug (and more):

--- complex.c~	Mon Aug 31 10:32:29 1998
+++ complex.c	Fri Sep 18 16:31:19 1998
@@ -91,6 +91,10 @@
 {
     double logr, logi, x, y;
 
+    if(b->i == 0.) {/* be fast (and more accurate)*/
+	if(b->r == 1.) { /* a^1 */ r->r = a->r; r->i = a->i; return;}
+	if(a->i == 0.) { r->r = pow(a->r, b->r); r->i = 0.; return;}
+    }
     logr = log(hypot(a->r, a->i) );
     logi = atan2(a->i, a->r);
 

Martin Maechler <maechler@stat.math.ethz.ch>	http://stat.ethz.ch/~maechler/
Seminar fuer Statistik, ETH-Zentrum SOL G1;	Sonneggstr.33
ETH (Federal Inst. Technology)	8092 Zurich	SWITZERLAND
phone: x-41-1-632-3408		fax: ...-1086			<><
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-devel-request@stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._