Recommendations for binomial system solver I am interested in solving binomial systems of the form
$$
  \begin{cases}
    a_1 x_1^{d_{11}} x_2^{d_{12}} \cdots x_n^{d_{1n}} +
    b_1 x_1^{d_{11}} x_2^{d_{12}} \cdots x_n^{d_{1n}} &= 0 \\\\
    \vdots &\vdots \\\\
    a_m x_1^{d_{m1}} x_2^{d_{m2}} \cdots x_n^{d_{mn}} +
    b_m x_1^{d_{m1}} x_2^{d_{m2}} \cdots x_n^{d_{mn}} &= 0
  \end{cases}
$$
where the exponents may be negative. I.e., each equation has exactly two terms. For example
$$  
  \begin{cases}
    3 x_1^{2} x_2^{-5} + 4 x_1^{-1} x_2^{6} &= 0 \\\\
    2 x_1^{3} x_2^{5}  - 7 x_1^{2}  x_2^{4} &= 0
  \end{cases}
$$
I am looking for some software packages that can solve these systems over the nonzero complex numbers. Here is what I mean by "solve":


*

*When the solution set consists of isolated points. Locate the points.

*When the solution set has positive dimension, find the number of the components and the dimension. I.e., an irreducible decomposition of the variety defined by the system. Ideally we should also get the parametrization of each component.

*Multiplicity information would be nice but not necessary.


Or more technically, I'd like to decompose a variety defined by a binomial system into (irreducible) toric varieties and enumerate the characters.
Binomials.m2, a Maclaulay2 package (http://thomas-kahle.de/bpd.html) seems to be the best.  But when applied to large systems with 50 or more variables, Binomials.m2 simply does not terminate within any reasonable amount of time.
I am familiar with the algorithms behind, but before I write my own, I would like to know if there is anything even better out there. Or maybe I'm using Binomials.m2 incorrectly?
 A: Let $I$ be the ideal generated by your binomial equations after clearing denominators.
As a general rule with binomial equations, disregard coefficients in the first run (see here for why and how: http://arxiv.org/abs/1107.4699 ).  I would first solve the system where all coefficients are set to be one, which means the equations should be of the form monomial minus monomial = 0. Once you have solved this, you can deform the result to a solution of your original problem.
Since you want to disregard any components contained in coordinate hyperplanes ("only nonzero solutions"), you are actually looking for the saturation $I : (\prod_i x_i)^\infty$.  If you have it, then you are basically done since the saturation is a lattice ideal.  In characteristic zero it is unmixed and has trivial primary decomposition given by extending characters (Corollary 2.2 in "Binomial ideals" [Eisenbud,Sturmfels, Duke Journal 84, 1996].
To compute the saturation $I : (\prod_i x_i)^\infty$ the fastest software is 4ti2 ( http://www.4ti2.de/ ). There the command markov does what you want to do.  You would have to specify the exponent vectors of your equations in a '.lat' file and then run markov.  Note that you can also do this from Macaulay2 using the FourTiTwo package.
To iterate, compared to computing generators of the saturation the rest is trivial. For instance finding parameterizations of the components of the saturation is basically just a Smith normal computation as in Chapter 7 of of the book "Combinatorial Commutative Algebra" [Miller/Sturmfels, Springer GTM].
If you like you can e-mail me a specific example.
