Hi all! I'm a PhD student at the University of Toronto, working on polynomials (optimization, algorithms, and applications).
As part of my work, I needed to solve systems of polynomials (inequalities and equalities), which SymPy currently lacks. I have therefore been working on such a solver for SymPy. The standard algorithm for this is *cylindrical algebraic decomposition, *which decomposes R^n into "cells" over which each polynomial is sign-invariant. You can then evaluate the system in each cell, which allows you to either get a satisfying point (or, a satisfying point from each true cell), or the algorithm returns that the system is inconsistent. Here is a demo of my code: [image: No description available.] The implementation here returns a single satisfying point if it exists, and otherwise returns an empty list if no satisfying point exists. The algorithm relies heavily on the *Hong projection operator*. This was really the main thing I had to code. The Hong operator relies on taking resultants, reducta, and derivatives -- all of which can already be done in SymPy, so it was just a matter of putting it together (after understanding Hong's paper). I previously posted this on my Twitter and was asked by the SymPy Twitter account to post here. I would love to make this part of SymPy. Note that thus far, the only two implementations of CAD are in Mathematica or a software called QEPCAD. For now, the algorithm works with sample points from each cell (these sample points are generated during the running of the algorithm). I'm thinking that the poly_solve() function could have a parameter like 'return_type' which can take values either 'one' or 'all' -- the former would return a satisfying point from a single cell, the latter would return a satisfying point from all valid cells. Another consideration is that CAD can also be used to construct algebraic formulas that describe all satisfying points. This is called the *solution formula*. Mathematica and QEPCAD implement this. Unfortunately, this step is quite challenging, and I have not fully understood it yet. Pretty much all the information about this step is given in a PhD thesis from the 90s by Christopher Brown (who actually wrote QEPCAD). I'm studying it right now. Overall, I have the following questions: 1. What would be the process to integrate this algorithm into SymPy? 2. Should I wait until I figure out the solution formula step to try to integrate into SymPy? 3. Are there any rules about me eventually writing a paper about my implementation in, for example, the INFORMS Journal of Computing (https://pubsonline.informs.org/page/ijoc/editorial-statement#Software%20Tools)? Thanks! -- You received this message because you are subscribed to the Google Groups "sympy" group. To unsubscribe from this group and stop receiving emails from it, send an email to sympy+unsubscr...@googlegroups.com. To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/d6688cc4-1e09-49fe-800e-a16e3f9e8ba4n%40googlegroups.com.