I am solving a polynomial which arises from plotting titration cures in
chemistry. The rule of signs suggests it has one positive root. Find_root
seems to find it. Solve with poly_solve=true does not. Instead it gives 4
complex roots, which don't appear to satisfy the equation. They do not even
form two conjugate pairs. Here is some code which illustrates the problem:
var('Kb','b0','K1','Ks','a0')
fx(x)=-Kb*x^4 - (b0*Kb + K1*Kb + Ks)*x^3 - (b0*K1*Kb - a0*K1*Kb + K1*Ks -
Kb*Ks)*x^2 + (a0*K1*Ks + K1*Kb*Ks + Ks^2)*x + K1*Ks^2
fx(x)=fx.subs(K1=10^-2.15,Ks=10^-14,Kb=10^5.0,a0=0.1,b0=0)
s=find_root(fx,0,0.025,1E-14)
show(s)
s1=solve(fx,x,to_poly_solve=true)
show('roots from solve:
',N(s1[0].rhs(),3));show(N(s1[1].rhs(),3));show(N(s1[2].rhs(),3));show(N(s1[3].rhs(),3))
show('Substituting first root from solve: ',N(fx(s1[0].rhs()),10))
show('Substituting first root from find_root: ',N(fx(s),10))
plot(fx,0.001,0.03)
Can someone shed light on this?
--
You received this message because you are subscribed to the Google Groups
"sage-support" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to sage-support+unsubscr...@googlegroups.com.
To post to this group, send email to sage-support@googlegroups.com.
Visit this group at https://groups.google.com/group/sage-support.
For more options, visit https://groups.google.com/d/optout.
