fix: treat floats like rationals in as_coeff_mul

oscarbenjamin commented 3 weeks ago

References to other Issues or PRs

Fixes gh-26571

Brief description of what is fixed or changed

Treat floats like rationals in as_coeff_mul by default. This mirrors the behaviour of as_coeff_Mul.

Other comments

Release Notes

core
- The as_coeff_mul method now treats floats like rationals by default. This fixes some integration bugs relating to integrals with floats.

sympy-bot commented 3 weeks ago

:white_check_mark:

Hi, I am the SymPy bot. I'm here to help you write a release notes entry. Please read the guide on how to write release notes.

Your release notes are in good order.

Here is what the release notes will look like:

core
- The as_coeff_mul method now treats floats like rationals by default. This fixes some integration bugs relating to integrals with floats. (#26573 by @oscarbenjamin)

This will be added to https://github.com/sympy/sympy/wiki/Release-Notes-for-1.13.

Click here to see the pull request description that was parsed.

#### References to other Issues or PRs Fixes gh-26571 #### Brief description of what is fixed or changed Treat floats like rationals in as_coeff_mul by default. This mirrors the behaviour of as_coeff_Mul. #### Other comments #### Release Notes * core * The `as_coeff_mul` method now treats floats like rationals by default. This fixes some integration bugs relating to integrals with floats.

github-actions[bot] commented 3 weeks ago

Benchmark results from GitHub Actions

Lower numbers are good, higher numbers are bad. A ratio less than 1 means a speed up and greater than 1 means a slowdown. Green lines beginning with + are slowdowns (the PR is slower then master or master is slower than the previous release). Red lines beginning with - are speedups.

Significantly changed benchmark results (PR vs master)

Significantly changed benchmark results (master vs previous release)

| Change   | Before [2487dbb5]    | After [26584616]    |   Ratio | Benchmark (Parameter)                                                |
|----------|----------------------|---------------------|---------|----------------------------------------------------------------------|
| -        | 70.1±0.9ms           | 44.2±0.3ms          |    0.63 | integrate.TimeIntegrationRisch02.time_doit(10)                       |
| -        | 68.3±1ms             | 42.9±0.2ms          |    0.63 | integrate.TimeIntegrationRisch02.time_doit_risch(10)                 |
| +        | 18.2±0.1μs           | 30.4±0.2μs          |    1.67 | integrate.TimeIntegrationRisch03.time_doit(1)                        |
| -        | 5.36±0.06ms          | 2.89±0.02ms         |    0.54 | logic.LogicSuite.time_load_file                                      |
| -        | 73.9±0.7ms           | 28.4±0.2ms          |    0.38 | polys.TimeGCD_GaussInt.time_op(1, 'dense')                           |
| -        | 25.6±0.1ms           | 16.8±0.03ms         |    0.66 | polys.TimeGCD_GaussInt.time_op(1, 'expr')                            |
| -        | 74.4±0.4ms           | 28.7±0.07ms         |    0.39 | polys.TimeGCD_GaussInt.time_op(1, 'sparse')                          |
| -        | 256±2ms              | 123±0.8ms           |    0.48 | polys.TimeGCD_GaussInt.time_op(2, 'dense')                           |
| -        | 256±1ms              | 124±0.8ms           |    0.49 | polys.TimeGCD_GaussInt.time_op(2, 'sparse')                          |
| -        | 652±7ms              | 373±0.9ms           |    0.57 | polys.TimeGCD_GaussInt.time_op(3, 'dense')                           |
| -        | 660±5ms              | 372±0.9ms           |    0.56 | polys.TimeGCD_GaussInt.time_op(3, 'sparse')                          |
| -        | 492±1μs              | 287±1μs             |    0.58 | polys.TimeGCD_LinearDenseQuadraticGCD.time_op(1, 'dense')            |
| -        | 1.80±0.04ms          | 1.05±0ms            |    0.58 | polys.TimeGCD_LinearDenseQuadraticGCD.time_op(2, 'dense')            |
| -        | 5.78±0.09ms          | 3.09±0.01ms         |    0.53 | polys.TimeGCD_LinearDenseQuadraticGCD.time_op(3, 'dense')            |
| -        | 452±2μs              | 231±2μs             |    0.51 | polys.TimeGCD_QuadraticNonMonicGCD.time_op(1, 'dense')               |
| -        | 1.47±0.01ms          | 670±3μs             |    0.46 | polys.TimeGCD_QuadraticNonMonicGCD.time_op(2, 'dense')               |
| -        | 4.88±0.03ms          | 1.63±0.02ms         |    0.33 | polys.TimeGCD_QuadraticNonMonicGCD.time_op(3, 'dense')               |
| -        | 376±3μs              | 205±1μs             |    0.55 | polys.TimeGCD_SparseGCDHighDegree.time_op(1, 'dense')                |
| -        | 2.44±0ms             | 1.25±0.01ms         |    0.51 | polys.TimeGCD_SparseGCDHighDegree.time_op(3, 'dense')                |
| -        | 10.2±0.1ms           | 4.47±0.02ms         |    0.44 | polys.TimeGCD_SparseGCDHighDegree.time_op(5, 'dense')                |
| -        | 363±0.4μs            | 173±1μs             |    0.48 | polys.TimeGCD_SparseNonMonicQuadratic.time_op(1, 'dense')            |
| -        | 2.52±0.02ms          | 912±8μs             |    0.36 | polys.TimeGCD_SparseNonMonicQuadratic.time_op(3, 'dense')            |
| -        | 9.66±0.02ms          | 2.66±0.02ms         |    0.28 | polys.TimeGCD_SparseNonMonicQuadratic.time_op(5, 'dense')            |
| -        | 1.04±0.01ms          | 426±3μs             |    0.41 | polys.TimePREM_LinearDenseQuadraticGCD.time_op(3, 'dense')           |
| -        | 1.73±0ms             | 510±3μs             |    0.29 | polys.TimePREM_LinearDenseQuadraticGCD.time_op(3, 'sparse')          |
| -        | 5.85±0.04ms          | 1.81±0.02ms         |    0.31 | polys.TimePREM_LinearDenseQuadraticGCD.time_op(5, 'dense')           |
| -        | 8.35±0.04ms          | 1.51±0.01ms         |    0.18 | polys.TimePREM_LinearDenseQuadraticGCD.time_op(5, 'sparse')          |
| -        | 289±1μs              | 65.9±0.3μs          |    0.23 | polys.TimePREM_QuadraticNonMonicGCD.time_op(1, 'sparse')             |
| -        | 3.48±0.04ms          | 391±3μs             |    0.11 | polys.TimePREM_QuadraticNonMonicGCD.time_op(3, 'dense')              |
| -        | 4.09±0.02ms          | 279±0.7μs           |    0.07 | polys.TimePREM_QuadraticNonMonicGCD.time_op(3, 'sparse')             |
| -        | 7.18±0.03ms          | 1.28±0.01ms         |    0.18 | polys.TimePREM_QuadraticNonMonicGCD.time_op(5, 'dense')              |
| -        | 8.87±0.04ms          | 844±7μs             |    0.1  | polys.TimePREM_QuadraticNonMonicGCD.time_op(5, 'sparse')             |
| -        | 5.24±0.1ms           | 2.98±0.02ms         |    0.57 | polys.TimeSUBRESULTANTS_LinearDenseQuadraticGCD.time_op(2, 'sparse') |
| -        | 12.0±0.1ms           | 6.59±0.06ms         |    0.55 | polys.TimeSUBRESULTANTS_LinearDenseQuadraticGCD.time_op(3, 'dense')  |
| -        | 22.4±0.1ms           | 9.09±0.02ms         |    0.41 | polys.TimeSUBRESULTANTS_LinearDenseQuadraticGCD.time_op(3, 'sparse') |
| -        | 5.24±0.02ms          | 872±7μs             |    0.17 | polys.TimeSUBRESULTANTS_QuadraticNonMonicGCD.time_op(1, 'sparse')    |
| -        | 12.6±0.04ms          | 7.07±0.04ms         |    0.56 | polys.TimeSUBRESULTANTS_QuadraticNonMonicGCD.time_op(2, 'sparse')    |
| -        | 103±1ms              | 25.7±0.2ms          |    0.25 | polys.TimeSUBRESULTANTS_QuadraticNonMonicGCD.time_op(3, 'dense')     |
| -        | 169±2ms              | 54.0±0.3ms          |    0.32 | polys.TimeSUBRESULTANTS_QuadraticNonMonicGCD.time_op(3, 'sparse')    |
| -        | 176±0.8μs            | 111±0.8μs           |    0.63 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(1, 'dense')      |
| -        | 364±1μs              | 216±3μs             |    0.59 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(1, 'sparse')     |
| -        | 4.28±0.02ms          | 843±5μs             |    0.2  | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(3, 'dense')      |
| -        | 5.36±0.1ms           | 381±2μs             |    0.07 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(3, 'sparse')     |
| -        | 20.5±0.05ms          | 2.78±0.01ms         |    0.14 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(5, 'dense')      |
| -        | 23.2±0.1ms           | 627±1μs             |    0.03 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(5, 'sparse')     |
| -        | 482±4μs              | 139±2μs             |    0.29 | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(1, 'sparse') |
| -        | 4.90±0.08ms          | 621±2μs             |    0.13 | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(3, 'dense')  |
| -        | 5.31±0.03ms          | 138±1μs             |    0.03 | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(3, 'sparse') |
| -        | 13.4±0.1ms           | 1.29±0ms            |    0.1  | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(5, 'dense')  |
| -        | 14.1±0.08ms          | 142±1μs             |    0.01 | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(5, 'sparse') |
| -        | 133±1μs              | 77.2±0.3μs          |    0.58 | solve.TimeMatrixOperations.time_rref(3, 0)                           |
| -        | 250±0.8μs            | 89.7±0.8μs          |    0.36 | solve.TimeMatrixOperations.time_rref(4, 0)                           |
| -        | 24.1±0.08ms          | 10.2±0.03ms         |    0.42 | solve.TimeSolveLinSys189x49.time_solve_lin_sys                       |
| -        | 28.0±0.1ms           | 15.5±0.09ms         |    0.55 | solve.TimeSparseSystem.time_linsolve_Aaug(20)                        |
| -        | 54.5±0.2ms           | 24.8±0.2ms          |    0.46 | solve.TimeSparseSystem.time_linsolve_Aaug(30)                        |
| -        | 27.9±0.2ms           | 15.0±0.08ms         |    0.54 | solve.TimeSparseSystem.time_linsolve_Ab(20)                          |
| -        | 54.1±0.3ms           | 24.4±0.2ms          |    0.45 | solve.TimeSparseSystem.time_linsolve_Ab(30)                          |

Full benchmark results can be found as artifacts in GitHub Actions (click on checks at the top of the PR).

smichr commented 3 weeks ago

I was distracted by the other changes the first time I looked at this. Now I am seeing that as_coeff_add doesn't have a flag to allow Rational discrimination while as_coeff_mul does. And all that you are doing is changing the default value of the flag. I would prefer to leave the flag unchanged since the notion of Rational coefficient (or slot 0 is for the Rational) is pretty well baked into SymPy (I think). All that needs to be done without breaking the behavior is to use the flag rational=True within the integration routines that are failing.

If you think the default should be changed, then as_coeff_add should be changed, too, to see how significant a change that will be in terms of codebase expectations. If it doesn't make a difference there, then perhaps my misgivings are unfounded. But if it does make a difference then we will have painted ourselves into an awkward corner, having as_coeff_mul behaving differently than as_coeff_add.

oscarbenjamin commented 3 weeks ago

But if it does make a difference then we will have painted ourselves into an awkward corner, having as_coeff_mul behaving differently than as_coeff_add.

We have already painted ourselves into an awkward corner by handling floats differently from rationals for no reason in such widely used basic routines.

There are endless float bugs in situations where a float should really be no different from a rational. A primary source is the fact that low-level functions like as_coeff_mul that are used throughout the codebase handle rationals differently for no reason. That is why I have changed the default value for the flag rather than overriding the default at particular call sites. Anywhere that uses as_coeff_mul should probably treat floats and rationals the same way.

oscarbenjamin commented 3 weeks ago

I just realised this needs to change the default in Mul.as_coeff_mul as well.

smichr commented 3 weeks ago

We have already painted ourselves into an awkward corner

The onus is on the user to select the correct flags for the function being used. We are not yet in a corner because the user can ask for the Rational coefficient if desired. Anyone writing a low level routine should be informed about the methods being used.

I suspect it will just be a matter of time before changing the default starts causing new bugs to arise because the routine thought it was handling a Rational coefficient when it was actually a Float. A case in point (as I review how as_coeff_mul is used) is in trigsimp.py where the following code appears:

trigterms = [(g.args[0].as_coeff_mul(), g.func) for g in gens
                     if g.func in allfuncs]

That is supposed to break up something like sin(3*x) into ((3,(x,)), sin) and igcd will be run on the coefficients. If those coefficients are Floats that are not equal to their int counterparts (e.g. sin(3.1 x)), that is going to raise an error. You did not run into this in tests, probably because test weren't written for cases that didn't apply -- i.e. thought being something like,

if as_coeff_mul()[0] only returns Rational then we don't have to test for Floats; and we already tested as_coeff_mul in test_expr to verify that that is the case

The seemingly innocuous change that this PR makes to test_expr is, in fact, a major change to the expectation of that low level routine. This has heretofore been an invariant: coeff will be Rational unless explicitly indicated otherwise. I do not see a good reason to change this when the routines that are failing can simply allow non-Rational to be returned instead.

sympy / sympy