Hello, On Wed, Jun 10, 2009 at 8:40 AM, faicel<faicelb...@gmail.com> wrote: > > hello > > The question is : > > How to substitute n floating values for the n variables of a > polynomial in Sage, where n is a variable ? > > For example, > > Let P > > P=x3*x6+x1*x6+x4*x5+x2*x5+x1*x4+x2*
Did you leave off something here? For the rest of the example, I'll just use sage: R.<x1,x2,x3,x4,x5,x6> = QQ[] sage: P=x3*x6+x1*x6+x4*x5+x2*x5+x1*x4+x2 > and > > f=x^6+2 sage: f=x^6+2 > with approximate roots > > rr=gp.polroots(f) > =[0.9720806486198328151427283823 + 0.5612310241546864907167665248*I, > 0.9720806486198328151427283823 - 0.5612310241546864907167665248*I, > 0.E-36 + 1.122462048309372981433533050*I, 0.E-36 - > 1.122462048309372981433533050*I, -0.9720806486198328151427283823 + > 0.5612310241546864907167665248*I, -0.9720806486198328151427283823 - > 0.5612310241546864907167665248*I] sage: rr = gp.polroots(f) > We want to substitute in P the n=6 values contained in rr for the > variables xi. > > P(rr[0],rr[1],rr[2],rr[3],rr[4],rr[5]) > is not convenient because we want > P(rr[0],rr[1],rr[2],...,rr[n]) > and to know how to construct the series rr[0],rr[1],rr[2],...,rr[n] ? We can use what's found at http://docs.python.org/tutorial/controlflow.html#unpacking-argument-lists sage: P(*rr) 1.6020411735672693975263336859507182059 - 1.6523546601264078942768391390293986394*I sage: _.parent() GP/PARI interpreter > The substitution may be performed as follows : > ************************ > n=6 > for i in range(1,n+1): > p2=P.subs({x[i]:rr[i]}) > > result: > (x5 + (0.9720806486198328151427283823 > +0.5612310241546864907167665248*I))*x4 + > ((x5 + x3)*x2 + (x3 +(0.9720806486198328151427283823 + > 0.5612310241546864907167665248*I))*x6) > ********************************** > But just the first variable is substituted. A type problem w.r.t. the > coefficients ? (I'll assume that you meant to have "p2 = P" first and that the line in the for was mean to be *p2 = p2.subs({x[i]:rr[i]})") Once you perform the first subs, the type changes to from a Sage polynomial to a GP element (since the things in rr are GP elements, and subs() doesn't work on GP elements: sage: gp(x1) x1 sage: gp(x1).subs({x1:4}) x1 > When the degree n of f is know, We can also substitute using the > following command > RESULT1= P.subs({x[1]:rr[0],x[2]:rr[1],x[3]:rr[2],x[4]:rr[3],x[5]:rr > [4],x[6]:rr[5]}) > > 1.259921049894873164767210607 - 2.524354897 E-29*I > > But in a program the degree of f is n and not necessarily 6. The > problem is to construct You can create the dictionary programmatically like this sage: dict((x[i], rr[i+1]) for i in range(6)) {x6: -0.97208064861983281514272838231160403058 - 0.56123102415468649071676652483958975811*I, x5: -0.97208064861983281514272838231160403058 + 0.56123102415468649071676652483958975811*I, x4: 0.E-55 - 1.1224620483093729814335330496791795162*I, x3: 0.E-55 + 1.1224620483093729814335330496791795162*I, x2: 0.97208064861983281514272838231160403058 - 0.56123102415468649071676652483958975811*I, x1: 0.97208064861983281514272838231160403058 + 0.56123102415468649071676652483958975811*I} dict can take a list (or iterable) of (key, value) tuples. > Moreover, We have used the map fonction of maxima in order to > construct a substitution list for the maxima.subst fonction: > ****************************** > n=6 > v=maxima.map('lambda([rac,i],concat(x,i)=rac)',rr,'makelist(i,i, > 1,n)');v > RESULT2 = gp(maxima.subst(v,P)); > > 1.259921049894873327847796133 + 2.524354897 E-29*I > ************************* > However, result RESULT2 is not the desired result RESULT1 computed > with the subs fonction of Sage > (derived from the subst fonction of maxima) because maxima troncates > the float: > > *********************** > rr[0] > 0.9720806486198328151427283823 + 0.5612310241546864907167665248*I > with maxima rr[0] is troncated: > maxima(rr[0]) > .5612310241546865*I+.9720806486198328 > ********************************************* > > The only solution is to extend the precision of the floats of maxima ? You'd have to use Maxima's bigfloats, but doing it via Maxima is probably not the best way to go about things. > And how in Sage (not in maxima) ? The best solution is to do this all natively in Sage rather than directly in the other systems. sage: R.<x1,x2,x3,x4,x5,x6> = QQ[] sage: S.<x> = QQ[] sage: C = ComplexField(53) sage: P=x3*x6+x1*x6+x4*x5+x2*x5+x1*x4+x2 sage: f = x^6 + 2 sage: rr = f.roots(C) sage: rr[0] (-0.972080648619833 - 0.561231024154686*I, 1) sage: rr = [root for root, multiplicity in rr] sage: P(*rr) -0.342120123672399 - 0.529892611817036*I sage: P.subs(dict(zip(R.gens(), rr))) -0.342120123672399 - 0.529892611817036*I Now at higher precision: sage: C = ComplexField(300) sage: rr = f.roots(C) sage: rr = [root for root, multiplicity in rr] sage: P(*rr) -0.342120123672396232759123078672489855294300914423752822089952855092800116533273243051072780 - 0.529892611817034912843306089350219123209667783433307563983436300295328107274293118489368470*I sage: P.subs(dict(zip(R.gens(), rr))) -0.342120123672396232759123078672489855294300914423752822089952855092800116533273243051072780 - 0.529892611817034912843306089350219123209667783433307563983436300295328107274293118489368470*I Or using exact roots via QQbar: sage: C = QQbar sage: rr = f.roots(C) sage: rr = [root for root, multiplicity in rr] sage: P(*rr) -0.3421201236723962? - 0.5298926118170349?*I sage: P.subs(dict(zip(R.gens(), rr))) -0.3421201236723962? - 0.5298926118170349?*I Let me know if this helps. --Mike --~--~---------~--~----~------------~-------~--~----~ To post to this group, send email to sage-support@googlegroups.com To unsubscribe from this group, send email to sage-support-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-support URLs: http://www.sagemath.org -~----------~----~----~----~------~----~------~--~---
