Intro
We have been working on a recent project and have been looking for a suitable system to calculate some values. SymPy was recommended as being a rich mathematical library. However, we have been unable to make it "work" with our project.
The issue we have been struggling with specifically is that many of the values we would be using have been rounded numerous times and are likely susceptible to float errors. To work around this on a previous project, we used Interval Arithmetic for JavaScript to fairly effective use. mpmath for Python appears to be similar, but SymPy not only uses mpmath, but also offers other potentially useful functions we may need in the future.
Problem
A sample equation that we have been working on lately is a = b * (1 + c * d * e) and we are looking to solve for e when all other variables are known. However, some of the variables need to be represented as a range of values as we don't know the exact value, but a small range.
Code
from sympy import * from sympy.sets.setexpr import SetExpr a, b, c, d, e = symbols('a b c d e') b = 40 c = 1 d = 0.1 a = SetExpr(Interval(45.995, 46.005)) equ = Eq(b * (1 + c * d * e), a) solveset(equ, e) ValueError: The argument '45.995*I' is not comparable. This was just the latest attempt, but I have tried setting domains, setting inequalities for symbols, using AccumBounds, and numerous other solutions but I can't help but think that we have completely overlooked something simple.
Solution
It appears that that using one interval is doable with the code provided by the selected answer, but it doesn't extend to multiple symbols requiring intervals or ranges of values. It appears we will be extending the mpmath library to support additional interval functions.
1 Answer
There is an intervalmath module in SymPy in the plotting module for some reason. It doesn't subclass Basic though so can't be used directly in an expression. We can however use lambdify to substitute it into an expression as
from sympy import * from sympy.plotting.intervalmath import interval b = 40 c = 1 d = 0.1 a, e = symbols('a, e', real=True) equ = Eq(b * (1 + c * d * e), a) sol_e, = solveset(equ, e) f_e = lambdify(a, sol_e) int_a = interval(45.995, 46.005) int_e = f_e(a=int_a) print(int_e) This gives
[1.498750, 1.501250] I don't think the intervalmath module is used much though so there's a good chance it might not fully work in your real problem.
Using sets is probably a better approach and it seems that imageset can do this:
In [16]: set_e = imageset(Lambda(a, sol_e), Interval(45.995, 46.005)) In [17]: set_e Out[17]: [1.49875, 1.50125] I'm not sure how well this works with more than one symbol/interval though.
EDIT: For completeness I'm showing how you would use intervalmath with more than one interval:
from sympy import * from sympy.plotting.intervalmath import interval b = 40 d = 0.1 a, c, e = symbols('a, c, e', real=True) equ = Eq(b * (1 + c * d * e), a) sol_e, = solveset(equ, e) f_e = lambdify((a, c), sol_e) int_a = interval(45.995, 46.005) int_c = interval(0.95, 1.05) int_e = f_e(a=int_a, c=int_c) print(int_e) That gives
[1.427381, 1.580263] 
