Ask Your Question

Solving a quartic equation

asked 2015-05-31 21:51:26 -0500

I'm attempting to rearrange an equation from an answer on the Mathematics StackExchange.

The answer given is this equation:

$$L^2 = (-ab(t)+p)^2-\left(\frac{(-ab(t)+p).(cd(t)-ab(t))}{(cd(t)-ab(t))^2}(cd(t)-ab(t))\right)^2$$

Where $a$, $b$, $c$, $d$, and $p$ are known 2D points, $L$ is a known length, and $t$ is an unknown scalar. $ab(t)$ indicates interpolation between $a$ and $b$.

I am interested in rearranging this to solve for $t$. Here's what I've tried in Sage:

def sqr(var): return var.dot_product(var)
var('ax bx cx dx px ay by cy dy py t L')
a = vector([ax, ay])
b = vector([bx, by])
c = vector([cx, cy])
d = vector([dx, dy])
p = vector([px, py])
g = a - at + bt
h = c - ct + dt
u = p - g
v = h - g
eq = L^2 == sqr(u) - sqr((u.dot_product(v)/sqr(v)) * v)

At the solve step I have observed it to sit for quite a while without producing a result. Two questions:

  1. Am I inputting the problem correctly?

  2. Is there any way to know if this is likely to terminate in a reasonable time? I have no idea what the solver looks like under the hood, and wouldn't want to wait for some O(n!) calculation to terminate :)

edit retag flag offensive close merge delete


The problem per se seems quite involved. The r.h.s when expanded displays 3 million characters:

sage: len(str(eq.rhs().expand()))
rws gravatar imagerws ( 2015-06-01 01:04:35 -0500 )edit

1 answer

Sort by ยป oldest newest most voted

answered 2015-06-02 08:02:15 -0500

updated 2015-06-02 08:32:42 -0500

There is a simpification possible. Instead of expanding the expression (which makes for a nice Pynac benchmark) we transform the equation into a fraction and get the numerator:

sage: g = a - a*t + b*t
sage: h = c - c*t + d*t
sage: u = p - g
sage: v = h - g
sage: eq = L^2 == sqr(u) - sqr((u.dot_product(v)/sqr(v)) * v)
sage: rhs = eq.rhs()
sage: num = (rhs.rational_simplify()-L^2).numerator()
sage: solve(num, t)

This is much smaller but unfortunately still takes forever to solve. Using collect and coefficient you can let Sage show you the coefficient of t^0, t, t^2, t^3, and t^4 to get an idea of what's involved and which equations determine the solution process of the quartic. For example this is the constant coefficient:

sage: num.coefficient(t,0)
-L^2*ax^2 - L^2*ay^2 + 2*L^2*ax*cx - L^2*cx^2 + ay^2*cx^2 + 2*L^2*ay*cy - 2*ax*ay*cx*cy - L^2*cy^2 + ax^2*cy^2 - 2*ay^2*cx*px + 2*ax*ay*cy*px + 2*ay*cx*cy*px - 2*ax*cy^2*px + ay^2*px^2 - 2*ay*cy*px^2 + cy^2*px^2 + 2*ax*ay*cx*py - 2*ay*cx^2*py - 2*ax^2*cy*py + 2*ax*cx*cy*py - 2*ax*ay*px*py + 2*ay*cx*px*py + 2*ax*cy*px*py - 2*cx*cy*px*py + ax^2*py^2 - 2*ax*cx*py^2 + cx^2*py^2

You could also get from Sage the general solution to the quartic via:

sage: (a,b,c,d,e)  = var('a,b,c,d,e')
sage: sol = solve(a*x^4+b*x^3+c*x^2+d*x+e,x)

This returns a vector of 4 solutions. You could then simply substitute your coefficients for a,b,c,d,e in one of the solutions:

sage: sol1 = sol[0]
sage: sol1 = sol1.subs(e==num.coefficient(t,0)).subs(d==num.coefficient(t,1)).subs(c==num.coefficient(t,2)).subs(b==num.coefficient(t,3)).subs(a==num.coefficient(t,4))
sage: len(str(sol1))

This solution however is still too big to handle. It needs 500 pages to print! Simplifying is out of the question.

edit flag offensive delete link more

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower


Asked: 2015-05-31 21:51:26 -0500

Seen: 125 times

Last updated: Jun 02 '15