Booki logo

Algebra: An Algorithmic Treatment

Appendix D. Exact Solutions of Polynomial Equations

Introduction

For thousands of years the solution of polynomial equations was one of the most important and productive problems in algebra, and indeed in all of mathematics. Apart from the question of solution methods, the study of quadratic, cubic, and quartic equations and their geometric equivalents led to great advances in the concept of number over a period of more than two millennia, including irrational, negative, and complex numbers, and the limitations of numbers constructible with ruler and compass. Later on, the problem of quintic and higher order equations led to what has been called higher algebra or abstract algebra, that is, to the study of structures such as groups, rings, fields, vector spaces, and algebras of many other kinds, further greatly expanding the concept of number, and many other concepts. All of these structures and many others were then gathered together into category theory, which studies all of the mappings within and between all such objects, providing insight not only into the algebraic structures themselves, but also into algebraic geometry, the study of geometric objects by means of algebraic structures defined on them.

Furthermore, the solution of quadratic equations turned out to be fundamental in elementary physics in the period from Galileo to Newton. (Going beyond the elementary level requires calculus, as explained, for example, in Ken Iverson's book Elementary Analysis, which can be used as a continuation from this study of algebra.) For example, in elastic collisions, where both momentum and energy are conserved, the conservation laws are respectively linear and quadratic in form, so that the solution which preserves both reduces to the two solutions of a quadratic equation. If m and v0 are two-element vectors of masses and initial velocities, we must find a vector v1 such that:

   (+/m*v0)=(+/m*v1)
   (+/m**:v0)=(+/m**:v1) 

One solution, v0, represents the state before the collision, and the other solution represents the state after the collision.

The solution of polynomial equations is fundamental to computer graphics, computer-aided design, and robotics.

The original problem of solving polynomials was to find solutions that could be expressed exactly in terms of the five basic algebraic operations, + - * % %:, starting from integers. The body of this textbook has been concerned only with numeric solutions,  but this Appendix reviews the classical results and methods for solving quadratic, cubic, and quartic equations.

General Considerations

We know that every polynomial is a product of linear terms, of the form */x-r possibly multiplied by a constant that has no effect on the roots. We can divide the coefficient vector c of a polynomial by its last element _1{c, giving a new coefficient vector for a polynomial with 1 for the coefficient of the highest power of x present. We call this a monic polynomial.

The monadic case of the polynomial verb p. allows us to convert between coefficient vectors and a representation of a polynomial with boxed multiplier and list of roots. For example:

   p. 1;2 3
6 _5 1
   p. 6 _5 1
┌─┬───┐
│1│3 2│
└─┴───┘
   6 _5 1 p.i.5
6 2 0 0 2

Each linear polynomial of the form x-r has a root r, so the product polynomial has a root for each of the r values. That is, at each of those values, at least one term has the value 0, so the product also has the value 0. These are the only roots, because at any other value of x, the factors are all non-zero and therefore the product is non-zero.

Furthermore, we know that the constant term in a polynomial */x-r is the product of all of the roots, by construction. The second-highest order term is the sum of all of the roots, again by construction. In the case of a quadratic polynomial, we saw in Chapter 14 that the coefficients of the polynomial are (r0*r1),(-r0+r1),1, possibly multiplied by a constant. Thus the polynomial with roots 2 3 has coefficients 6 _5 1, as above. A cubic polynomial has a term composed of the sum of products of pairs of roots, that is, (r1*r2)+(r0*r2)+(r0*r1). A quartic polynomial has a term composed of sums of products of triples of roots, and so on.

These terms represent functions of the roots with the special property that interchanging roots does not change the value of the function, unlike, say r0-r1, which is negated by interchanging r0 and r1. Functions with this property are called symmetric. Their study was fundamental in the proof that quintic and higher-order polynomials do not have solutions using only the functions + - * % %:, known as solutions in radicals.

Quadratic Equations

Although many mathematicians in Babylonia, Egypt, Greece, India, and China developed partial methods for solving quadratic equations, the general method to be described here is due to Muhammad ibn Musa al-Khwarizmi of Persia in the 9th century. It was introduced into Europe in the 12th century, first appearing in a Hebrew work by Abraham bar Ḥiyya ha-Nasi, then in Latin translations of al-Khwarizmi and ha-Nasi. However, the full solution including negative and complex roots was not recognized until the 19th century.

Given a three-element coefficient vector c from which we can make the polynomial p=.c&p., where the last element of c is not 0, how can we go about solving the equation 0=p x for its two roots, r0 and r1? We can begin by dividing by the last element of c to give an equivalent monic polynomial.

Graph

None of these methods is entirely satisfactory. What we would like is an expression for the roots in terms of the coefficients of the polynomial. There is a particular case where this is very easy, the case where the two roots are equal. We have

   (x-r)*(x-r)
   (x^2)+(-2*r*x)+(r^2)

Setting 0=.*:(x-r) and solving gives us

   0=*:x-r
   (%:0)=x-r
   r+0=x
   r=x

as desired. Graphically, we can see that this is the case where the curve just touches the x axis at the point (r,0), as in the following graph:

Graph

Usually the graph of a quadratic polynomial in real numbers intersects the x axis in two distinct points, or not at all. But in those cases, adding a constant to the polynomial would move the graph up (for a positive constant) or down (for a negative constant), so there must be an intermediate value that would move the graph to just touch the x axis. Let us call that number n. Then adding n,0 0 to c, or setting p1=.n+p, has this result in the graph.

Graph

We now know that p1 is the square of x-r for some still unknown r, but we also know that

   (1{c)=2*r
   ((1{c)%2)=r

Using this r, we can proceed as follows:

   0=p1 x
   n=n+p x
   (%:n)=x-r
   (r+%:n)=x

Now in fact n has two square roots, one of which is returned by the %: function. We call this the principal square root, or the principal value of the square root function. For positive y, %:y is the positive square root, and there is also a negative root. For negative y, *:y is 0j1*%:-y. Its negation is also a square root of y. So to get both square roots of y, it suffices to write 1 _1*%:y. (For complex numbers in general other than real numbers, %:z is the square root with positive imaginary part. The other square root is the negation of this principal square root.) Inserting this into the solution above above gives:

    (r+1 _1*%:n)=x

It now remains to find n. We know that the square of x-r is (x^2)+(_2*r*x)+r^2. So if we have (x^2)+(_2*r*x)+k we want to add s=.(r^2)-k, or equivalently, in terms of values we know (*:b%2)-k, which results in the square of the monomial. We then know how to take the square roots of (*:b%2)-k, and it is easy to arrange the results in a convenient form for calculation. This method of solution is called completing the square.

Given a polynomial in the form (a*x^2)+(b*x)+c, we therefore proceed as follows:

   0=(a*x^2)+(b*x)+c
   0=(x^2)+((b%a)*x)+c%a
   (-c%a)=(x^2)+(b%a)*x
   ((-c%a)+(*:b%2*a))=(x^2)+((b%a)*x)+*:(b%2*a)
   ((-c%a)+((*:b)%4**:a))=(x^2)+((b%a)*x)+*:(b%(2*a))
   (((-4*a*c)%4**:a)+(*:b)%4**:a)=(x^2)+((b%a)*x)+*:(b%(2*a))
   (((-4*a*c)+*:b)%4**:a)=(x^2)+((b%a)*x)+*:(b%(2*a))
   (((*:b)-4*a*c)%4**:a)=(x^2)+((b%a)*x)+*:(b%(2*a))
   *./((1 _1*%:(*:b)-4*a*c)%2*a)=x+b%2*a
   *./((-b%2*a)+(1 _1*%:(*:b)-4*a*c)%2*a)=x
   *./(((-b)+(1 _1*%:(*:b)-4*a*c))%2*a)=x

For example:

   ('c';'b';'a')=.1 2 1
   1 2 1 p.i.5
1 4 9 16 25
   ]x=.((-b)+(1 _1*%:(*:b)-4*a*c))%2*a
_1 _1
   (a*x^2)+(b*x)+c
0 0
   (c,b,a) p. x
0 0
   ('c';'b';'a')=._1 _1 1
   _1 _1 1 p._4+i.10
19 11 5 1 _1 _1 1 5 11
   ]x=.((-b)+(1 _1*%:(*:b)-4*a*c))%2*a
1.61803 _0.618034
   (c,b,a) p. x
2.22045e_16 2.22045e_16

This result is correct within the limit of precision of computer arithmetic. Examining the formula for the solutions shows that these roots are of the form:

   -b%a
1
   (*:b)-4*a*c)
5
   (1+1 _1*%:5)%2
1.61803 _0.618034

The first of these is known as the Golden Ratio. It plays an important role in art, mathematics, biology and other areas.

We can now use the expression derived above to create an explicit function definition for solving quadratic equations given the coefficient vector c, as follows:

   qe=.3 : 0
('c';'b';'a')=.y
((-b)+(1 _1*%:(*:b)-4*a*c))%2*a
)

For example:

   c=._12 1 1
   qe c
3 _4
   c p. 3 _4
0 0

For an equation with real coefficients, the character of the solution depends on the sign of the expression (*:b)-4*a*c, the argument of the square root function. If it is positive, there are two real roots. If it is 0, there are two equal real roots. If it is negative, there are two complex roots. This expression is therefore called the discriminant, because it discriminates among these cases. For example, if we first define disc as the discriminant function:

   disc=.3 : 0
('c';'b';'a')=.y
(*:b)-4*a*c
)
   c=.1 2 1    NB. 1+(2*x)+x^2
   disc c
0
   qe c        NB. Equal roots
_1 _1
   c=.1 0 1    NB. 1+x^2
   disc c
_4
   qe c        NB. Complex roots
0j1 0j_1
   c=._1 0 1   NB. _1+x^2
   disc c
4
   qe c        NB. Distinct real roots
1 _1

Cubic Equations

The cube of a positive real number is positive, and the cube of a negative real number is negative. It follows that any cubic polynomial takes on both positive and values for real arguments of sufficient magnitude, where the cubic term is much greater in magnitude than the other terms. The polynomial thus must be zero for some real argument, having at least one real root. It may have three, as shown in the following graph.

Graph

The simplest cubic equation is 1=x^3. It has the three roots

1 _0.5j0.866025 _0.5j_0.866025

If we set w1=._0.5j0.866025 and w2=._0.5j_0.866025, we get the following relationships:

   w1=w2^2
   w2=w1^2
   1=w1*w2
   _1=w1+w2

Thus w1 and w2 are the roots of the quadratic equation 0=_1 0 1 p.x . Using the quadratic formula above we get the expression 0.5*_1+1 _1*%:_3, which yields these complex cube roots of 1.

There is no method for completing a cube, that is, for rearranging a cubic into the form (x-r)^3. What we have to do instead is to find a quadratic relationship among its coefficients (and thus among its roots) that we can solve and substitute back into the cubic. This method was worked out by Scipione del Ferro and Tartaglia, and published by Gerolamo Cardano in 1545. The full solution for negative and complex roots was not recognized until much later. However, the appearance of complex numbers in expressions for real roots was essential to the development of complex arithmetic, raising questions in the 16th century that were not fully resolved until the 19th.

To solve a cubic equation 0=cv p. x, first we want to simplify the form of the coefficient vector cv=.d,c,b,a. As before, we divide through by the last element of cv to produce a monic polynomial. Then we make the substitution t-(2{cv)%3 for x in this monic polynomial, which has the result that the new coefficient vector has the form d,c,0 1. This form is called a depressed cubic. Note that:

   (t-(2{cv)%3)=x
   t=x+(2{cv)%3

For example, starting from the monic coefficient vector cv, let us define:

   ('d';'c';'b';'a')=.cv
   p=.((3*c)-(b^2))%3
   q=.((2*b^3)+(_9*b*c)+27*d)%27

Then we can make the substitution:

   0=c p. x
   0=c p. t-b%3
   0=d+(c*t-b%3)+(b*(t-b%3)^2)+d*(t-b%3)^3
   0=q+(p*t)+t^3
   0=(q,p,0 1) p. t

A significant amount of routine algebraic manipulation has been omitted above. It can be treated as an exercise, or the student can verify the equivalence numerically by assigning appropriate values to the variables and executing each expression. For example:

   cv=._12 5 6 1
   b=.2{cv
6
   c=.1{cv
5
   d=.0{cv
_12
   p=.((3*c)-(b^2))%3
_7
   q=.((2*b^3)+(_9*b*c)+27*d)%27
_6

This yields the depressed cubic with coefficients cvd=._6 _7 0 1. We can verify this as follows:

   ]x=._5+i.11
_5 _4 _3 _2 _1 0 1 2 3 4 5
   ]t=.x+(2{cv)%3
_3 _2 _1 0 1 2 3 4 5 6 7
   poly=._12 5 6 1&p.
   dpoly=._6 _7 0 1&p.
   poly x
_12 0 0 _6 _12 _12 0 30 84 168 288
   dpoly t
_12 0 0 _6 _12 _12 0 30 84 168 288

Note that in this example, both forms have three real, integer roots. This will not be the case in general.

If we now introduce two new variables linked by the conditions (u+v)=t and 0=p+3*u*v, and substitute u+v for t in the depressed cubic, we get

   0=((u+v)^3)+(p*(u+v))+q
   0=(u^3)+(3*(u^2)*v)+(3*u*v^2)+v^3+(p*(u+v))+q
   0=(u^3)+(v^3)+(3*u*v*(u+v))+(p*(u+v))+q
   0=(u^3)+(v^3)+((p+3*u*v)*(u+v))+q
   0=(u^3)+(v^3)+q
   (-q)=(u^3)+(v^3)

 We also know that

   0=p+3*u*v
   (-p)=3*u*v
   (-p%3)=u*v
   (-(p^3)%27)=(u^3)*(v^3)

Knowing the sum and product of u^3 and v^3, we can construct a quadratic equation with u^3 and v^3 as roots, namely

   0=((-(p^3)%27),q,1) p. x

The solutions for u^3 and v^3, using the quadratic formula above, are

   (-q%2)+1 _1*%:((q^2)%4)+(p^3)%27)
In the case where the discriminant of the depressed cubic
((q^2)%4)+(p^3)%27
  is positive or 0, we can take its real principal square root and get
   u0=.3%:(-q%2)+%:((q^2)%4)+(p^3)%27
   v0=.3%:(-q%2)-%:((q^2)%4)+(p^3)%27
   x0=.u0+v0

The other two solutions of the cubic equation require finding the other two cube roots of u^3 and v^3, which turn out to be w1 and w2 times the primary cube roots. We must be careful to match the resulting four values in pairs, observing the condition 0=p+3*u*v. The result is

   x1=.(w1*u)+(w2*v)
   x2=.(w2*u)+(w1*v)

These pairings work because 1=w1*w2, so (u*v)=(w1*u)*(w2*v) and (u*v)=(w2*u)*(w1*v).

However, if the discriminant is negative, using the cube root function 3&^: does not in general match the appropriate cube roots. We must use the relationship

   p+3*u*v=0
   3*u*v=(-p)
   v=(-p)%(3*u)

There is one other case to consider, if u is 0, in which case p is 0. Then we cannot use the relationship above, but must use

   ((u^3)+(v^3))=-q
   (v^3)=-q
   v=3%:-q

Continuing the previous example, we have

   cv=._12 5 6 1
   poly=.cv&p.   NB. In x
   cvd=._6 _7 0 1
   dpoly=.cvd&p.  NB. In t=.x+(2{cv)%3
   a=.1
   b=.0
   c=.1{cvd
   d=.0{cvd
   ]p=.((3*c)-(b^2))%3
_7
   ]q=.((2*b^3)+(_9*b*c)+27*d)%27
_6
   ]u=.3%:(-q%2)+%:((q^2)%4)+(p^3)%27
1.5j0.288675
   ]v=.(-p)%(3*u)
1.5j_0.288675
   ]t0=.u+v
3
   ]t1=.(w1*u)+(w2*v)
_1
   ]t2=.(w2*u)+(w1*v)
_2
   dpoly t=.t0,t1,t2  NB. 3 _1 _2
0 1.77636e_15 _2.66454e_15
   x=.(t-(2{cv)%3)
1 _3 _4
   poly x
0 0 0

These results are thus excellent approximations of the exact roots, with errors less than 1e_14.

Using the method described above, we can write a function to solve cubic equations as follows:

   cubic=.3 : 0
cvm=.y%3{y
('d';'c';'b';'a')=.cvm
p=.((3*c)-(b^2))%3
q=.((2*b^3)+(_9*b*c)+27*d)%27
NB. cvd=.q,p,0 1
u=.3%:(-q%2)+%:((q^2)%4)+(p^3)%27
if. u=0
do. v=3:-q
else. v=.(-p)%(3*u)
end.
w1=._0.5 j. _0.5*%:3
w2=.%w1
(u+v),((w1*u)+(w2*v)),((w2*u)+(w1*v))
)

For example:

   cv=._1 3 _3 1
   ]x=.cubic cv
1 1 1
   cv p. x
0 0 0
   cv=.1 2 3 4
   ]x=.cubic cv
_0.60583 _0.0720852j0.638327 _0.0720852j_0.638327
   cv p. x
_2.22045e_16 _4.44089e_16j_1.80411e_16 0j_1.80411e_16

Again, these are excellent approximations of the roots, with errors less than 1e_15.

Quartic Equations

Gerolamo Cardano's student Lodovico Ferrari discovered a method for solving quartic equations by means of an auxiliary cubic equation in 1540, before Cardano published the solution of the cubic. Ferrari's method uses the concept of completing the square twice. The condition for completing the second square results in the auxiliary cubic. Several other methods have been discovered since using techniques such as factorization of the polynomial, Galois theory, and algebraic geometry.

The first two steps are as for the cubic equation, to reduce the quartic to monic form (divide by the coefficient of x^4) and then to depressed form (substitute t=.x-b%4, where b is the coefficient of x^3). The new coefficient vector will have the form cvd=.ed,dd,cd,0 1. In terms of the original coefficient vector cv=.e,d,c,b,a, these have the values

   ed=.((_3*b^4)%(256*a^4)+((c*b^2)%16*a^3)+((b*d)%4*a^2)+e%a
   dd=.((b^3)%8*a^3)+((b*c)%(2*a^2))+d%a
   cd=.(-3*(b^2)%8*a^2)+c%a

If dd=0, the polynomial can be rewritten as (ed,cd,1)p. x^2, that is, as an easily-solved quadratic in x^2. If ed=0, one of the roots is 0, and we can write an equation for the other roots as (dd,cd,0 1) p. x, that is, as a depressed cubic solvable by methods discussed above.

At this point we are given that 0=(ed,dd,cd,0 1) p.x for certain unknown values of x, or that

   (ed+(dd*x)+(cd*x^2)+x^4)=0

It is a standard identity that:

   (((x^2)+cd)^2)=(x^4)+(2*cd*x^2)+cd^2
   ((((x^2)+cd)^2)+(-x^4)+(-2*cd*x^2))=cd^2

Adding this equation to the previous one yields:

   ((ed+(dd*x)+(cd*x^2)+x^4)+(((x^2)+cd)^2)+(-x^4)+(-2*cd*x^2))=cd^2
   (((((x^2)+cd)^2)+((x^4)+-x^4)+(dd*x)+(cd*x^2)-2*cd*x^2)+ed)=cd^2
   ((((x^2)+cd)^2)+(dd*x)+(-cd*x^2)+ed)=cd^2
   ((((x^2)+cd)^2)+(dd*x)+ed)=(cd*x^2)+cd^2
   (((x^2)+cd)^2)=(cd*x^2)+(cd^2)+(-dd*x)-ed

The left side of this equation is now a perfect square. Now we can add another identity to introduce a variable y, while maintaining the lefthand term as a perfect square. We derive it as follows:

   (((x^2)+cd+y)^2)=(((x^2)+cd)^2)+(2*((x^2)+cd)*y)+y^2
   ((((x^2)+cd+y)^2)-((x^2)+cd)^2)=(2*((x^2)+cd)*y)+y^2
   ((((x^2)+cd+y)^2)-(((x^2)+cd)^2))=(2*(x^2)*y)+(2*cd*y)+y^2

Adding the identity

   0=((cd+2*y)*(x^2))+(-2*y*x^2)+(-cd*x^2)

to this yields:

   (((x^2)+cd+y)^2)-((x^2)+cd)^2
   (2*(x^2)*y)+(2*cd*y)+(y^2)+((cd+2*y)*(x^2))+(_2*y*x^2)+(-cd*x^2)
   (2*(x^2)*y)+(2*cd*y)+(y^2)

Adding this result to the earlier equation:

   (((x^2)+cd)^2)=(cd*x^2)+(cd^2)+(-dd*x)-ed

now yields:

   (((x^2)+cd)^2)+(((x^2)+cd+y)^2)-((x^2)+cd)^2
   ((x^2)+cd+y)^2
   ((2*(x^2)*y)+(2*cd*y)+(y^2))+(cd*x^2)+(cd^2)+(-dd*x)-ed

which is equivalent to:

   (((x^2)+cd+y)^2)=((cd+2*y)*(x^2))+(-dd*x)+(y^2)+(2*cd*y)+(cd^2)-ed

We would like to make the right-hand side a perfect square by choosing a suitable value of y. We know that a quadratic polynomial (c,b,a) p. x is a perfect square when the discriminant is 0. That is:

   0=(b^2)-4*a*c
   4*a*c=(b^2)
   c=(b^2)%4*a

The polynomial is thus

   (a*x^2)+(b*x)+(b^2)%4*a
   a*(x^2)+((b%a)*x)+(b^2)%⁴*a^2
   a*(x+b%2*a)^2

Setting

   aq=.cd+2*y
   bq=.-dd
   cq=.(y^2)+(2*cd*y)+(cd^2)-ed

gives us

   0
   (bq^2)-4*aq*cq
   ((bq^2)%_4)+aq*cq
   ((dd^2)%_4)+(cd+2*y)*((y^2)+(2*cd*y)+(cd^2)-ed)
   ((dd^2)%_4)+(cd*(y^2)+(2*cd*y)+(cd^2)-ed)+2*y*((y^2)+(2*cd*y)+(cd^2)-ed)
   ((dd^2)%_4)+((cd*y^2)+(2*(cd^2)*y)+(cd^3)-cd*ed)+(2*y^3)+(4*cd*y^2)+(2*(cd^2)*y)-2*ed*y
   (2*y^3)+((5*cd)*y^2)+(((4*cd^2)-2*ed)*y)+((cd^3)+(-(dd^2)%4)-cd*ed)

 This is a cubic equation in y, which we know how to solve. Set

A=.2
B=.5*cd
C=.(4*cd^2)-2*ed
D=.(cd^3)+(-(dd^2)%4)-cd*ed
r=.cubic D,C,B,A

Any of the three solutions will serve our purpose. Each solution of the cubic allows us to complete the second square, resulting in a pair of quadratic equations. Each such pair matches the roots of the quartic in pairs, in one of three possible arrangements:

r0 r1  r2 r3
r0 r2  r1 r3
r0 r3  r1 r2

Substituting 0{r for y, our equation now takes the form:

   (((s^2)*u^2)+((2*s*t)*u)+t^2)
   (s^2)*(u+2*s*t)%2*s^2))^2
   (s*(u+2*s*t)%2*s^2)^2

   (((x^2)+cd+y)^2)=((cd+2*y)*(x^2))+(-dd*x)+(y^2)+(2*cd*y)+(cd^2)-ed
   (((x^2)+cd+y)^2)=(((%:cd+2*y)*x)-bd%(2*%:cd+2*y))^2

Then we can take square roots and collect like powers of x:

   ((x^2)+cd+y)=((%:cd+2*y)*x)-bd%(2*%:cd+2*y)
   0=(x^2)+(-(%:cd+2*y)*x)+((-(cd+y))+(-bd%(2*%:cd+2*y)))

Or, with the opposite sign for the square root:

   ((x^2)+cd+y)=-((%:cd+2*y)*x)-bd%(2*%:cd+2*y)
   ((x^2)+cd+y)=(-(%:cd+2*y)*x)+bd%(2*%:cd+2*y)
   0=(x^2)+((%:cd+2*y)*x)+(-(cd+y))+bd%(2*%:cd+2*y)

Now we have two quadratics in x, with coefficients

   aq=.1
   bq=._1 1*(%:cd+2*y)
   cq=.-((cd+y)+1 _1*bd%(2*%:cd+2*y))

where the signs of b and c must be correctly matched, as above. The four solutions of these two quadratics are the roots of the depressed quartic, as follows:

   ((-bq)+,1 _1*/%:(bq^2)-4*aq*cq))%2*aq
   ((-bq)+,1 _1*/%:(bq^2)-4*cq))%2

We transform these results to roots of the original quartic by adding -b%4.

Collecting just the expressions that we need to evaluate out of the derivation above, and from the solutions to cubic and quadratic equations, we can define a quartic equation solver as follows:

   quartic=:3 : 0
('e';'d';'c';'b';'a')=:y
NB. Depressed cubic
alpha=:(_3r8*(b^2)%a^2)+c%a
beta=:(1r8*(b^3)%a^3)+(-(b*c)%2*a^2)+d%a
gamma=:((_3r256*b^4)%a^4)+((c*b^2)%16*a^3)+(-(b*d)%4*a^2)+e%a
NB. Biquadratic special case
if. 0=beta
do. roots0=:(-b%4*a)+,1 _1*/%:((-alpha)+1 _1*%:(alpha^2)-4*gamma)%2
else.
NB. Solve auxiliary cubic
P=:((-alpha^2)%12)-gamma
Q=:((-alpha^3)%108)+(alpha*gamma%3)-(beta^2)%8
R=:(-Q%2)+%:((Q^2)%4)+(P^3)%27
U=:3%:R
if. 0=U
do. r3=:(_5r6*alpha)-3%:Q
else. r3=:(_5r6*alpha)+U-P%3*U
end.
NB. The second square is completed; use quadratic formula
w=:%:alpha+2*r3
v=:%:-((3*alpha)+(2*r3)+1 _1*2*beta%w)
(-b%4)+((1 _1 1 _1*w)+,1 _1*/v)%2
end.
)

For example:

   quartic 1 4 6 4 1   NB. 0=(x+1)^4
_1 _1 _1 _1
   1 4 6 4 1 p. _1
0
   quartic _1 0 0 0 1  NB. 1=x^4
1 0j1 _1 0j_1
   _1 0 0 0 1 p. 1 0j1 _1 0j_1
0 0 0 0
   quartic 30 _61 41 _11 1
5 2 3 1
   30 _61 41 _11 1 p. i.6
30 0 0 0 _6 0

Equations of Higher Degree

The methods described above fail to work for quintic (fifth-degree) equations and for all higher degrees. When one tries Ferrari's method on quintic polynomials, one does indeed get an auxiliary equation, but it is of the sixth degree, and so is of no help. The Abel-Ruffini theorem established that equations of higher degree cannot be solved in radicals, and Galois theory greatly clarified the relationship between the coefficients and the roots of any polynomial. However, these topics are far beyond the scope of this textbook.


EDIT

There has been error in communication with booki server. Not sure right now where is the problem.

You should refresh this page.