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.

Reply via email to