Issues with upgrading from SymPy 0.7.1 to 0.7.2, 0.7.3, or 0.7.4

[petarmaric]

Hi, I'm the author of beam_integrals, a project heavily reliant on SymPy. Our extensive test suite is demonstrating weird issues when upgrading from SymPy 0.7.1 to 0.7.2 or 0.7.3.

What's causing these issues?

Testing environments

Ubuntu 10.04 32bit @ our Jenkins CI

SymPy version	gmpy version	SYMPY USE CACHE	Tests pass?	Test run time	Details
0.7.1	not installed	undefined	YES	17840.769s	http://ci.petarmaric.com/job/beam_integrals/49/
0.7.2	not installed	undefined	NO, stopped	~15 hours	http://ci.petarmaric.com/job/beam_integrals/50/
0.7.1	1.10-1 via Ubuntu package	undefined	YES	15928.379s	http://ci.petarmaric.com/job/beam_integrals/51/

Notes:

1. These tests were performed before I was made aware of the memory leak caused by the SymPy internal caching issues, soon afterwards I've patched our test suite to avoid the memory leak
2. http://ci.petarmaric.com/job/beam_integrals/50/ has been stopped manually as it ate our entire RAM and most of the swap (memory leak due to SymPy internal caching issues, see previous note)

Ubuntu 12.04 64bit @ DigitalOcean Cloud

SymPy version	gmpy version	SYMPY USE CACHE	Tests pass?	Test run time	Test error details
0.7.1	not installed	no	YES	10068.788s
0.7.1	1.14-3build1 via Ubuntu package	no	YES	7787.417s
0.7.2	not installed	no	YES	8278.762s
0.7.2	1.14-3build1 via Ubuntu package	no	NO	368.534s	https://gist.github.com/petarmaric/6346863
0.7.3	not installed	no	NO	47966.926s	https://gist.github.com/petarmaric/6350717
0.7.3	1.14-3build1 via Ubuntu package	no	NO	763.135s	https://gist.github.com/petarmaric/6351831
0.7.3	2.0.2 via source code install	no	NO	810.033s	https://gist.github.com/petarmaric/6352787

Notes:

1. 3 droplets were created (sympy-0-7-1, sympy-0-7-2 and sympy-0-7-3) which were setup identically, with the exception of pinning the SymPy version, by editing beam_integrals/requirements.txt locally
2. gmpy2 was installed from source, based on the official instructions

How I've tested beam_integrals on DigitalOcean Cloud

1. Log into DigitalOcean Control Panel at https://www.digitalocean.com/login

2. Create a new Ubuntu 12.04 64bit Droplet with 2048 MB of RAM

   • Wait until the server is up and login as root

3. Update the system:

   aptitude update
   aptitude safe-upgrade
   

4. Install various system monitoring utilities:

   aptitude install htop iftop
   

5. Install the Ubuntu packages required to build and test beam_integrals:

   aptitude install build-essential python-dev mercurial python-pip
   

6. Install beam_integrals package in the "development" mode, along with the rest of its requirements:

   hg clone https://bitbucket.org/petar/beam_integrals
   cd beam_integrals
   python setup.py develop
   pip install -r requirements-test.txt
   

7. [OPTIONAL] It's strongly recommended to install gmpy. Without it code will still run correctly, but much slower at high precision. Install it as an Ubuntu package, to avoid compiling it manually:

   aptitude install python-gmpy
   

8. Reboot the server, to apply kernel updates and clear disk cache from RAM (to start fresh)

   • Wait until the server is up and login as root

9. Run unit tests with SymPy internal caching disabled (see https://bitbucket.org/petar/beam_integrals/commits/5bcca1e for details):

   cd beam_integrals
   SYMPY_USE_CACHE=no python setup.py nosetests
   

[asmeurer]

We are planning on replacing the global cache with more localized caches, which should fix the memory issues.

In the meanwhile, I would recommend clearing the cache periodically (such as between each test) rather than disabling it outright, as some functionality is orders of magnitude slower without it.

[asmeurer]

As to the test failures, I'm not sure what to make of them. What sort of operations is it performing?

[petarmaric]

Re: Caching

Yeah, the global cache is a bit of a bother right now, it's great we're moving towards local/function level caching.

I wanted to avoid disabling the global cache, however I was unable to find an easy way to do that with nose after each test - that is without monkey-patching nose/unittest itself (too much magic == debugging nightmare), or polluting the test code with a new code snippet (requiring us to remember to add it for each and every test teardown). However I've added an option to run our tests without disabling the global cache, by setting the environment variable IGNORE_SYMPY_MEMORY_LEAK to yes.

As you can see from our build time trend the testing performance didn't suffer much from disabling the SymPy internal cache (done at build #56).

[asmeurer]

Some operations suffer performance hits from disabling the cache more than others. If you want to get an idea of what can happen, try running the SymPy tests with and without the cache (use bin/test -C).

[asmeurer]

See #2448 about caching.

[petarmaric]

Still more weird issues, now with SymPy 0.7.4.1.

The correct behavior, with SymPy 0.7.1 (no gmpy/gmpy2 installed):

>>> from beam_integrals.characteristic_equation_solvers import find_best_root
d:\_fax\rnd\beam_integrals\beam_integrals\__init__.py:45: PerformanceWarning: gmpy is not available for speedups (code will still run correctly, but much slower at high precision)
  category=PerformanceWarning
>>> from beam_integrals.beam_types import BaseBeamType

>>> beam_type = BaseBeamType.coerce(2)

# With `decimal_precision=100` the error is at approx. 6.24e-87
>>> find_best_root(beam_type, mode=10, use_cache=False, include_error=True, decimal_precision=100)
(32.98672286269281956154719686483756311040117568283227892102769304821772658683437565346748656460046447,
 6.239916657399792542797394081821281856785313748862743716779932350352203130542726844127910805693946236e-87)

# With `decimal_precision=200` the error has, as expected, substantially decreased (to approx. 2.21e-187)
>>> find_best_root(beam_type, mode=10, use_cache=False, include_error=True, decimal_precision=200)
(32.986722862692819561547196864837563110401175682832278921027693048217726586834375653467486564600464466799729091688477730536938182713305381298734791755624232559116615389154222342533608728820773693414151,
 2.2060364672001921377137344126949511259376484093355774469862228563431583786670362636014501447438971581007617378278309520871891472197878372228096180832647648284978834368664148014992478935353203898044366e-187)

# With `decimal_precision=400` the error has, as expected, substantially decreased (to approx. 1.60e-386)
>>> find_best_root(beam_type, mode=10, use_cache=False, include_error=True, decimal_precision=400)
(32.98672286269281956154719686483756311040117568283227892102769304821772658683437565346748656460046446679972909168847773053693818271330538129873479175562423255911661538915422234253360872882077369341415082291083560682483998739270166380947899449185251590416935767850811732133893906742954901472996336387571722543432233907421403398716945791523183080542174552979588066467890474734232966702033580464543028702,
 1.598114239518829991502353652997735406710979438725371440914289990920001389094212796806840915477033217673894589544405863334100539159882416670486196556368495660065616376920381544471058648907971103611824394940967457940153690105639247114047090691687802601044838845951151379553178701526533808969875159425009134471427915194837844156055107347237831434868342437079182927310965982449386203711667325766283501291e-386)

The incorrect behavior, with SymPy 0.7.4.1 (no gmpy/gmpy2 installed):

>>> from beam_integrals.characteristic_equation_solvers import find_best_root
d:\_fax\rnd\beam_integrals\beam_integrals\__init__.py:45: PerformanceWarning: gmpy is not available for speedups (code will still run correctly, but much slower at high precision)
  category=PerformanceWarning
>>> from beam_integrals.beam_types import BaseBeamType

>>> beam_type = BaseBeamType.coerce(2)

# UNEXPECTED: With `decimal_precision=100` the error is at approx. 0.453, which is way to large for this working precision
>>> find_best_root(beam_type, mode=10, use_cache=False, include_error=True, decimal_precision=100)
(32.98672286269281528348074061796069145202636718750000000000000000000000000000000000000000000000000000,
 0.4530741103747513228739408406112422142400775469596284469083916899851644037039050342437212276106206340)

# UNEXPECTED: With `decimal_precision=200` the error has remained at pretty much the same level
>>> find_best_root(beam_type, mode=10, use_cache=False, include_error=True, decimal_precision=200)
(32.986722862692815283480740617960691452026367187500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000,
 0.45307411037475132287394084061124221424007754695962844690839168998516440370390503424372122761062063403665594391208613321068949713141126978226162659393402874350758617387952575737565964184528698448198416)

# SERIOUS ERROR: With `decimal_precision=400` the resulting root is nonsense
>>> find_best_root(beam_type, mode=10, use_cache=False, include_error=True, decimal_precision=400)
(0.0, 0.e-808)

[asmeurer]

Can you check the git master?

[petarmaric]

Sure, just tested the above code with SymPy master and I've got exactly the same behaviour and numbers as with SymPy 0.7.4.1 (no gmpy/gmpy2 installed).

[asmeurer]

So is the behavior different if gmpy is installed?

[petarmaric]

When running the above test snippet with SymPy 0.7.1 I get exactly the same behaviour and numbers in both cases: without gmpy and when gmpy (gmpy-1.15.win-amd64-py2.7.exe) is installed.

With SymPy 0.7.4.1 and gmpy (gmpy-1.15.win-amd64-py2.7.exe) I get:

>>> from beam_integrals.characteristic_equation_solvers import find_best_root
>>> from beam_integrals.beam_types import BaseBeamType

>>> beam_type = BaseBeamType.coerce(2)

# UNEXPECTED: With `decimal_precision=100` the error is at approx. 0.453, which is way to large for this working precision
>>> find_best_root(beam_type, mode=10, use_cache=False, include_error=True, decimal_precision=100)
(32.98672286269281528348074061796069145202636718750000000000000000000000000000000000000000000000000000,
 0.4530741103747513228739408406112422142400775469596284469083916899851644037039050342437212276106206340)

# UNEXPECTED: With `decimal_precision=200` the error has remained at pretty much the same level
>>> find_best_root(beam_type, mode=10, use_cache=False, include_error=True, decimal_precision=200)
(32.986722862692815283480740617960691452026367187500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000,
 0.45307411037475132287394084061124221424007754695962844690839168998516440370390503424372122761062063403665594391208613321068949713141126978226162659393402874350758617387952575737565964184528698448198416)

# SERIOUS ERROR: With `decimal_precision=400` the resulting root is nonsense
>>> find_best_root(beam_type, mode=10, use_cache=False, include_error=True, decimal_precision=400)
(0.0, 0.e-808)

which is exactly the same behaviour and numbers as with SymPy 0.7.4.1 when no gmpy/gmpy2 is installed.

Seems some kind of a regression has been introduced in SymPy after 0.7.1, as I'm also getting hit by #2724 when running SYMPY_USE_CACHE=no python setup.py nosetests for the beam_integrals project (SymPy 0.7.4.1 with gmpy-1.15.win-amd64-py2.7.exe):

======================================================================
ERROR: tests.characteristic_equation_solvers.test_cache.TestBestRootsCache.test_mode_list_ordering
----------------------------------------------------------------------
Traceback (most recent call last):
  File "C:\_dev\Python27\lib\site-packages\nose\case.py", line 197, in runTest
    self.test(*self.arg)
  File "D:\_FAX\RnD\beam_integrals\tests\characteristic_equation_solvers\test_cache.py", line 88, in test_mode_list_orde
ring
    eq_(mode_list, sorted(mode_list))
  File "C:\_dev\Python27\lib\site-packages\sympy\core\numbers.py", line 889, in __lt__
    other = other.evalf()
  File "C:\_dev\Python27\lib\site-packages\sympy\core\evalf.py", line 1300, in evalf
    re = C.Float._new(re, p)
  File "C:\_dev\Python27\lib\site-packages\sympy\core\numbers.py", line 655, in _new
    obj._mpf_ = mpf_norm(_mpf_, _prec)
  File "C:\_dev\Python27\lib\site-packages\sympy\core\numbers.py", line 55, in mpf_norm
    rv = mpf_normalize(sign, man, expt, bc, prec, rnd)
TypeError: argument is not an mpz

----------------------------------------------------------------------

[asmeurer]

Are any of these issues resolved if you enable the cache?

[petarmaric]

Not really, when running IGNORE_SYMPY_MEMORY_LEAK=yes python setup.py nosetests -x for the beam_integrals project (SymPy 0.7.4.1 with gmpy-1.15.win-amd64-py2.7.exe) I get:

======================================================================
ERROR: tests.characteristic_equation_solvers.test_cache.TestBestRootsCache.test_mode_list_ordering
----------------------------------------------------------------------
Traceback (most recent call last):
  File "C:\_dev\Python27\lib\site-packages\nose\case.py", line 197, in runTest
    self.test(*self.arg)
  File "D:\_FAX\RnD\beam_integrals\tests\characteristic_equation_solvers\test_cache.py", line 88, in test_mode_list_orde
ring
    eq_(mode_list, sorted(mode_list))
  File "C:\_dev\Python27\lib\site-packages\sympy\core\numbers.py", line 889, in __lt__
    other = other.evalf()
  File "C:\_dev\Python27\lib\site-packages\sympy\core\evalf.py", line 1300, in evalf
    re = C.Float._new(re, p)
  File "C:\_dev\Python27\lib\site-packages\sympy\core\numbers.py", line 655, in _new
    obj._mpf_ = mpf_norm(_mpf_, _prec)
  File "C:\_dev\Python27\lib\site-packages\sympy\core\numbers.py", line 55, in mpf_norm
    rv = mpf_normalize(sign, man, expt, bc, prec, rnd)
TypeError: argument is not an mpz

----------------------------------------------------------------------

Judging from the tracebacks it seems Float.__lt__ or some code it uses is causing these issues. Seems like a regression to me.

Would love to help in drilling down to the root cause I'm very busy with writing my PhD thesis for the next month or so.

[asmeurer]

Can you print what mode_list is?

[hargup]

@petarmaric @asmeurer is this issue still valid?

[petarmaric]

Very much so. you can try it out youself. I'll have more time to analyze the underlying issue once the semester ends.

[asmeurer]

If you could narrow the issue down to a simple example that could be reproduced without the use of beam_integrals, that would help a lot.

[petarmaric]

I'm sick (yay?), so I finally have some time to try and figure out which SymPy commit introduced a regression that caused beam_integrals tests to fail.

Please note, this is all done without gmpy installed and without disabling the SymPy internal cache.
If the variables:

• gmpy_installed (no, 1.x via Ubuntu package, 2.x via source code install) and
• sympy_cache (disabled, enabled but without the memory leak warning)

were introduced it would result in way too many testing variations for me to handle right now.

How to test beam_integrals on DigitalOcean Cloud

1. Log into DigitalOcean Control Panel at https://www.digitalocean.com/login

2. Create a new Ubuntu 14.04 64bit Droplet with 2048 MB of RAM

   • Wait until the server is up and login as root // I know, I know... It's easy and I'm lazy and I'm very sick right now

3. Update the system:

   aptitude update
   aptitude safe-upgrade
   

4. Install the Ubuntu packages required to build and test beam_integrals and SymPy:

   aptitude install build-essential python-dev mercurial git python-pip
   

5. Install beam_integrals package in the "development" mode, along with the rest of its requirements:

   hg clone https://bitbucket.org/petar/beam_integrals ~/beam_integrals
   cd ~/beam_integrals
   python setup.py develop
   pip install -r requirements-test.txt
   

6. Remove the SymPy 0.7.1 package that beam_integrals installed as its dependency:

   pip uninstall sympy
   

7. Clone the SymPy repository:

   git clone git://github.com/sympy/sympy.git ~/sympy
   

8. Switch over to SymPy 0.7.1 and install it in the "development" mode:

   cd ~/sympy
   git checkout sympy-0.7.1
   pip install -e .
   

9. Make sure SymPy 0.7.1 works correctly:

   cd ~/sympy
   ./bin/test
   ./bin/doctest
   

10. Reboot the server, to apply kernel updates and clear disk cache from RAM (to start fresh)

    • Wait until the server is up and login as root

11. Run the entire beam_integrals test suite for SymPy 0.7.1 (merely a sanity test), without disabling the SymPy internal cache
    (see https://bitbucket.org/petar/beam_integrals/commits/5bcca1e and Issues with upgrading from SymPy 0.7.1 to 0.7.2, 0.7.3, or 0.7.4 #2439 (comment) for details):

    cd ~/beam_integrals
    IGNORE_SYMPY_MEMORY_LEAK=yes python setup.py nosetests
    

Using beam_integrals to find a regression in SymPy > 0.7.1

From now on we'll only run tests within tests.characteristic_equation_solvers.test_solvers
(as the full beam_integrals test suite takes almost 2 hours to complete), in stop-on-failure mode and with generally lowered verbosity:

1. Save this in ~/test.sh:

   #!/bin/bash
   (cd ~/beam_integrals && IGNORE_SYMPY_MEMORY_LEAK=yes MPMATH_NOGMPY=yes python setup.py -q nosetests -x --tests tests.characteristic_equation_solvers.test_solvers)
   

2. And make the ~/test.sh executable:

   chmod +x ~/test.sh
   

3. To further decrease the reporting verbosity and also disable coverage plugins please remove the ~/beam_integrals/setup.cfg file:

   rm ~/beam_integrals/setup.cfg
   

Verify that SymPy 0.7.1 is "good" and SymPy 0.7.6 is "bad":

1. Switch over to SymPy 0.7.6 to make sure tests fail:

   cd ~/sympy
   git checkout sympy-0.7.6
   ~/test.sh
   

2. Switch over to SymPy 0.7.1 to make sure tests still pass:

   cd ~/sympy
   git checkout sympy-0.7.1
   ~/test.sh
   

Now let's use git bisect to figure out which commit introduced a regression between SymPy 0.7.1 and 0.7.6:

cd ~/sympy
git bisect start sympy-0.7.6 sympy-0.7.1
git bisect run ~/test.sh
git bisect reset

The results

It seems b84674d is the commit that introduced this regression.
Full git bisect output can be found at https://gist.github.com/petarmaric/5985ad8674247e0de78a

I've also ran ~/test.sh against the current master (50b4f3d ATM)

cd ~/sympy
git checkout master
pip install mpmath
~/test.sh

and, unfortunately, the tests still fail:

root@sympy-beam-ingrals-test:~/sympy# ~/test.sh
F
======================================================================
FAIL: tests.characteristic_equation_solvers.test_solvers.test_find_best_root(1, 1)
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/root/beam_integrals/nose-1.3.4-py2.7.egg/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/root/beam_integrals/tests/characteristic_equation_solvers/test_solvers.py", line 31, in check_find_best_root
    assert_less_equal(mu_m_error, t.MAX_ERROR_TOLERANCE)
AssertionError: 0.00000000000000032162452993532729844684607400882247424678888277467224259784769572226917549479169303989821176790032008305968824678769927238903329590530768150168932 not less than or equal to 1e-16

----------------------------------------------------------------------
Ran 1 test in 0.392s

FAILED (failures=1)

If it's worth anything, the current master generates the AssertionError with the number:

0.00000000000000032162452993532729844684607400882247...932

while all of the git bisect "bad" commits (https://gist.github.com/petarmaric/5985ad8674247e0de78a) return:

0.00000000000000012246467991473531772260659322749979...412

I've tried @cbm755 suggestion from #8821 (comment) to try his merge request #8826, but it resulted in the exact same failure as for the current master (50b4f3d ATM):

root@sympy-beam-ingrals-test:~/sympy# ~/test.sh
F
======================================================================
FAIL: tests.characteristic_equation_solvers.test_solvers.test_find_best_root(1, 1)
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/root/beam_integrals/nose-1.3.4-py2.7.egg/nose/case.py", line 197, in runTest
    self.test(*self.arg)
  File "/root/beam_integrals/tests/characteristic_equation_solvers/test_solvers.py", line 31, in check_find_best_root
    assert_less_equal(mu_m_error, t.MAX_ERROR_TOLERANCE)
AssertionError: 0.00000000000000032162452993532729844684607400882247424678888277467224259784769572226917549479169303989821176790032008305968824678769927238903329590530768150168932 not less than or equal to 1e-16

----------------------------------------------------------------------
Ran 1 test in 0.375s

FAILED (failures=1)

[petarmaric]

Can't prove it right now, but I have a "feeling" that this line b84674d#diff-6cd5274897d73197c25bc41e9afe63e4R565 (L565-L566) is causing issues with https://bitbucket.org/petar/beam_integrals/src/9a74e91d5a5f50e2639a7e3aeb6de9a486d6e59d/beam_integrals/characteristic_equation_solvers.py?at=default#cl-46 (see the comment on L50-L51)

[asmeurer]

@smichr that was your commit.

[cbm755]

@petamaric: if you're right, my change in #8826 won't help because the new path still runs that same code (as you've already discovered).

[cbm755]

Can't prove it right now, but I have a "feeling" that ...

@petarmaric one way to "prove" it would be a minimum (non)-working example---several lines of code that demonstrate what does wrong. I expect you're in a better position that we are to do this, as you (obviously) understand the beam_integrals code better.

This would go a long way to getting this resolved!

[petarmaric]

@cbm755 ask and ye shall receive

I've created a gist with a minimal working example showcasing which SymPy versions/commits are causing beam_integrals tests to fail. As an added bonus, the code has no external dependencies besides SymPy itself.

The testing environment was based on instructions in #2439 (comment)

The results are very interesting:

1. SymPy v0.7.1 - PASS. All mpmath_root and sympy_root pairs are the same.
2. SymPy v0.7.6 - FAIL. All mpmath_root's are equal to their values from SymPy v0.7.1. However sympy_root's get zero-clipped at about the 39th decimal, independent of decimal_precision used.
3. SymPy b84674d - FAIL. Same problem as with SymPy v0.7.6
4. SymPy 829bc1b (parent commit of b84674d) - PASS. Same results as with SymPy v0.7.1

IMO, this pretty much proves that b84674d introduced a regression.

[petarmaric]

If we remove the context, we can make the minimal working example even shorter:

>>> from sympy import mpmath, Float

>>> text = '3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068'
>>> decimal_precision = 100

>>> mpmath.mp.dps = decimal_precision

>>> mpmath_num = mpmath.mpf(text)
>>> sympy_num = Float(mpmath_num, decimal_precision)

>>> print "mpmath_num = %s\nsympy_num  = %s" % (mpmath_num, sympy_num)
mpmath_num = 3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068
sympy_num  = 3.141592653589793115997963468544185161590576171875000000000000000000000000000000000000000000000000000

This test was done with SymPy 0.7.5 on Python 2.7.8 win-amd64 |Anaconda 2.1.0 (64-bit)|

[petarmaric]

In fact, after closer inspection of the results discussed in my previous 2 comments I noticed something even more interesting:

>>> from sympy import mpmath, Float

>>> text = '3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068'
>>> decimal_precision = 100

>>> mpmath.mp.dps = decimal_precision

>>> mpmath_num = mpmath.mpf(text)
>>> sympy_num = Float(mpmath_num, decimal_precision)

>>> print "mpmath_num = %s\nsympy_num  = %s" % (mpmath_num, sympy_num)
mpmath_num = 3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068
sympy_num  = 3.141592653589793115997963468544185161590576171875000000000000000000000000000000000000000000000000000

>>> diff = mpmath_num - sympy_num
>>> print diff
0.0000000000000001224646799147353177226065932275001058209749445923078164062862089986280348253421170680078538784935795

>>> print mpmath.log10(diff)
-15.91198914827844864114174208533867293255232912230592779427060130266807642186107421296808278366048327

It seems not only are the SymPy values getting zero-clipped at about the 39th decimal, but there are also digit-level inaccuracies long before that.

Judging from the diff I'd say they are introduced at about the 16th decimal, and log10 of the diff is about 15.9 which looks suspiciously familiar... Something in sympy.Float seems to be lowering it's numeric precision to the level of float64, independent of high of a decimal precision is used.

[cbm755]

Thanks @petarmaric "context-free" is exactly what I meant. This does look like a bug in creating a Float from an mpf object: probably intermediate result converted to double precision. I'll take a look into this.

Note that sympify(mpmath_num) should work if that helps you work around.

[cbm755]

Have to do something else now, here's what I've found so far...

IMO, this pretty much proves that b84674d introduced a regression.

Could be. As I think you mentioned earlier, ln 565 causes mpf input to eventually get to this line:

                      _mpf_ = mpmath.mpf(
                        S.NegativeOne**num[0]*num[1]*2**num[2])._mpf_

This creates a 15 digit Float (that is it calls Float again makes a Python float while doing the num[1]*2**num[2] bit).

Deleting this bit:

elif isinstance(num, mpmath.mpf):
    num = num._mpf_

sort of helps with your problem:

mpmath.mp.dps = 100
a = mpmath.mpf(str(pi.evalf(128)))
sin(a)    % good, 2.8689e-102
b = Float(a, 100)
sin(b)    % good, 2.86889280175964e-102

Note however:

c = Float(a)
sin(c)   % also 2.86889280175964e-102
srepr(c)   % but above is much more accurate than expected from:
"Float('3.1415926535897932', prec=15)"

I tihnk this is because c just copied the mpf object... and technically it does have (at least) 15 digits of precision so it follows spec. @smichr may have been trying to avoid this "problem", based on his tests.

Perhaps the right fix is S.NegativeOne**num[0]*num[1]*2**_sympif\y(num[2]))._mpf_, that is to _sympify() the exponent to prevent Python float from coming from 2**<int>. It bothers me that this line makes a Float internally.

More later (couple days at least), but if @smichr has ideas go for it.

[updated comment to clarify a bit]

[petarmaric]

@cbm755 this is turning into a quite interesting discussion, thanks for the effort so far.

It bothers me that this line makes a Float internally.

Sorry, but I'm not sure if you meant Python float or SymPy Float here?

I tihnk this is because c just copied the mpf object... and technically it does have (at least) 15 digits of precision so it follows spec. @smichr may have been trying to avoid this "problem", based on his tests.

I solve most of these issues in beam_integrals by forcing a decimal_precision kwarg onto all functions performing numerical computations, to control the decimal precision explicitly between any API calls.
Is there a reason that this mpmath call _mpf_ = mpmath.mpf(S.NegativeOne**num[0]*num[1]*2**num[2])._mpf_ is:

1. done with the mixed SymPy/Python/mpmath eval and not pure mpmath code, i.e. using mpmath.power() instead of 2**num[2]?
2. doing computations in the mpmath default context precision, and not using the prec arg already provided within its owner Float.__new__(cls, num, prec=15) to create a local mpmath context via mpmath.workdps() before any mpmath related computations commence? I ask because mpmath.mpf() gives different results depending on the mpmath context:

>>> from sympy import mpmath

>>> mpmath.mpf('1/3')
mpf('0.33333333333333331')

# Changing the mpmath context gives different results
>>> mpmath.mp.dps = 50
>>> mpmath.mpf('1/3')
mpf('0.33333333333333333333333333333333333333333333333333311')

# Forcing a Python float here is not a very good idea
>>> mpmath.mpf(1/3.0)
mpf('0.33333333333333331482961625624739099293947219848632813')

[cbm755]

Well it does both: 2**num[2] makes a Python float and then S.NegativeOne stuff makes a SymPy Float. Then it converts the SymPy Float to an mpf, then finally to a SymPy Float. At very least, this is sub-optimal.

Your point 1. sounds reasonable. Maybe try to patch it and see if it fixes anything? Sorry I don't understand your point 2., probably because I've not dug deep enough into this code yet.

---

We have the minimal nonworking example, so just need someone to dig it and fix it now.

This is a long thread. For the benefit of others who might get here, here the minimal nonworking example (or #2439 (comment))

mpmath.mp.dps = 100
a = mpmath.mpf(str(pi.evalf(128)))  # a is 100 digit pi in mpf format
sin(a)    # good, 2.8689e-102
b = Float(a, 100)   # bad, b has only 15 digits precision
sin(b)
0.000000000000000122464....

I think you should change the title to "loss in precision when converting mpf to Float".

[cbm755]

@petarmaric a further thought: if you'd like to help but don't have time to really dive into the code, a PR for an "@xfail" test of the minimal nonworking example would probably be well-received!

[petarmaric]

RE point 2: currently all mpmath related code within Float.__new__() is executed in the mpmath default context. I might be in a situation where I want to keep the mpmath default context precision low (say 25 decimals) yet I wish to create a new sympy.Float at high precision (say 300 decimals).
As the code is right now any mpmath calls within Float.__new__(cls, num, prec) would be executed at 25 decimals instead of the desired 300 decimals. To amend this we should maybe build a new, local, mpmath context (i.e. via mpmath.workdps()) at the start of Float.__new__() and set it's decimal precision to prec. Anything else is bound to create accuracy issues, like the one shown in the code sample above.

RE changing the title: I may be in favor of creating a new issue instead whilst keeping this one open, as my original tests on top of this issue have shown there may be even more problems with SymPy > 0.7.1 (i.e. https://gist.github.com/petarmaric/6346863).

RE @xfail test: Sure, any documentation on how this is done?

[cbm755]

Sure new issue fine with me: it is nicely self contained. CC me.

Re: @xfail documentation, I'm not sure. But for example, look in sympy/core/tests/test_evalf.py, sympy/core/tests/test_numbers.py. Basically write the test that should pass and then put @XFAIL before it:

@XFAIL
def test_float_from_mpf_issue_xxxx()
    ....

where xxxx is the new issue number.

[cbm755]

... where I want to keep the mpmath default context precision low (say 25 decimals) yet I wish to create a new sympy.Float at high precision (say 300 decimals).

Minimal nonworking example please?

[cbm755]

RE @xfail test: Sure, any documentation on how this is done?

Oops, maybe I only answered the @XFAIL part of this, more generally:
https://github.com/sympy/sympy/wiki/Development-workflow

[jcrist]

ping @cbm755 , @petarmaric : what is the status of this?

[cbm755]

Maybe there are more issues here but I summarized my view here: #2439 (comment)

@petarmaric do you have time to work on this? If not, let us know. Re-reading the above, my suggestions were:
1. File new issue for #2439 (comment).
2. Make PR for XFAIL test.

[petarmaric]

@jcrist, @cbm755: I'd love to work on this, but I'm afraid I won't have any spare time till October, as I'm split between spending time with my newborn and working the PhD thesis.

I think it may be for the best to open a new issue, and once it is fixed I'll use beam_integrals tests to check again if there are more regressions in SymPy > 0.7.1

[oscarbenjamin]

The minimum working example described above now works fine:

In [13]: import mpmath                                                                                                            

In [14]: mpmath.mp.dps = 100 
    ...: a = mpmath.mpf(str(pi.evalf(128)))  # a is 100 digit pi in mpf format 
    ...: sin(a)    # good, 2.8689e-102 
    ...: b = Float(a, 100)   # bad, b has only 15 digits precision 
    ...: sin(b)                                                                                                                   
Out[14]: 2.8689e-102

Can anyone summarise what exactly is still to fix here?

[smichr]

I think you should change the title to "loss in precision when converting mpf to Float".
Can anyone summarise what exactly is still to fix here?

As you have confirmed, the loss of precision when converting mpf to Float has been resolved.