Pre and post treatment of expression when cse=True in lambdify #26463

mleila1312 · 2024-04-05T10:10:25Z

Fixes #26404

These additions to the lambdify file fixes the problems seen in the issue mentionned above.

What is changed

When lambdify is called and the optionnal parameter cse==True, then we check if there are Derivatives in the arguments of the function. If so, they are replaced by variable names to prevent them from being ignored later in the cse process (see issue page here to see the problem). After this and the cse process is applied, we do a post process : we change back the Derivatives replaced at the beginning to their original values in the expression and in the cse replacements (ie the Derivatives).

To do this, I have added two functions in lambdify.py : pre_treatment_cse(args, expr) (line 181 in lambdify.py) and post_treatment_cse(dictionnary, args, expr, cses) (line 247 in lambdify.py)

I have also added new tests in the test_lambdify.py file. (end of the file)

If you have any questions, would like more comments in the code.... , please let me know!

Release Notes

functions
- Fixed a bug with cse in lambdify when there are Derivatives in the arguments and the expression.

sympy-bot · 2024-04-05T10:10:32Z

✅

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:

functions
- Fixed a bug with cse in lambdify when there are Derivatives in the arguments and the expression. (#26463 by @mleila1312)

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.

<!-- Your title above should be a short description of what
was changed. Do not include the issue number in the title. -->

### Fixes #26404
<!-- If this pull request fixes an issue, write "Fixes #NNNN" in that exact
format, e.g. "Fixes #1234" (see
https://tinyurl.com/auto-closing for more information). Also, please
write a comment on that issue linking back to this pull request once it is
open. -->
These additions to the lambdify file fixes the problems seen in the issue mentionned above.

### What is changed
When lambdify is called and the optionnal parameter cse==True, then we check if there are Derivatives in the arguments of the function. If so, they are replaced by variable names to prevent them from being ignored later in the cse process (see issue page [here](https://github.com/sympy/sympy/issues/26404) to see the problem). After this and the cse process is applied, we do a post process : we change back the Derivatives replaced at the beginning to their original values in the expression and in the cse replacements (ie the Derivatives).  

To do this, I have added two functions in lambdify.py : pre_treatment_cse(args, expr) (line 181 in lambdify.py) and post_treatment_cse(dictionnary, args, expr, cses)   (line 247 in lambdify.py) 

I have also added new tests in the test_lambdify.py file. (end of the file)

If you have any questions, would like more comments in the code.... , please let me know!


#### Release Notes

<!-- Write the release notes for this release below between the BEGIN and END
statements. The basic format is a bulleted list with the name of the subpackage
and the release note for this PR. For example:

* solvers
  * Added a new solver for logarithmic equations.

* functions
  * Fixed a bug with log of integers. Formerly, `log(-x)` incorrectly gave `-log(x)`.

* physics.units
  * Corrected a semantical error in the conversion between volt and statvolt which
    reported the volt as being larger than the statvolt.

or if no release note(s) should be included use:

NO ENTRY

See https://github.com/sympy/sympy/wiki/Writing-Release-Notes for more
information on how to write release notes. The bot will check your release
notes automatically to see if they are formatted correctly. -->

<!-- BEGIN RELEASE NOTES -->
* functions
  * Fixed a bug with cse in lambdify when there are Derivatives in the arguments and the expression.

<!-- END RELEASE NOTES -->

moorepants · 2024-04-05T12:39:47Z

Hi, thanks for preparing this.

I have a few questions:

Is there a reason not to use the optimizations kwarg in cse for this (as I suggested in the issue)? The machinery is already there to to pre and post work on expression that are processed.
The pre and post functions seem elaborate, yet there is only one simple test. Do you need to test other possibilities given the complexity (I assume it is there there to handle other cases).
Also given the complexity mentioned in point 2., how much will this slow things down? You can use the n_link_pendulum model in the mechanics models module to generate some really long expressions with derivatives to test, if you need a test case.

mleila1312 · 2024-04-05T14:26:37Z

@moorepants

concerning the optimization kwarg, I didn't use them here because I have some problems passing optionnal arguments to the function cse (you can see it in the last message on the issue, but basically, I noticed that I couldn't pass a flag that said cse is used in lambdify, it still retained the default value in the function, same for the optionnal argument list). So, by measure of security, I didn't want to implement it by passing the pre and post operations in the optionnal argument optimization. Although, with the way I implemented this, it would not be possible since a dictionnary containing the Derivative- newname associations is returned by the pre cse operation to be used in the post cse operation.
concerning the test, I could add others. I can say that even if a Derivative (A) of a Derivative (B) and B are passed as arguments, it doesn't cause an issue, same if the Derivative is in a function but I could add them as a test also. (I will add it in my next commit, if you have other ideas for tests, don't hesitate)
For the time complexity, if I test :

import sympy as sm
import sympy.physics.mechanics.models as models
from sympy import (cos, sin, Matrix, symbols, zeros)                     
import sympy.physics.mechanics.tests as tests
from sympy.physics.mechanics import (dynamicsymbols)

def timing(func, a, b, c, d, e, f,g, h, i , j, k, l, m, n):
    import time
    s= time.time()
    print(func(a, b, c, d, e, f,g, h, i , j, k, l, m, n))
    print(time.time()- s)


l0, m0 = sm.symbols("l0 m0")
l1, m1 = sm.symbols("l1 m1")
m2 = sm.symbols("m2")
g = sm.symbols("g")
q0, q1, q2 = dynamicsymbols("q0 q1 q2")
u0, u1, u2 = dynamicsymbols("u0 u1 u2")
F, T1 = dynamicsymbols("F T1")
kane1 = models.n_link_pendulum_on_cart(2)
massmatrix1 = Matrix([[m0 + m1 + m2, -l0*m1*cos(q1) - l0*m2*cos(q1),
                           -l1*m2*cos(q2)],
                          [-l0*m1*cos(q1) - l0*m2*cos(q1), l0**2*m1 + l0**2*m2,
                           l0*l1*m2*(sin(q1)*sin(q2) + cos(q1)*cos(q2))],
                          [-l1*m2*cos(q2),
                           l0*l1*m2*(sin(q1)*sin(q2) + cos(q1)*cos(q2)),
                           l1**2*m2]])
timing(sm.lambdify((l0, m0 ,l1, m1, m2, g, q0, q1, q2 , u0, u1, u2, F, T1), massmatrix1 - \
            kane1.mass_matrix), 0, 0 ,1, 1, 2, 9, 0, 1, 2 , 0, 1, 2, 13, 1)
#and
timing(sm.lambdify((l0, m0 ,l1, m1, m2, g, q0, q1, q2 , u0, u1, u2, F, T1), massmatrix1 -\
            kane1.mass_matrix, cse=True), 0, 0 ,1, 1, 2, 9, 0, 1, 2 , 0, 1, 2, 13, 1)

Average on a 100 executions :
For the not cse : 9.9945068359375 μs
For the one with cse= True : 9.99927520751953e μs

When there are Derivatives, I'm having trouble testing it, I tried :

import sympy as sm
import sympy.physics.mechanics.models as models
from sympy import (cos, sin, Matrix, symbols, zeros)                     
import sympy.physics.mechanics.tests as tests
from sympy.physics.mechanics import (dynamicsymbols)
def timing(func, a, b, c, d, e, f,g, h, i , j, k, l, m, n):
    import time
    s= time.time()
    print(func(a, b, c, d, e, f,g, h, i , j, k, l, m, n))
    print(time.time()- s)

q0, q1, q2 = dynamicsymbols("q0 q1 q2")
u0, u1, u2 = dynamicsymbols("u0 u1 u2")
l0, m0 = sm.symbols("l0 m0")
l1, m1 = sm.Function('l1')(l0),sm.Function('m1')(m0)
m2 = l1.diff(l0)
Derivative(l1(l0), l0)
g = m2.diff(l0)
F, T1 = dynamicsymbols("F T1")
kane1 = models.n_link_pendulum_on_cart(2)
massmatrix1 = Matrix([[m0 + m1 + m2, -l0*m1*cos(q1) - l0*m2*cos(q1),
                           -l1*m2*cos(q2)],
                          [-l0*m1*cos(q1) - l0*m2*cos(q1), l0**2*m1 + l0**2*m2,
                           l0*l1*m2*(sin(q1)*sin(q2) + cos(q1)*cos(q2))],
                          [-l1*m2*cos(q2),
                           l0*l1*m2*(sin(q1)*sin(q2) + cos(q1)*cos(q2)),
                           l1**2*m2]])
timing(sm.lambdify((l0, m0 ,l1, m1, m2, g, q0, q1, q2 , u0, u1, u2, F, T1), massmatrix1 -\
                   kane1.mass_matrix), 0, 0 ,1, 1, 2, 9, 0, 1, 2 , 0, 1, 2, 13, 1)

I get the error :

Traceback (most recent call last):
  File "<pyshell#139>", line 1, in <module>
    timing(sm.lambdify((l0, m0 ,l1, m1, m2, g, q0, q1, q2 , u0, u1, u2, F, T1), massmatrix1 - kane1.mass_matrix), 0, 0 ,1, 1, 2, 9, 0, 1, 2 , 0, 1, 2, 13, 1)
  File "<pyshell#127>", line 4, in timing
    print(func(a, b, c, d, e, f,g, h, i , j, k, l, m, n))
  File "<lambdifygenerated-1>", line 2, in _lambdifygenerated
    return ImmutableDenseMatrix([[_Dummy_56 + _Dummy_63 - m1 - m2, -_Dummy_56*l0*cos(_Dummy_61) - \
_Dummy_63*l0*cos(_Dummy_61) + l0*m1*cos(_Dummy_61) + l0*m2*cos(_Dummy_61), -\
_Dummy_56*_Dummy_64*cos(_Dummy_60) + l1*m2*cos(_Dummy_60)], [-_Dummy_56*l0*cos(_Dummy_61) - \
_Dummy_63*l0*cos(_Dummy_61) + l0*m1*cos(_Dummy_61) + l0*m2*cos(_Dummy_61), _Dummy_56*l0**2 + _Dummy_63*l0**2 -\
 l0**2*m1 - l0**2*m2, _Dummy_56*_Dummy_64*l0*(sin(_Dummy_60)*sin(_Dummy_61) + cos(_Dummy_60)*cos(_Dummy_61)) -\
 l0*l1*m2*(sin(_Dummy_60)*sin(_Dummy_61) + cos(_Dummy_60)*cos(_Dummy_61))], [-\
_Dummy_56*_Dummy_64*cos(_Dummy_60) + l1*m2*cos(_Dummy_60), _Dummy_56*_Dummy_64*l0*\
(sin(_Dummy_60)*sin(_Dummy_61) + cos(_Dummy_60)*cos(_Dummy_61)) - l0*l1*m2*(sin(_Dummy_60)*sin(_Dummy_61) +\
 cos(_Dummy_60)*cos(_Dummy_61)), _Dummy_64**2*_Dummy_56 - l1**2*m2]])
  File "C:\Users\Utilisateur\AppData\Local\Programs\Python\Python312\Lib\site-packages\mpmath\matrices\matrices.py", line 302, in __init__
    self[i, j] = a
  File "C:\Users\Utilisateur\AppData\Local\Programs\Python\Python312\Lib\site-packages\mpmath\matrices\matrices.py", line 560, in __setitem__
    value = self.ctx.convert(value)
  File "C:\Users\Utilisateur\AppData\Local\Programs\Python\Python312\Lib\site-packages\mpmath\ctx_mp_python.py", line 671, in convert
    return ctx._convert_fallback(x, strings)
  File "C:\Users\Utilisateur\AppData\Local\Programs\Python\Python312\Lib\site-packages\mpmath\ctx_mp.py", line 634, in _convert_fallback
    raise TypeError("cannot create mpf from " + repr(x))
TypeError: cannot create mpf from -m1 - m2 + 3

It seems that the derivatives aren't dummyfied in the matrix, this isn't a result of my changes since my changes are only applied when cse is True. Maybe I'm not using this right, otherwise, this could be another issue.

github-actions · 2024-04-05T15:49:18Z

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 [e89ee937]    |   Ratio | Benchmark (Parameter)                                                |
|----------|----------------------|---------------------|---------|----------------------------------------------------------------------|
| -        | 71.9±3ms             | 44.1±0.4ms          |    0.61 | integrate.TimeIntegrationRisch02.time_doit(10)                       |
| -        | 67.3±1ms             | 43.0±0.3ms          |    0.64 | integrate.TimeIntegrationRisch02.time_doit_risch(10)                 |
| +        | 18.8±0.2μs           | 30.1±0.4μs          |    1.6  | integrate.TimeIntegrationRisch03.time_doit(1)                        |
| -        | 5.49±0.05ms          | 2.87±0.04ms         |    0.52 | logic.LogicSuite.time_load_file                                      |
| -        | 74.2±0.5ms           | 29.1±0.2ms          |    0.39 | polys.TimeGCD_GaussInt.time_op(1, 'dense')                           |
| -        | 26.0±0.1ms           | 17.1±0.06ms         |    0.66 | polys.TimeGCD_GaussInt.time_op(1, 'expr')                            |
| -        | 73.8±0.5ms           | 29.1±0.2ms          |    0.39 | polys.TimeGCD_GaussInt.time_op(1, 'sparse')                          |
| -        | 260±1ms              | 126±0.2ms           |    0.48 | polys.TimeGCD_GaussInt.time_op(2, 'dense')                           |
| -        | 258±1ms              | 126±0.3ms           |    0.49 | polys.TimeGCD_GaussInt.time_op(2, 'sparse')                          |
| -        | 665±4ms              | 376±2ms             |    0.57 | polys.TimeGCD_GaussInt.time_op(3, 'dense')                           |
| -        | 664±4ms              | 377±2ms             |    0.57 | polys.TimeGCD_GaussInt.time_op(3, 'sparse')                          |
| -        | 505±4μs              | 285±0.8μs           |    0.56 | polys.TimeGCD_LinearDenseQuadraticGCD.time_op(1, 'dense')            |
| -        | 1.77±0.02ms          | 1.07±0ms            |    0.6  | polys.TimeGCD_LinearDenseQuadraticGCD.time_op(2, 'dense')            |
| -        | 5.83±0.03ms          | 3.05±0.01ms         |    0.52 | polys.TimeGCD_LinearDenseQuadraticGCD.time_op(3, 'dense')            |
| -        | 457±2μs              | 230±2μs             |    0.5  | polys.TimeGCD_QuadraticNonMonicGCD.time_op(1, 'dense')               |
| -        | 1.48±0ms             | 670±8μs             |    0.45 | polys.TimeGCD_QuadraticNonMonicGCD.time_op(2, 'dense')               |
| -        | 4.98±0.04ms          | 1.67±0.01ms         |    0.34 | polys.TimeGCD_QuadraticNonMonicGCD.time_op(3, 'dense')               |
| -        | 382±1μs              | 208±1μs             |    0.54 | polys.TimeGCD_SparseGCDHighDegree.time_op(1, 'dense')                |
| -        | 2.44±0.01ms          | 1.25±0.01ms         |    0.51 | polys.TimeGCD_SparseGCDHighDegree.time_op(3, 'dense')                |
| -        | 10.4±0.09ms          | 4.44±0.03ms         |    0.43 | polys.TimeGCD_SparseGCDHighDegree.time_op(5, 'dense')                |
| -        | 366±4μs              | 170±0.3μs           |    0.46 | polys.TimeGCD_SparseNonMonicQuadratic.time_op(1, 'dense')            |
| -        | 2.53±0.01ms          | 923±3μs             |    0.36 | polys.TimeGCD_SparseNonMonicQuadratic.time_op(3, 'dense')            |
| -        | 9.75±0.03ms          | 2.67±0.02ms         |    0.27 | polys.TimeGCD_SparseNonMonicQuadratic.time_op(5, 'dense')            |
| -        | 1.04±0.01ms          | 428±5μs             |    0.41 | polys.TimePREM_LinearDenseQuadraticGCD.time_op(3, 'dense')           |
| -        | 1.74±0.01ms          | 512±6μs             |    0.29 | polys.TimePREM_LinearDenseQuadraticGCD.time_op(3, 'sparse')          |
| -        | 5.87±0.06ms          | 1.79±0.01ms         |    0.31 | polys.TimePREM_LinearDenseQuadraticGCD.time_op(5, 'dense')           |
| -        | 8.50±0.1ms           | 1.48±0.01ms         |    0.17 | polys.TimePREM_LinearDenseQuadraticGCD.time_op(5, 'sparse')          |
| -        | 286±2μs              | 64.7±0.6μs          |    0.23 | polys.TimePREM_QuadraticNonMonicGCD.time_op(1, 'sparse')             |
| -        | 3.45±0.05ms          | 396±3μs             |    0.11 | polys.TimePREM_QuadraticNonMonicGCD.time_op(3, 'dense')              |
| -        | 4.01±0.01ms          | 279±2μs             |    0.07 | polys.TimePREM_QuadraticNonMonicGCD.time_op(3, 'sparse')             |
| -        | 7.28±0.08ms          | 1.26±0.01ms         |    0.17 | polys.TimePREM_QuadraticNonMonicGCD.time_op(5, 'dense')              |
| -        | 8.95±0.07ms          | 830±3μs             |    0.09 | polys.TimePREM_QuadraticNonMonicGCD.time_op(5, 'sparse')             |
| -        | 5.01±0.04ms          | 2.95±0.01ms         |    0.59 | polys.TimeSUBRESULTANTS_LinearDenseQuadraticGCD.time_op(2, 'sparse') |
| -        | 12.2±0.07ms          | 6.60±0.03ms         |    0.54 | polys.TimeSUBRESULTANTS_LinearDenseQuadraticGCD.time_op(3, 'dense')  |
| -        | 22.6±0.2ms           | 8.94±0.02ms         |    0.4  | polys.TimeSUBRESULTANTS_LinearDenseQuadraticGCD.time_op(3, 'sparse') |
| -        | 5.26±0.02ms          | 862±3μs             |    0.16 | polys.TimeSUBRESULTANTS_QuadraticNonMonicGCD.time_op(1, 'sparse')    |
| -        | 12.8±0.1ms           | 6.97±0.01ms         |    0.54 | polys.TimeSUBRESULTANTS_QuadraticNonMonicGCD.time_op(2, 'sparse')    |
| -        | 104±1ms              | 26.1±0.1ms          |    0.25 | polys.TimeSUBRESULTANTS_QuadraticNonMonicGCD.time_op(3, 'dense')     |
| -        | 168±0.9ms            | 53.0±0.3ms          |    0.32 | polys.TimeSUBRESULTANTS_QuadraticNonMonicGCD.time_op(3, 'sparse')    |
| -        | 173±2μs              | 111±1μs             |    0.64 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(1, 'dense')      |
| -        | 363±6μs              | 216±1μs             |    0.59 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(1, 'sparse')     |
| -        | 4.32±0.02ms          | 833±4μs             |    0.19 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(3, 'dense')      |
| -        | 5.31±0.04ms          | 381±2μs             |    0.07 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(3, 'sparse')     |
| -        | 19.9±0.05ms          | 2.82±0.01ms         |    0.14 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(5, 'dense')      |
| -        | 23.2±0.3ms           | 617±3μs             |    0.03 | polys.TimeSUBRESULTANTS_SparseGCDHighDegree.time_op(5, 'sparse')     |
| -        | 482±4μs              | 135±2μs             |    0.28 | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(1, 'sparse') |
| -        | 4.84±0.02ms          | 619±3μs             |    0.13 | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(3, 'dense')  |
| -        | 5.48±0.05ms          | 137±1μs             |    0.03 | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(3, 'sparse') |
| -        | 13.7±0.1ms           | 1.29±0.01ms         |    0.09 | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(5, 'dense')  |
| -        | 14.3±0.09ms          | 145±10μs            |    0.01 | polys.TimeSUBRESULTANTS_SparseNonMonicQuadratic.time_op(5, 'sparse') |
| -        | 136±0.5μs            | 75.1±0.2μs          |    0.55 | solve.TimeMatrixOperations.time_rref(3, 0)                           |
| -        | 254±1μs              | 87.8±0.2μs          |    0.35 | solve.TimeMatrixOperations.time_rref(4, 0)                           |
| -        | 24.5±0.1ms           | 10.2±0.03ms         |    0.42 | solve.TimeSolveLinSys189x49.time_solve_lin_sys                       |
| -        | 28.2±0.3ms           | 15.3±0.06ms         |    0.54 | solve.TimeSparseSystem.time_linsolve_Aaug(20)                        |
| -        | 55.0±0.4ms           | 25.1±0.2ms          |    0.46 | solve.TimeSparseSystem.time_linsolve_Aaug(30)                        |
| -        | 28.0±0.2ms           | 15.2±0.08ms         |    0.54 | solve.TimeSparseSystem.time_linsolve_Ab(20)                          |
| -        | 54.0±0.3ms           | 24.7±0.05ms         |    0.46 | 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 · 2024-04-05T16:08:53Z

You might consider using _masked from solveset.py.

mleila1312 · 2024-04-05T17:41:32Z

You might consider using _masked from solveset.py.

I looked into it, but I'm not sure if it would be compatible for dummyfying different derivatives, from what I've seen, it can only replace atoms like x, y, ... And a derivative has atoms as arguments, so I don't think it's an atom, but please tell me if I'm wrong, it would be more practical

I've detected another problem, cses can also replace combinations of the replaced Derivatives in the pre-treatment, I've solved this issue, I'm in the process of building a test with some matrix, I will push the updates soon.

mleila1312 · 2024-04-05T18:29:33Z

About the time of execution, in this new case(with more arguments ans more instances of Derivatives to replace in the expression) :

l0, m0 = sm.symbols("l0 m0")
l1, m1 = sm.symbols("l1 m1")
m2 = sm.symbols("m2")
g = symbols("g")
q0, q1, q2 = sm.Function("q0")(l0),sm.Function("q1")(l1),sm.Function("q2")(m0)
u1, u2 =q1.diff(l1), q2.diff(m0)
F, T1 = dynamicsymbols("F T1")
massmatrix1 = Matrix([[m0 + m1 + m2, -l0*m1*cos(q1) - l0*m2*cos(q1),
                           -l1*m2*cos(q2)],
                          [-l0*m1*cos(q1) - l0*m2*cos(q1), l0**2*m1 + l0**2*m2,
                           l0*l1*m2*(sin(q1)*sin(q2) + cos(q1)*cos(q2))],
                          [-l1*m2*cos(q2),
                           l0*l1*m2*(sin(q1)*sin(q2) + cos(q1)*cos(q2)),
                           l1**2*m2]])

forcing1 = Matrix([[-l0*m1*u1**2*sin(q1) - l0*m2*u1**2*sin(q1) -
                        l1*m2*u2**2*sin(q2) + F,
                       g*l0*m1*sin(q1) + g*l0*m2*sin(q1) -
                        l0*l1*m2*(sin(q1)*cos(q2) - sin(q2)*cos(q1))*u2**2,
                       g*l1*m2*sin(q2) - l0*l1*m2*(-sin(q1)*cos(q2) +
                                                    sin(q2)*cos(q1))*u1**2],
                   [-l0*m1*u1**2*sin(q1) - l0*m2*u1**2*sin(q1) -
                        l1*m2*u2**2*sin(q2) + F,
                       g*l0*m1*sin(q1) + g*l0*m2*sin(q1) -
                        l0*l1*m2*(sin(q1)*cos(q2) - sin(q2)*cos(q1))*u2**2,
                       g*l1*m2*sin(q2) - l0*l1*m2*(-sin(q1)*cos(q2) +
                                                    sin(q2)*cos(q1))*u1**2],
                   [-l0*m1*u1**2*sin(q1) - l0*m2*u1**2*sin(q1) -
                        l1*m2*u2**2*sin(q2) + F,
                       g*l0*m1*sin(q1) + g*l0*m2*sin(q1) -
                        l0*l1*m2*(sin(q1)*cos(q2) - sin(q2)*cos(q1))*u2**2,
                       g*l1*m2*sin(q2) - l0*l1*m2*(-sin(q1)*cos(q2) +
                                                    sin(q2)*cos(q1))*u1**2]])
# time of this =t1
(lambdify((l0, m0 ,l1, m1, m2, g, q0, q1, q2 ,  u1, u2, F, T1), massmatrix1 -\
           forcing1)( 0, 0 ,1, 1, 2, 9, 0, 1, 2 , 1, 2, 13, 1))
#time of this = t2
(lambdify((l0, m0 ,l1, m1, m2, g, q0, q1, q2 ,  u1, u2, F, T1), massmatrix1 -forcing1, \
                  cse=True)( 0, 0 ,1, 1, 2, 9, 0, 1, 2 , 1, 2, 13, 1))

(replaced time of execution, the formers were false, it was longer because of print)
Average on 100 tries :
t1 = 9.9945068359375 μs
t2 = 9.95159149169922 μs

Another example (that wasn't added as a test in the test file) :

import sympy as sm
from sympy import (cos, sin, Matrix, symbols, zeros)

def moy_timing(func, a, b, c, d, e, f,g, h, i , j, k, l, m, n, o, p, q, r, s, t, u, v, w, x, y, z):
    import time
    sum_time = 0
    for i in range(100) :
        t = time.time()
        _ = func(a, b, c, d, e, f,g, h, i , j, k, l, m,n, o, p, q, r, s, t, u, v, w, x, y, z)
        sum_time += time.time()-t
    return sum_time/100

#Variables
a,b,c,d,e,f,g,h,i,j,k,l,m =sm.symbols("a b c d e f g h i j k l m")
n, o, p,q = sm.Function("n")(a),sm.Function("o")(b),sm.Function("p")(c),sm.Function("q")(d)
r,s,t,u =n.diff(a), o.diff(b), p.diff(c), q.diff(d)
v,w,x,y =r.diff(a), s.diff(b), t.diff(c), u.diff(d)
z = v.diff(a)

massmatrix1 = Matrix([[a*b*c*d*w*r - cos(p*q-u*x) +z*sin(e+f+g+h+i+x+w) -cos(z*g*5-l*m+p*q)+\
       j*k*l-m+sin(z*cos(n-o+p*s-t)+v*y)+\
       sin( e+sin( sin(t*r*e*r-t-e+o+r*d)*z+e+r*f-s*k)*d+s*g-d*e*y-c-a*e),
                       a*f*d*e*-e+r-t*c-w+z*sin(z*cos(z*z)),
                           d*c*s*9+q-e*r+sin(z*z*z)],
                          [d*c*s*9+q-e*r+sin(z*z*z), a*f*d*e*-e+r-t*c-w+z*sin(z*cos(z*z)),
                           a*b*c*d*w*r - cos(p*q-u*x) +z*sin(e+f+g+h+i+x+w) -cos(z*g*5-l*m+p*q)+\
       j*k*l-m+sin(z*cos(n-o+p*s-t)+v*y)+\
       sin( e+sin( sin(t*r*e*r-t-e+o+r*d)*z+e+r*f-s*k)*d+s*g-d*e*y-c-a*e)],
                          [d*c*s*9+q-e*r+sin(w*w*w),
                           a*b*c*d*w*r - cos(p*q-u*x) +z*sin(e+f+g+h+i+x+w) -cos(z*g*5-l*m+p*q)+\
       j*k*l-m+sin(z*cos(n-o+p*s-t)+v*y)+\
       sin( e+sin( sin(t*r*e*r-t-e+o+r*d)*z+e+r*f-s*k)*d+s*g-d*e*y-c-a*e),
                           a*f*d*e*-e+r-t*c-w+z*sin(z*cos(h*h))]])

massmatrix2 = Matrix([[sin(a*v*z-e-f+t)+cos(s+d+f-e-t+sin(z*z*z*z)),
                       a+d+e-f-g*c+v*b-s*cos(cos(cos(cos(z*w)))),
                           r*e+c*k-q+m*l+y+u*sin(sin(sin(u*u*u*f)))],
                          [r*e+c*k-q+m*l+y+u*sin(sin(sin(u*u*u*f))),
                           sin(a*v*z-e-f+t)+cos(s+d+f-e-t+sin(z*z*z*z)),
                           a+d+e-f-g*c+v*b-s*cos(cos(cos(cos(z*w))))],
                          [a+d+e-f-g*c+v*b-s*cos(cos(cos(cos(z*w)))),
                           r*e+c*k-q+m*l+y+u*sin(sin(sin(u*u*u*f))),
                           sin(a*v*z-e-f+t)+cos(s+d+f-e-t+sin(z*z*z*z))]])
expr = massmatrix1-massmatrix2

t1 = moy_timing(sm.lambdify((a, b, c, d, e, f,g, h, i , j, k, l, m,n, o, p, q, r, s, t, u, v, w, x, y, z), expr),\
 0,12,54,6,2,3,4,5,8,9,7,6,-5,-7,23,9,5,6,1,2,0,-5,6,4,9,8)

t2 = moy_timing(sm.lambdify((a, b, c, d, e, f,g, h, i , j, k, l, m,n, o, p, q, r, s, t, u, v, w, x, y, z), expr,\
 cse=True), 0,12,54,6,2,3,4,5,8,9,7,6,-5,-7,23,9,5,6,1,2,0,-5,6,4,9,8)

Here :

t1= 19.96755599975586 μs
t2=10.018348693847657 μs

If you want me to do more complicated tests or if you have any other question, please let me know

.mailmap

sympy/utilities/lambdify.py

sympy/utilities/tests/test_lambdify.py

sympy/utilities/lambdify.py

moorepants · 2024-04-06T08:02:04Z

To test the effects of slowdown you need something like kane1 = models.n_link_pendulum_on_cart(20) at least.

mleila1312 · 2024-04-06T10:30:22Z

To test the effects of slowdown you need something like kane1 = models.n_link_pendulum_on_cart(20) at least.

I've done several tests to see the duration of execution.
For kane = models.n_link_pendulum_on_cart(20) , average for a 100 executions :

without cse : 1.0021281242370604 ms
with cse : 0.5107760429382324 ms
Before implementation with cse enabled : 0.4900932312011719 ms

Here, activating cse reduces the time of execution by half with the implementation.
Before the implementation, execution time was 0.039% shorter.

For kane = models.n_link_pendulum_on_cart(40) , average for a 100 executions :

without cse : 5.916361808776856 ms
with cse : 2.1598124504089355 ms
Before the implementation with cse enabled : 2.75665283203125 ms

Here, activating cse reduces the time of execution by 63% with the implementation.
Before the implementation, execution time was 0.21% longer.

For kane = models.n_link_pendulum_on_cart(100) , average for a 100 executions :

without cse : 81.51222705841064 ms
with cse : 22.987971305847166 ms
Before the implementation with cse enabled : 23.05973768234253 ms

Here, activating cse reduces the time of execution by 72% with the implementation.
Before the implementation, execution time was 0.003% longer.

If you want to try it for youselves, I've put the tests I've used here (the execution time should vary depending on the user's computer and sorry Ifor the length of the tes's code , I could have made it shorter)

Do you want me to test this with a larger n for kane = models.n_link_pendulum_on_cart(n) to see if the gap widens?

moorepants · 2024-04-06T11:52:19Z

Sorry if I misstated. The slowdown I'm interested in is the timing of lambdify(). If we do pre and post substitutions, there is a slowdown associated with executing the lambdify function. It will only be seen on very large expressions that are full of Derivative() terms.

moorepants · 2024-04-06T11:53:13Z

Are you saying that your modification makes lambdify faster?

sympy/utilities/lambdify.py

moorepants · 2024-04-06T11:56:23Z

sympy/utilities/tests/test_lambdify.py

+                            x0*l1*m2*(sin(q1)*sin(q2) + cos(q1)*cos(q2))],
+                            [-l1*m2*cos(q2),
+                            x0*l1*m2*(sin(q1)*sin(q2) + cos(q1)*cos(q2)),
+                            l1**2*m2]])


These expressions would be more appropriate if the u values were all swapped with derivatives: u1 = Derivative(q1(t) t), for example.

u values are : u1, u2 =q1.diff(l1), q2.diff(m0), so already derivatives, did you mean putting more of them in the expression of the matrix? Or do you mean putting Derivatives directly in the matrix? If so, I haven't handled the case were the expression of a Derivative is directly in the expression and not also as an argument, I will add it

sympy/utilities/lambdify.py

mleila1312 · 2024-04-06T12:02:53Z

Sorry if I misstated. The slowdown I'm interested in is the timing of lambdify(). If we do pre and post substitutions, there is a slowdown associated with executing the lambdify function. It will only be seen on very large expressions that are full of Derivative() terms.

No problem, I wasn't sure sorry! I'll do that today(so max in 8 hours) and post the new tests results, but I won't be able to compare to before the implementation(since it didn't work, although I could compare to the execution time of expressions without Derivative terms)

Are you saying that your modification makes lambdify faster?

It seems counterintuitive, I still wanted to add the results since the goal was to compare, but the fact that the execution is faster with the implementation might be due to other processes that were running on my computer (I didn't execute anything else, but I don't see any other explanation since the kane matrix didn't have any Derivative so didn't need modifications by the pre treatment)

moorepants · 2024-04-06T12:05:33Z

Ideally a fix would work with any types of SymPy objects where this may occur. Derivative() is one, but this may also occur with things like Integral(f(t), t) or even things like x(t) + sin(t). Is your solution only for Derivative()?

moorepants · 2024-04-06T12:07:30Z

No problem, I wasn't sure sorry! I'll do that today(so max in 8 hours)

No apologies are needed because you have not done anything wrong and there are no time limits so please take whatever time you need/desire.

mleila1312 · 2024-04-06T13:52:48Z

Derivative() is one, but this may also occur with things like Integral(f(t), t) or even things like x(t) + sin(t). Is your solution only for Derivative()?

l0=symbols('l0')
x = Function('x')(sin(l0)+l0)
y=Function('y')(x.diff(l0)*x.diff(l0))
sm.lambdify((y, x, l0), y+x*y*l0 +x*x +l0*x*x, cse=True)(1,1,1)

And

In[1] : print(y+x*y*lo +x*x +lo*x*x)
                                         
Out[1] : lo*x(lo + sin(lo))**2 + 
             lo*x(lo + sin(lo))*y((cos(lo) + 1)**2*Subs(Derivative(x(_xi_1), _xi_1), _xi_1, lo + sin(lo))**2) +
             x(lo + sin(lo))**2 + y((cos(lo) + 1)**2*Subs(Derivative(x(_xi_1), _xi_1), _xi_1, lo + sin(lo))**2)

With this, I get the original error, so ideally, a solution for this case should be implemented as well. I checked and the pre and post operations did not change the expression. And

sm.lambdify((x, lo), x*lo +x*x +lo*x*x, cse=True)(1,1)

Worked correctly

For other classes of objects, I think my solution would work, we would just need to add to the isinstance line 238 in lambdify.py all the other classes we want to mask, but it should be tested.

Although, concerning the integral, I tried to pass it as an argument (without cse activated) and I got an error, but maybe I didn't write the command lines correctly.

No apologies are needed because you have not done anything wrong and there are no time limits so please take whatever time you need/desire.

Thank you

mleila1312 · 2024-04-06T17:30:14Z

I did the new tests, here are the results on average for a 100 executions :

for kane 20 :
+ with cse enabled : 10.004043579101563 μs
+ without cse : 10.006427764892578 μs

Execution with cse is almost the same as the one without.

for kane 40 :
+ with cse enabled : 60.02426147460937 μs
+ without cse : 60.048103332519535 μs

Execution with cse is 0.0003% shorter.

for kane 100 :
+ with cse enabled : 0.9108424186706543 ms
+ without cse : 0.8992886543273925 ms

Execution with cse is 0.01% longer.

for kane 200 : (added on 07/04 to see if the difference with and without cse increased)
+ with cse enabled : 4.0581893920898435 ms
+ without cse : 4.071612358093262 ms

Execution with cse is 0.003% shorter.

The percentage (updated with the new changes to accomodate lists) is really small so I think on average the execution time with and without cse are equivalent but if you want to execute yourselves the tests I've done and measure the average time spent per execution, they are here

I changed and measured the time of execution with different Matrices(or list of lists here) similar to the kane built-in(the built-in function didn't take Derivatives as arguments so I made one similar in dimensions)

Also, one of the test before the last commit failed because some tests sent tuples as expr to lambdify, so I updated the description on the function : it also takes tuple type of expressions as expr

Before, when there were Derivatives in expr ans args given to lambdify with cse enabled, there was an error because the cse treatment changed the arguments of the Derivative object. With this implementation, the expression is pre-treated by a function to mask the instances of Derivative objects, then the cse process is applied and finally we do a post-treatment to put back the Derivative expressions in expr.

mleila1312 · 2024-04-10T10:01:57Z

Rebased with more explicit commit message

abhiphile suggested changes Apr 6, 2024

View reviewed changes

abhiphile reviewed Apr 6, 2024

View reviewed changes

sympy/utilities/lambdify.py Outdated Show resolved Hide resolved

moorepants reviewed Apr 6, 2024

View reviewed changes

sympy/utilities/lambdify.py Outdated Show resolved Hide resolved

moorepants reviewed Apr 6, 2024

View reviewed changes

sympy/utilities/lambdify.py Outdated Show resolved Hide resolved

moorepants reviewed Apr 6, 2024

View reviewed changes

sympy/utilities/lambdify.py Show resolved Hide resolved

moorepants reviewed Apr 6, 2024

View reviewed changes

sympy/utilities/lambdify.py Outdated Show resolved Hide resolved

sympy deleted a comment Apr 6, 2024

mleila1312 force-pushed the master branch from c17ad92 to 3fd1ffc Compare April 10, 2024 09:57

mleila1312 force-pushed the master branch from 3fd1ffc to 5e97b9d Compare April 10, 2024 10:00

mleila1312 changed the title ~~Issue #26404 : pre and post treatment of expression when cse=True in lambdify~~ Pre and post treatment of expression when cse=True in lambdify Apr 16, 2024

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Pre and post treatment of expression when cse=True in lambdify #26463

Pre and post treatment of expression when cse=True in lambdify #26463

mleila1312 commented Apr 5, 2024 •

edited

sympy-bot commented Apr 5, 2024 •

edited

moorepants commented Apr 5, 2024 •

edited

mleila1312 commented Apr 5, 2024 •

edited

github-actions bot commented Apr 5, 2024 •

edited

smichr commented Apr 5, 2024

mleila1312 commented Apr 5, 2024

mleila1312 commented Apr 5, 2024 •

edited

moorepants commented Apr 6, 2024

mleila1312 commented Apr 6, 2024 •

edited

moorepants commented Apr 6, 2024

moorepants commented Apr 6, 2024

moorepants Apr 6, 2024

mleila1312 Apr 6, 2024

mleila1312 commented Apr 6, 2024

moorepants commented Apr 6, 2024

moorepants commented Apr 6, 2024

mleila1312 commented Apr 6, 2024

mleila1312 commented Apr 6, 2024 •

edited

mleila1312 commented Apr 10, 2024

Pre and post treatment of expression when cse=True in lambdify #26463

Are you sure you want to change the base?

Pre and post treatment of expression when cse=True in lambdify #26463

Conversation

mleila1312 commented Apr 5, 2024 • edited

Fixes #26404

What is changed

Release Notes

sympy-bot commented Apr 5, 2024 • edited

moorepants commented Apr 5, 2024 • edited

mleila1312 commented Apr 5, 2024 • edited

github-actions bot commented Apr 5, 2024 • edited

smichr commented Apr 5, 2024

mleila1312 commented Apr 5, 2024

mleila1312 commented Apr 5, 2024 • edited

moorepants commented Apr 6, 2024

mleila1312 commented Apr 6, 2024 • edited

moorepants commented Apr 6, 2024

moorepants commented Apr 6, 2024

moorepants Apr 6, 2024

Choose a reason for hiding this comment

mleila1312 Apr 6, 2024

Choose a reason for hiding this comment

mleila1312 commented Apr 6, 2024

moorepants commented Apr 6, 2024

moorepants commented Apr 6, 2024

mleila1312 commented Apr 6, 2024

mleila1312 commented Apr 6, 2024 • edited

mleila1312 commented Apr 10, 2024

mleila1312 commented Apr 5, 2024 •

edited

sympy-bot commented Apr 5, 2024 •

edited

moorepants commented Apr 5, 2024 •

edited

mleila1312 commented Apr 5, 2024 •

edited

github-actions bot commented Apr 5, 2024 •

edited

mleila1312 commented Apr 5, 2024 •

edited

mleila1312 commented Apr 6, 2024 •

edited

mleila1312 commented Apr 6, 2024 •

edited