complex numbers multiplication in double precision

49 次查看（过去 30 天）
FastCar
FastCar 2017-7-1
评论： James Tursa 2017-7-3
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

请先登录，再回答此问题。

回答（3 个）

Walter Roberson
Walter Roberson 2017-7-1
编辑：Walter Roberson 2017-7-1
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.
  2 个评论
FastCar
FastCar 2017-7-1
Ok, I have missed an important detail. The information I have is the phase. So the first equation uses cos and sin to get the real part and the imaginary part. I can tell you that doing a 10^7 addition on the multiplications the two methods give a different result. My idea was that the second one is more precise because does less computations and thus less approximations, but this is far from being a demonstration that it is more accurate. Any idea how to demonstrate that it should be more accurate?

请先登录，再进行评论。

James Tursa
James Tursa 2017-7-2
编辑：James Tursa 2017-7-2
"... 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
FastCar 2017-7-2
编辑：FastCar 2017-7-2
I try to clarify my question: the phase is given by the distance of a satellite from the Earth multiplied by
4 x pi / lambda
The distance is big (around 40000 km) and thus gives a huge phase. So far I used these expressions:
distance = norm ( position vector )
phi = ( 4 x pi x distance ) / lambda
method one
c1= cos(phi1) + j x sin(phi1)
c2= cos(phi2) + j x sin(phi2)
c1 x c2* = cos(phi1) x cos(phi2) + sin(phi1) x sin(phi2) + j x (sin(phi1) x cos (phi2) - cos(phi2) x sin(c2)
method 2
c1 x c2* = cos(phi1-phi2) + j x sin(phi1-phi2)
For what I know every computation in matlab should be done in double precision and without the symbolic packed pi should have only 15 figures (don't blame me if I am wrong). Mind that 15 figures (and thus 14 sure figures) is a good output in my point of view.
I hope I have clarified the question.
I haven't answered to the single question because I think both of you needed the full picture.
  1 个评论
James Tursa
James Tursa 2017-7-3
I don't quite understand your physical application yet, but keep in mind what I have already written in my Answer. Because phi is so large, you have effectively already lost 9 significant digits from your answer. So the c1 x c2* result will only have about 6 meaningful digits of accuracy at most.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Logical 的更多信息

标签

产品

Community Treasure Hunt

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

Start Hunting!

Translated by