86 views (last 30 days)

Hello,

I have to multiply couple of complex numbers and then I have to add all the product.

Let's say the two complex numbers (with unitary module) are c1 and c2. What I need is prod=c1 x c2* (but I can substitute c2 with its conjiugate in our example). c1 has phase phi1, c2 has phase phi2. I can do it in chartesian coordinates

c1 x c2* = real(c1) x real (c2) + imag(c1) x imag(c2) + j x (imag(c1) x real (c2) - real(c1) x imag(c2)

(if I have used the right formula)

or I can use

c1 x c2* = cos(phi1-phi2) + j x sin(phi1-phi2)

Mind that the phases are of the order of magnitude of 10^10 radians.

The second way is of course faster, but my question is which one is the more precise given that the computation is done in double precision?

Mind also that I have to sum 10^7 products and I can tell you that there is a difference, a not neglectable difference.

Many thanks in advance

Carlo

Walter Roberson
on 1 Jul 2017

Edited: Walter Roberson
on 1 Jul 2017

On my system, even excluding the time it takes to calculate the angle, the cos/sin version is about 1.8 to 2.4 times slower than the other version. This is to be expected as cos and sin require a many more calculations. There are hardware instructions, yes, but they are slower than plain multiply or addition.

The second version appears to have marginally better numerical properties eps(1) compared to 2*eps(1). That is not really a meaningful difference compared to the rough-off error in the rest of the calculations.

James Tursa
on 2 Jul 2017

Edited: James Tursa
on 2 Jul 2017

"... Mind that the phases are of the order of magnitude of 10^10 radians. ..."

When you have angles that large, I have to start questioning where they came from and how exactly you are using them in your code. Realize that 10^10 is more than 9 orders of magnitude larger than 2pi, so this will come into play in the "accuracy" you expect for calculations. E.g., a naive example:

>> a = 10^10

a =

1.000000000000000e+010

>> cos(a)+sin(a)*i

ans =

0.873119622676856 - 0.487506025087511i

>> cos(mod(a,2*pi))+sin(mod(a,2*pi))*i

ans =

0.873119809656583 - 0.487505690208075i

So you can see a difference in the 6th digit. What is going on? It turns out that the trig functions have a very sophisticated method of dealing with the argument. I don't know the exact method that is used behind the scenes, but effectively it treats the IEEE double precision argument as if it has the same absolute precision as a IEEE double precision pi. E.g.,

>> va = vpa(a)

va =

10000000000.0

>> vp = vpa(pi)

vp =

3.1415926535897932384626433832795

>> n = floor(va/(2*vp))

n =

1591549430

>> mod(a,2*pi)

ans =

5.773954618557402

>> double(va-n*2*vp)

ans =

5.773954235013852

>> cos(ans)+sin(ans)*i

ans =

0.873119622676856 - 0.487506025087511i

So you can see by treating the angle as if it had the same absolute precision as a double precision pi value we can do the mod in such a way that we can uncover how cos() and sin() are behaving behind the scenes. It is effectively a different angle that is used compared to the naive mod(a,2*pi) angle result.

Bottom line is that with 10^10 radian angles you need to be very careful that your problem has real meaning at these values, and then be extra careful how you are doing calculations with them to come up with meaningful results. For your phi1-phi2 method, realize that the cos() and sin() functions are going to treat that phi1-phi2 double precision result as if it had absolute precision the same as a double precision pi. This effect may be creeping into your results.

FastCar
on 2 Jul 2017

Edited: FastCar
on 2 Jul 2017

James Tursa
on 3 Jul 2017

Opportunities for recent engineering grads.

Apply TodayFind the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!
## 0 Comments

Sign in to comment.