Opened 7 years ago

Closed 4 years ago

#12807 closed defect (fixed)

Taking the real part of a sum of exponentials with imaginary exponents gives wrong result

Reported by: inaki Owned by: burcin
Priority: major Milestone: sage-6.7
Component: symbolics Keywords:
Cc: Merged in:
Authors: Ralf Stephan Reviewers: Travis Scrimshaw
Report Upstream: Fixed upstream, in a later stable release. Work issues:
Branch: 410bab7 (Commits) Commit: 410bab7e7016f06220e1a34335abbaa8f5e3e876
Dependencies: #18362 Stopgaps:

Description

In sage 4.8 the following computation gives 1/2, while the right result is -3/2:

sage: a = exp(i*2*pi/3)
sage: b = exp(i*pi/3)
sage: c = 2*a-b
sage: c.real()
1/2

Doing the computation numerically gives the right result:

sage: N(c).real()
-1.50000000000000

and asking maxima directly also works:

sage: (c.maxima_methods().rectform()).real()
-3/2

This is sage 4.8 running on Fedora 16, linux 3.1.0-7.fc16.x86_64 SMP.

Change History (11)

comment:1 Changed 7 years ago by burcin

Thanks for the report! This is a bug in Pynac. I opened a ticket here:

http://hg.pynac.org/pynac/issue/15/wrong-result-from-real-part

comment:2 Changed 6 years ago by jdemeyer

  • Milestone changed from sage-5.11 to sage-5.12

comment:3 Changed 5 years ago by vbraun_spam

  • Milestone changed from sage-6.1 to sage-6.2

comment:4 Changed 5 years ago by vbraun_spam

  • Milestone changed from sage-6.2 to sage-6.3

comment:5 Changed 5 years ago by vbraun_spam

  • Milestone changed from sage-6.3 to sage-6.4

comment:6 Changed 5 years ago by kcrisman

This gets really weird, actually.

sage: a.real()
-1/2
sage: b.real()
1/2
sage: (a-b).real() # good
-1
sage: (a-2*b).real() # good
-3/2
sage: (a-3*b).real()  # ho-hum
-2
sage: (2*a-b).real() # yikes
1/2
sage: (3*a-b).real()  # what?
3/4
sage: (4*a-b).real()  # consistent...
1
sage: (5*a-b).real()  # ?!?!?
5/4
sage: (-a-b).real()  # seeing but not believing
-1/4
sage: (-2*a-b).real()
-1/2
sage: (-3*a-b).real()
-3/4
sage: (2*a).real()  # okay
-1
sage: (2*a+b).real()  #  okay
-1/2
sage: (3*a+b).real()
-1
sage: (4*a+b).real()
-3/2

This is consistent with Pynac somehow thinking that a.real()==1/4 when preceded by a positive constant and -1/4 when given a negative constant. But only when part of a sum including -b, with +b it's fine.

I tried to figure out what it was but so far no luck.

comment:7 Changed 4 years ago by rws

  • Dependencies set to pynac-0.3.7
  • Milestone changed from sage-6.4 to sage-6.7
  • Report Upstream changed from N/A to Fixed upstream, in a later stable release.

Actually and surprisingly what's happening is that the simplified values of the exp expressions (1/2 and -1/2) and their coefficients are all multiplied together on final evaluation:

sage: (4*exp(i*pi/3)-4*exp(i*2*pi/3)).real_part()
4
sage: (6*exp(i*pi/3)-6*exp(i*2*pi/3)).real_part()
9

This is evident in the GiNaC source, this line(https://github.com/pynac/pynac/blob/master/ginac/add.cpp#L432):

    oc = oc.mul(ex_to<numeric>(j->rest)).mul(ex_to<numeric>(j->coeff));

is so rarely executed that e.g. only three doctests in symbolic/ touch it but don't trigger the bug because the expressions are too short.

The fix will be in Pynac-0.3.7.

comment:8 Changed 4 years ago by rws

  • Branch set to u/rws/taking_the_real_part_of_a_sum_of_exponentials_with_imaginary_exponents_gives_wrong_result

comment:9 Changed 4 years ago by rws

  • Authors set to Ralf Stephan
  • Commit set to 410bab7e7016f06220e1a34335abbaa8f5e3e876
  • Dependencies changed from pynac-0.3.7 to #18362
  • Status changed from new to needs_review

New commits:

23892cb18362: upgrade to pynac-0.3.7; fix one changed doctest result
410bab712807: doctest

comment:10 Changed 4 years ago by tscrim

  • Reviewers set to Travis Scrimshaw
  • Status changed from needs_review to positive_review

LGTM.

comment:11 Changed 4 years ago by vbraun

  • Branch changed from u/rws/taking_the_real_part_of_a_sum_of_exponentials_with_imaginary_exponents_gives_wrong_result to 410bab7e7016f06220e1a34335abbaa8f5e3e876
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.