Hacker News Board Hacker News, browsed like a classic message board
Live data from Hacker News

Doing Symbolic Math with SymPy

lwn.net

Post reply on HN
Print view View original

61–70 of 81 posts

Re: Doing Symbolic Math with SymPy

#61

Post by mehrdadn »

ogogmad wrote: post #54
BeetleB wrote: post #48

Earlier quoted context omitted.

Yikes those are bad. As a non-user interested in it, I must confess that seeing those issues (especially the 2nd one) makes me not want to use it! I do get that it is a hard problem, but I would then recommend they not have an option for "real" if they can't get it right. Users will expect it to work (unless the docs point out it is unreliable).

> I would then recommend they not have an option for "real" if they can't get it right. Users will expect it to work (unless the docs point out it is unreliable). You should be suspicious of using floats and expecting the solver to detect if a solution is real. Consider the polynomial $x^2 + C$ for $C \approx 0$. Are all of its roots real or complex? The problem is that even a very small change in C, maybe caused by…

> You should be suspicious of using floats and expecting the solver to detect if a solution is real. Consider the polynomial $x^2 + C$ for $C \approx 0$. Are all of its roots real or complex? The problem is that even a very small change in C, maybe caused by rounding errors, can easily change a root from real to not real.

C \approx 0?

The imaginary parts in the first bug are nowhere close to zero.

Re: Doing Symbolic Math with SymPy

#62

Post by ogogmad »

mehrdadn wrote: post #61
ogogmad wrote: post #54

Earlier quoted context omitted.

> I would then recommend they not have an option for "real" if they can't get it right. Users will expect it to work (unless the docs point out it is unreliable). You should be suspicious of using floats and expecting the solver to detect if a solution is real. Consider the polynomial $x^2 + C$ for $C \approx 0$. Are all of its roots real or complex? The problem is that even a very small change in C, maybe caused by…

> You should be suspicious of using floats and expecting the solver to detect if a solution is real. Consider the polynomial $x^2 + C$ for $C \approx 0$. Are all of its roots real or complex? The problem is that even a very small change in C, maybe caused by rounding errors, can easily change a root from real to not real. C \approx 0? The imaginary parts in the first bug are nowhere close to zero.

I was actually responding to the following claim:

> I would then recommend they not have an option for "real" if they can't get it right. Users will expect it to work (unless the docs point out it is unreliable).

It can't work reliably if the coefficients of the polynomial consist of floats. At least not in general. This point is a more fundamental one than your example.

Also, regarding your particular example: When solving cubics using the cubic formula, it is often the case that you will end up computing with complex numbers even if the result ends up real. If you introduce floats into this, it's impossible to guarantee that the imaginary parts will cancel out completely. Instead, they may cancel up to 10^-21 or something like that, which is precisely what happens in your example.

Re: Doing Symbolic Math with SymPy

#63

Post by mehrdadn »

ogogmad wrote: post #62
mehrdadn wrote: post #61

Earlier quoted context omitted.

> You should be suspicious of using floats and expecting the solver to detect if a solution is real. Consider the polynomial $x^2 + C$ for $C \approx 0$. Are all of its roots real or complex? The problem is that even a very small change in C, maybe caused by rounding errors, can easily change a root from real to not real. C \approx 0? The imaginary parts in the first bug are nowhere close to zero.

I was actually responding to the following claim: > I would then recommend they not have an option for "real" if they can't get it right. Users will expect it to work (unless the docs point out it is unreliable). It can't work reliably if the coefficients of the polynomial consist of floats. At least not in general. This point is a more fundamental one than your example. Also, regarding your particular example: When…

You're missing my point entirely. They don't need to solve this by keeping the variable real at every step or anything like that. What they need to do is to ensure the final solutions are real. Which means they need to solve it in the complex domain entirely and then filter out non-real solutions at the end, where they have the final magnitudes of the real and imaginary components. The only place that would run into trouble is where the imaginary components are small but nonzero, where you'd just expect it to use some tolerance threshold as the cutoff, but that obviously isn't the case in the first example (hence my comment).

There's no reason that the algorithm must fail on univariate real polynomials if it can't work equally well on arbitrary functions with arbitrary domains. It's trivial for a symbolic engine to recognize a univariate polynomial. And practically speaking it's the job of a computer algebra system to recognize different well-known classes of functions and treat them as well as it can. What I suspect likely happened here is they've been just so busy dealing with the more general problems (which are far more difficult) that they haven't gotten around to implementing many better methods for the special cases (like polynomials). Which I sympathize with, as I'd have virtually no clue how to solve some of the more general problems to begin with.

Re: Doing Symbolic Math with SymPy

#64

Post by dang »

If curious see also

a few weeks ago https://news.ycombinator.com/item?id=25254648

2020 https://news.ycombinator.com/item?id=23747133

2019 https://news.ycombinator.com/item?id=20287486

Lots in 2014:

https://news.ycombinator.com/item?id=8243452

https://news.ycombinator.com/item?id=8109265

https://news.ycombinator.com/item?id=7214783

https://news.ycombinator.com/item?id=7145219

2012 https://news.ycombinator.com/item?id=4463471

2011 https://news.ycombinator.com/item?id=3141699

Re: Doing Symbolic Math with SymPy

#65

Post by mehrdadn »

bee_rider wrote: post #56
mehrdadn wrote: post #16

The trouble with SymPy is it's, well, buggy. I tried it a few years ago, and as soon as I got serious, I quite quickly ran into problems that I reported, some of which I now see they apparently still haven't gotten around to addressing. [1] [2] Symbolic math is hard ; they have my sympathies. I don't think I could do better. But as long as bugs like these exist, it's going to be hard to convince people to switch away…

Have you tried out SageMath? It seemed to work pretty well when I tried it out, but I never got serious really. My needs for symbolic math programming are pretty limited.

I haven't, though I did download it. I was trying to work on an algorithm rather than solve specific problem instances, so I ended up just working around it by using different problem instances that didn't result in these errors.

Re: Doing Symbolic Math with SymPy

#66

Post by jupdarter »

I’ve had similar problems as others posting here have. To be honest when I’ve had to do symbolic math at work (usually some combination of computational geometry and PDEs), I typically fall back to maxima or Mathematica (if my collaborators also have licenses). Those systems are quite solid and very mature when it comes to symbolic work, and I never saw any value in using SymPy since I can get to maxima from Python via Sage. I will admit though: I tend to not do that since I’m usually pretty happy in maxima or Mathematica directly and don’t see a reason to complicate things with a layer of Python. I taught a course a couple years back where I tried to use sympy since we were using a bunch of other Python tools (numpy, scipy, etc...), and the sympy part felt very awkward. I think I remember showing the students maxima and Mathematica to convince them that computers really could do complex symbolic work and not to get too turned off by the Python tools.

Re: Doing Symbolic Math with SymPy

#67

Post by NoThisIsMe »

fantod wrote: post #38
leephillips wrote: post #18

Earlier quoted context omitted.

This is problem I have with Python in general; that even if you prefer to use it in a functional style, most libraries are written by real Python programmers, and using them will force you to grapple with the object-oriented inversion of common sense ( https://lee-phillips.org/pythonhate/ ). My favorite example is: ','.join(['a', 'b'])

Yeah, the Python language favors an imperative style due to the absence of features necessary for writing more functionally, like lambdas that can do more than return an expression, tail call optimization, and lexical scope. This last one, additionally, tends to make people use classes for essentially no other reason than that closures don't work natively. > My favorite example is: ','.join(['a', 'b']) I've always fo…

Closures do work natively in Python though. Lexical scope does not prevent that. Most languages, including functional languages, use lexical scope these days.

Re: Doing Symbolic Math with SymPy

#68

Post by carapace »

I did a fun thing with Sympy for Advent of Code 2017 December 3rd, here's a Jupyter notebook if anyone's interested: https://osdn.net/projects/joypy/scm/hg/Joypy/blobs/tip/docs/...

The puzzle involves a weird "spiral" RAM:

> You come across an experimental new kind of memory stored on an infinite two-dimensional grid.

> Each square on the grid is allocated in a spiral pattern starting at a location marked 1 and then counting up while spiraling outward.

The problem is to compute "the Manhattan Distance between the location of the data and square 1" for any given square.

While working on this I came up with the following equation:

    from sympy import floor, lambdify, solve, symbols

    k = symbols('k')

    E = 2 + 8 * k * (k + 1) / 2
Sympy figures out that this reduces to E = 4k(k + 1) + 2

I needed a function to solve for k given some n... Sympy can do that. Take a new symbol `n`, subtract it from the equation `E`, 0 = 4k(k + 1) + 2 - n and solve for `k`. There are two solutions because the equation is quadratic so it has two roots.

    n = symbols('n')

    g, f = solve(E - n, k)
In the context of the puzzle we only care about the larger root:

    (sqrt(n - 1) / 2 - 0.5) + 1
For reasons, I need to take the floor and add 1. Then Sympy can lambdify it and create a fast Python function to compute `k`, given `n`:

    F = lambdify(n, floor(f) + 1)

    for n in (9, 10, 25, 26, 49, 50):
        print(n, int(F(n)))

    9 1
    10 2
    25 2
    26 3
    49 3
    50 4
I'm sorry that it's hard to follow, I'm just excited to share how cool it was to let the computer do the symbolic math for me.

Re: Doing Symbolic Math with SymPy

#69

Post by snicker7 »

ChrisRackauckas wrote: post #34
thatcherc wrote: post #32

Worth noting that Julia's SymPy binding [1] is pretty pretty nice to work with too. If anyone's looking for big Julia project, I think a symbolic math package written fully in Julia would be a really exciting development. As far as I know, there isn't one yet. The better-known symbolic math packages for Julia still use bindings to C++ (SymEngine.jl [2]) or Python (SymPy.jl, Symata.jl [3]). [1] - https://github.com/Ju…

ModelingToolkit.jl is a symbolic math package written fully in Julia (with a bunch of extra symbolic-numerics features) https://github.com/SciML/ModelingToolkit.jl https://mtk.sciml.ai/dev/tutorials/symbolic_functions/ It's more like SymEngine in terms of completeness right now, though there's a good amount of simplification and equation solving built in. It's still growing, it's not at SymPy yet, but it's moving fas…

What is great about ModelingToolkit.jl is how its used in practical ways for other packages. E.g. NeuralPDE.jl.

Compared to SymPy, I feel that it is less of a "how do I integrate this function" package and more about "how can I build this DSL" framework.

https://github.com/SciML/NeuralPDE.jl

Re: Doing Symbolic Math with SymPy

#70

Post by BeetleB »

mehrdadn wrote: post #58
BeetleB wrote: post #55

Earlier quoted context omitted.

That one seems a bit more reasonable, though. When doing numerical computation (as opposed to symbolic), you should never assume the algorithm "just works". You need to know the quirks of the underlying method, etc. And I believe most root finding numerical algorithms will be sensitive to the initial guess. Surprised they haven't added the warning and closed the bug.

Not at all. I've worked on numerical stuff myself; unlike the symbolic bugs, I have little sympathy for these particular numerical errors. It would be far more reasonable if it was convoluted composition of functions, or some kind of bizarre numerical edge case. However: (a) Univariate quadratics are extremely well-understood and extremely important. For example, the convergences of optimization algorithms (this is r…

It's been a while since I did root finding in SciPy, but in principle I do not agree, and it doesn't seem like you're proposing an actual solution. You can make the case that that root finding algorithm shouldn't be the default (if it is). You can argue for a warning, as mentioned in the bug, but what else? I personally prefer less magic in the algorithms (e.g. detect if it is a quadratic and pick a suitable algorithm), and more information in the docs about the method being used so the user is aware.

Actually, I can see the argument both ways. To each his own. You're not fundamentally wrong, but you're also not fundamentally right either :-)

I personally should think their returned result should display the value of the function evaluated at their root so the programmer can check it against his own tolerance.

Post reply on HN