🚀 go-pugleaf

In article <CBws90.K8B@cid.aes.doe.ca>, ojastej@wfh7546.ice.ncr.ca (James Ojaste) writes:
|> In article <1993Aug16.194800.301@cs.uah.edu>, vchagant@ (Venkata Chaganti) writes:
|> |> I would like to know is there any public domain routines or library
|> |> to do Complex arthimetic in C language.
|>
|> |> I appreciate for pointers to ftp sites.
|> |> 
|> |> The main question is 
|> |> Z1=a1+ib1
|> |> Z2=a2+ib2
|> |> 
|> |> as and bs are real numbers. some times fractions too.
|> 
|> Just to be picky, all fractions are real numbers (assuming that the i is not
|> included except as a coefficient)...
|> 
|> I don't know about C (you certainly cannot use that syntax in C - you'd have
|> to use complex_add (a1, b1); ), but in C++ it shouldn't be too hard - just
|> overload the + operator and define a couple of functions...  If you don't 
|> know C++ - sorry (I'm just learning, and so I'm not quite capable of posting
|> source).
|> 
|> -- 
|> /*  0F 90 3E 44 F9 13 E7 CC	(jojaste@descartes.uwaterloo.ca)
|> **  10 20 44 C8 42 21 08 98	 
|> **  38 40 F9 50 87 C2 1F 30	"Woof bloody woof" - Gaspode the Wonder Dog	 
|> **  40 81 12 61 08 84 22 00	     (Terry Pratchett, Moving Pictures)
|> */  F9 FA 24 42 11 3F 44 C0
|> 

I found a book on C++ by Scott Robert Ladd helpful, it has a PC disk with
a C++ complex class.

Here is Ladd's equation:

LComplex pow(const LComplex & c, const LComplex & power)
    {
[...]
    result = exp(log(c) * power);
[...]
where exp and log are complex functions defined else where.  Then he
defines sqrt in terms of power.

LComplex sqrt(const LComplex & c)
    {
    return pow(c,LComplex(0.5,0.0));
}

Now I object to using exp and log (complex) to do a square root!!!  I
replace this with

LComplex sqrt(const LComplex & c)
{
  long double a=abs(c);
  if(c.Imag>0.0)
    return LComplex(sqrtl( (a+c.Real)*0.5 ),sqrtl( (a-c.Real)*0.5 ));
  else
    return LComplex(sqrtl( (a+c.Real)*0.5 ),-sqrtl( (a-c.Real)*0.5 ));
}

but I suspect this is inaccurate for some range of arg(c), near 0 and pi.  
I use 4 multiplies, 3 adds/subtracts and 3 square roots.  Do any of you 
good folk have an improvement for this?

Warmly,
Jay 'Not tonight honey, I'm on the edge of an hyperbolic component' Hill
-- 
{ hilljr@jupiter.saic.com } begin writeln(3*ln(640320)/sqrt(163):17:15) end.
void main(){double sqrt(), y=1/sqrt(2.), a=.5, m=1,z; int n=0;
for(;m*=2,z=sqrt(1-y*y),y=(1-z)/(1+z),a=a*(1+y)*(1+y)-m*y,n<4;n++);
printf("%17.15lf\n",1/a);}