On Feb 16, 2011, at 4:35 PM, Mateusz Paprocki wrote: > Hi, > > On 16 February 2011 12:49, Aaron S. Meurer <[email protected]> wrote: > > On Feb 16, 2011, at 10:43 AM, Mateusz Paprocki wrote: > >> Hi, >> >> On 16 February 2011 04:49, Vinzent Steinberg >> <[email protected]> wrote: >> Hi, >> >> On 16 Feb., 10:37, jeal2 <[email protected]> wrote: >> > Hello, >> > >> > I'm trying to use Python to solve (symbolically, for x, y and a0) the >> > following set of three simultaneous equations: >> > lam+2*y-a0*(1 - x/2)*x-0.005*x/2*x==0 >> > a0*(1 - x/2)*x-1*y-0.743436700916726*y==0 >> > x+y-conc==0 >> > >> > I don't see any reason it should be mathematically difficult - just >> > complicated and tedious - and Mathematica can do it without any >> > trouble. But when I give it to Python it just starts working and >> > doesn't finish. Do you have any idea why? >> > >> > >>> from sympy import * >> > >>> lam=Symbol('lam') >> > >>> a0=Symbol('a0') >> > >>> conc=Symbol('conc') >> > >>> x=Symbol('x') >> > >>> y=Symbol('y') >> > >> > # this works ok:>>> solve([lam+2*y-a0*(1 - x/2)*x-0.005*x/2*x], [x]) >> > >> > [(-(-4*(-4*y - 2*lam)/(0.005 - a0) + 4*a0**2/(0.005 - a0)**2)**(1/2)/2 >> > - a0/(0.005 - a0),), ((-4*(-4*y - 2*lam)/(0.005 - a0) + 4*a0**2/(0.005 >> > - a0)**2)**(1/2)/2 - a0/(0.005 - a0),)] >> > >> > # but this doesn't >> > >> > >>> solve([lam+2*y-a0*(1 - x/2)*x-0.005*x/2*x, a0*(1 - >> > >>> x/2)*x-1*y-0.743436700916726*y, x+y-conc], [x, y,a0]) >> > >> > If you have any ideas I would be very grateful! >> > Thanks. >> >> With current master I get this traceback: >> >> In [4]: eqs = [lam+2*y-a0*(1 - x/2)*x-0.005*x/2*x, a0*(1 - x/ >> 2)*x-1*y-0.743436700916726*y, x+y-conc] >> >> In [5]: sym = [x,y,a0] >> >> In [6]: solve(eqs, sym) >> --------------------------------------------------------------------------- >> ExactQuotientFailed Traceback (most recent call >> last) >> >> /home/one/src/sympy/<ipython console> in <module>() >> >> /home/one/src/sympy/sympy/solvers/solvers.pyc in solve(f, *symbols, >> **flags) >> 379 soln = solve_linear_system(matrix, *symbols, >> **flags) >> 380 else: >> --> 381 soln = solve_poly_system(polys) >> 382 >> 383 # Use swap_dict to ensure we return the same type >> as what was >> >> >> /home/one/src/sympy/sympy/solvers/polysys.pyc in >> solve_poly_system(system, *gens) >> 104 raise TypeError("expected iterable container, got %s" >> % system) >> 105 >> --> 106 solutions = solve_reduced_system(system, gens, entry=True) >> 107 >> 108 if solutions is None: >> >> /home/one/src/sympy/sympy/solvers/polysys.pyc in >> solve_reduced_system(system, gens, entry) >> 55 def solve_reduced_system(system, gens, entry=False): >> 56 """Recursively solves reduced polynomial systems. """ >> ---> 57 basis = groebner(system, gens, polys=True) >> 58 >> 59 if len(basis) == 1 and basis[0].is_ground: >> >> /home/one/src/sympy/sympy/polys/polytools.pyc in groebner(F, *gens, >> **args) >> 2373 F[i] = sdp_from_dict(f.rep.to_dict(), order) >> 2374 >> -> 2375 G = sdp_groebner(F, lev, order, dom) >> 2376 >> 2377 G = [ Poly(DMP(dict(g), dom, lev), *gens) for g in G ] >> >> /home/one/src/sympy/sympy/polys/groebnertools.pyc in sdp_groebner(F, >> u, O, K) >> 620 >> 621 p = sdp_mul_term(p, (monomial_div(M, p_LM), >> K.quo(K.one, sdp_LC(p, K))), u, O, K) >> --> 622 q = sdp_mul_term(q, (monomial_div(M, q_LM), >> K.quo(K.one, sdp_LC(q, K))), u, O, K) >> 623 >> 624 h = normal(sdp_sub(p, q, u, O, K), G) >> >> /home/one/src/sympy/sympy/polys/algebratools.pyc in quo(self, a, b) >> 462 """Quotient of `a` and `b`, implies `__floordiv__`. >> """ >> 463 if a % b: >> --> 464 raise ExactQuotientFailed('%s does not divide %s >> in %s' % (b, a, self)) >> 465 else: >> 466 return a // b >> >> ExactQuotientFailed: DMP([[Fraction(-1471256443825, >> 2309382115472961)], [Fraction(294251288765000, 2309382115472961), >> Fraction(0, 1)]], QQ) does not divide DMP([[Fraction(1, 1)]], QQ) in >> QQ[conc,lam] >> >> So it seems related to polynomials. Mateusz, any idea why this fails? >> Maybe the operation should rather be done over RR[conc,lam]? >> >> Yes, recently I added support for domains of this kind, so now (polys12) >> solve() gives: >> >> DomainError: can't compute a Groebner basis over RR >> >> for this example. > > Can this ever be implemented? Or are there numerical algorithms that can > solve systems of polynomial equations (like an extension of nsolve())? > > This can be implemented in two ways (still using Groebner bases): > > 1) Use "approximate arithmetics". In this case the domain would be RR(eps), > not simply RR, where eps would be the precision of computations. This way > Groebner basis algorithm would terminate, because we would be able to solve > the zero equivalence problem in the new type of domain. I'm not sure how > complex this change would be: surely we would have to change how RR works, > and possibly robustify groebnertools.py. Anyway, this would work in a similar > way to how CoefficientDomain->InexactNumbers[n] work in Mathematica. > > 2) Implement a Groebner bases-based solver using interval/rectangle > arithmetics. Some primitive elements of such an algorithm were already > implemented, but it's still a lot of work (I guess).
The first one sounds like it would be a nice thing to have in general for RR, not just for Groebner bases. The Buchberger algorithm is not the only algorithm that relies on zero equivalence checking (for example, the division algorithm, for which the Groebner basis can be considered a generalization of, relies one it). Perhaps the second one has advantages too (is it maybe faster?). If so, maybe we should aim for both. Aaron Meurer > > > Aaron Meurer > >> >> >> Vinzent >> >> Mateusz >> >> -- >> You received this message because you are subscribed to the Google Groups >> "sympy" group. >> To post to this group, send email to [email protected]. >> To unsubscribe from this group, send email to >> [email protected]. >> For more options, visit this group at >> http://groups.google.com/group/sympy?hl=en. > > > -- > You received this message because you are subscribed to the Google Groups > "sympy" group. > To post to this group, send email to [email protected]. > To unsubscribe from this group, send email to > [email protected]. > For more options, visit this group at > http://groups.google.com/group/sympy?hl=en. > > Mateusz > > -- > You received this message because you are subscribed to the Google Groups > "sympy" group. > To post to this group, send email to [email protected]. > To unsubscribe from this group, send email to > [email protected]. > For more options, visit this group at > http://groups.google.com/group/sympy?hl=en. -- You received this message because you are subscribed to the Google Groups "sympy" group. To post to this group, send email to [email protected]. To unsubscribe from this group, send email to [email protected]. For more options, visit this group at http://groups.google.com/group/sympy?hl=en.
