views:

200

answers:

2

Hi all,

I'm writing some coordinate transformations (more specifically the Joukoswky Transform, Wikipedia Joukowsky Transform), and I'm interested in performance, but of course precision. I'm trying to do the coordinate transformations in two ways:

1) Calculating the real and complex parts in separate, using double precision, as below:

double r2 = chi.x*chi.x + chi.y*chi.y;

//double sq = pow(r2,-0.5*n) + pow(r2,0.5*n); //slow!!!
double sq = sqrt(r2); //way faster!
double co = cos(atan2(chi.y,chi.x));
double si = sin(atan2(chi.y,chi.x));

Z.x = 0.5*(co*sq + co/sq);
Z.y = 0.5*si*sq;

where chi and Z are simple structures with double x and y as members.

2) Using complex :

Z = 0.5 * (chi + (1.0 / chi));

Where Z and chi are complex . There interesting part is that indeed the case 1) is faster (about 20%), but the precision is bad, giving error in the third decimal number after the comma after the inverse transform, while the complex gives back the exact number. So, the problem is on the cos(atan2), sin(atan2)? But if it is, how the complex handles that?

EDIT: Just figured out that this was not exactly the question that I had in mind. I have to do the general transformation, as

Z = 1/2*(chi^n + (1/chi)^n), and so far the code above was the way I've figured to do it. More precisely,

    double sq = pow(sqrt(r2),n); //way faster!
double co = cos(n*atan2(chi.y,chi.x));
double si = sin(n*atan2(chi.y,chi.x));

Z.x = 0.5*(co*sq + co/sq);
Z.y = 0.5*(si*sq - sq/si);

Also correcting the bug on Z.y.

Thanks!

+3  A: 

I think that in 1) it should be

Z.y = 0.5*(si*sq - si/sq);

If you want really good performance you may want to go back to first principles and observe that

1/(a+ib) = (a-ib)/(a*a+b*b)

No sqrt(), atan2() or cos() or sin().

brainjam
Thanks for the answer! It was really a stupid bug on the si. I use this notation because I have to raise this expression to power n, asZ = 0.5(chi^n + (1/chi)^n)so it becomessq = pow(r2,n)co = cos(n*atan2(chi.y,chi.x))and so far this is the only way I've found to do this.
Ivan
+4  A: 

Given r = sqrt(x*x+y*y):

cos(atan2(y,x)) == x/r
sin(atan2(y,x)) == y/r

Calculating it this way should be more accurate and faster.

When you plug these values into the formulas for Z.x and Z.y, the square root will cancel out as well, so you'll be left with only basic arithmetic operations.

interjay
Just what I was thinking.
Mike Dunlavey
In the explicit case, yes, but not in the general case z = 1/5(chi^n + (1/chi)^n). I've should be more precise, sorry about that - I've corrected the original post. Thanks for your help!
Ivan