• Back to MAIN

  • In the year 1834 Nikolaj Ivanovich Lobachevsky (1.Dec. 1792 - 24.Febr.1856) found a way to approximate the
    roots of algebraic equations. This method has been discovered independently in the year 1837 by the German
    mathematician Karl Heinrich Graeffe (7.Nov.1799 - 2.Dec.1873) and the French mathematician Germinal
    Pierre Dandelin (12.Apr.1794 - 15.Febr.1847)
    . We can see it today as a very convenient numerical method.

    Before we start here, we should be aware of the fact that the equation:P(x) = xn+a1xn-1+a2xn-2+ ... +an=0,
    with the roots r1,r2, ... ,rn and a given symetrical polynomial F[x1,x2,...,xn], we can determine the number
    s = F[r1,r2,...,rn] without solving our original equation P(x) = 0.

    As a consequence - in case of a given equation: P1(x)=a0xn+a1xn-1+a2xn-2+ ... +an=0, and its roots
    r1,r2, ... ,rn we can find(without solving P1) an other polynomial Pk(x)=ak0(x-rk1)(x-rk2)...(x-rkn), eg which has
    roots which are exactly the squares, cubes etc. of the roots of P1(x).

    Just as an example, if the roots of the (above given)polynomial P(x) are x1,x2,x3,x4,x5 then we can express the
    sum of their 5th powers as:(under the given index set {1,2,3,4,5})

    S[x5i] := x51+x52+x53+x54+x55 ; following this notation, we can now deduct:

    S[x4i]S[xi] = S[x5i] + S[[x4ixj];i < j] ==> S[x5i] = S[x4i]S[xi]-S[x4ixj;i < j]= -a1S[x4i]-S[x4ixj;i < j], but
    S[x3i]S[xixj;i < j] = S[x4ixj;i < j] + S[x3ixjxk;i < j < k] ==> S[x4ixj;i < j] = a2S[x3i] - S[x3ixjxk;i < j < k], but
    S[x2i]S[xixjxk;i < j < k] = S[x3ixjxk;i < j < k] + S[x2ixjxkxl;i < j < k < l] .. which means that
    S[x3ixjxk;i < j < k] = -a3S[x2i] - S[x2ixjxkxl;i < j < k < l] .. and finally
    S[xi]S[xixjxkxk;i < j < k < l] = S[x2ixjxkxl;i < j < j < l] + S[xixjxkxlxm;i < j < k < l < m] ..and so
    S[x2xixjxkxl;i < j < j < l] = a4S[xi] - 5a5 summarized as :

    S[x5i] = -a1S[x4i]-a2S[x3i]-a3S[x2i]-a4S[xi]-5a5 and by setting qk:=ak1+ak2+...+akn we can reformulate it
    as a general statement:

    If the ai are elementary symetrical functions of the variables x1,x2,...,xn, then for all k < n+1 holds:
    qk+a1qk-1+...+ak-1q1+kak = 0

    ..and in the same way we would be able to show that if k > n

    qk+a1qk-1+a2qk-2...+anqk-n= 0

    These both equations we know as Newton's recursion formula's . In the case n=4 we would obtain:

    q1+a1=0, q2+a1q1+2a2=0, q3+a1q2+a2q1+3a3=0, q4+a1q3+a2q2+a3q1+4a4=0, q5+a1q4+a2q3+a3q2+a4q1=0
    q6+a1q5+a2q4+a3q3+a4q2=0, q7+a1q6+a2q5+a3q4+a4q3=0, e.t.c. and consequently..

    q1=-a1, q2= a21-2a2, q3=-a31+3a1a2-3a3 q4= a41-4a21a2+4a1a3+2a22-4a4
    q5=-a51+5a31a2-5a21a3-5a1a22+5a1a4+5a2a3 e.t.c.

    A small sparkling spin-off from all of this:

    if we have P1(x)=a0xn+a1xn-1+a2xn-2+ ... +an=0, and our roots r1,r2,...,rn, then r21+r22+...+r2n=a21-2a2
    which means that if a21 <= 2a2, then at least one pair of roots must be complex.(we assume also an not 0,
    and n > 2)

    One important remaining point here:

    Lets once look again at our P1(x)=a0xn+a1xn-1+a2xn-2+ ... +an=0, with the roots r1,r2,...,rn, and that means :
    P1(x)=a0(x-r1)(x-r2) ... (x-rn), and therefore P1(-x)=(-1)na0(x+r1)(x+r2) ... (x+rn). What we look for is
    a polynomial Q2(x) having r21,r22,...,r2n as roots - without solving P1(x)=0, of course.
    But (-1)nP1(x)P1(-x)=a20(x2-r21)(x2-r22) ... (x2-r2n)=Q2(z) with z=x2, and by decomposing P1(x)=H1(x)+xH2(x),
    where H1(x) consists of even exponents of x, and xH2(x) of the odd ones, we can conclude that
    P1(-x)=H1(x)-xH2(x), eg:


    But now we are better prepared for our little rendez-vous with an idea which now celebrates its 170th birthday.

    As a nature of the thing, we can ask us the question, if the polynomials Q2,Q4,Q8,Q16..(having the 2nth power
    of the original polynomial as roots) can help us to find r1,r2,...,rn in a numerical way. The answer to that
    question is yes. Lets assume, we know all ri and lets for a moment also assume they are all real and different,
    so after re-indexing we have: |r1| > |r2| > |r3| > ... > |rn-1| > |rn| .. so if we take Q2k = a0,kxn+a1,kxn-1+...+an,k=0

    and the known relations we see that:

    -a1,k/a0,k = r2k1+r2k2+ ... +r2kn = r2k1[1+(r2k2/r2k1)+...+(r2kn/r2k1)]
    +a2,k/a0,k = r2k1r2k2+r2k1r2k3 + ... + r2kn-1r2kn = r2k1r2k2(1+r2k3/r2k2 + ... + r2kn-1r2kn/r2k1r2k2)
    -a3,k/a0,k = r2k1r2k2r2k3+r2k1r2k2r2k4 + ... + r2kn-2r2kn-1r2kn = r2k1r2k2r2k3(1+(r2k4/r2k3) + ... + (r2kn-2r2kn-1r2kn/ r2k1r2k2r2k3))


    (-1)nan,k/a0,k = r2k1r2k2 ... r2kn-1r2kn

    .. with other words...

    -a2,k/(a1,kr2k2)=(1+r2k3/r2k2 + ... + r2kn-1r2kn/r2k1r2k2) / 1+(r2k2/r2k1)+...+(r2kn/r2k1)
    -a3,k/(a2,kr2k3)=(1+(r2k4/r2k3) + ... + (r2kn-2r2kn-1r2kn/ r2k1r2k2r2k3))/1+r2k3/r2k2 + ... + r2kn-1r2kn/r2k1r2k2



    But as |r2ki+1/r2ki| < 1, all -ai+1,k/(ai,kr2ki+1) must converge to 1 if i --> infinity. That means for "larger" values k,
    we have

    r2ki+1 ~ -ai+1,k/ai,k

    It is absolutely clear, that we could have taken any exponent - not just only powers of 2, we took it for practical

    One possible solution using k4 can be:

    /some preliminaries

    mod:{x-y xbar x}
    dix:{x@&~(!#x) in y}

    xexp:{exp y*log x}

    q:{+/rot'[-!#x;(,/'(s#'x) *\:+,x),\:(s+1-2*s:#x)#0]}

    /lobachevsky .. just the correct sign needs to be verified after run .. a trivial K exercise
    lb:{xexp[(abs@ - %/(1,-1) _\: x {f@ (q@ 1f*x*~b) - q@ 1f*x*b:mod . ((!#x);2)}/1f*y);%xe2 x]}

    Sample 1: x5-2x4-19x3+13x2+36x-5 .. lb . (4;1 -2 -19 13 36 -5) returns:
    x1=5.001294, x2=3.548031, x3= 1.675143, x4= 1.258001, x5=0.133711
    ..the exakt solutions are (6 digits)
    x1=5.000000, x2= -3.548948, x3=1.674056, x4= -1.258819, x5=0.133711

    ..we see here that the polynomials having the 16th power as roots of the original polynomial(eg.. we handle
    here polynomials of degree 80) bring us actually very close to the real world..as we said, just ignoring the sign.
    The precision (although very convincing in a normal case) is not the real strength of this method. The advantage
    of the brilliant Lobachevsky approach is clearly, that we get all solutions in one go.

    By definition, the assumption of strictly decreasing absolute values of the roots excludes the search for complex

    If we are interested in the search for the complex roots - we slightly change our approach .. Just as an
    illustrative display for a cubic polynomial(..and real coefficients), having one pair of complex solutions:

    We assume that, r1=Reit and r2=Re-it and |r1|=|r2|>|r3|,
    so: rk1+rk2+rk3=Rk(eikt+e-ikt)+rk3=Rk(2cos kt + rk3R-k)=-a1,k/a0,k
    and -a2,k/a0,k=rk1rk2+rk1rk3+rk2rk3=R2k[1+2(r3/R)kcos kt] and -a3,k/a0,k = rk1rk2rk3 = R2krk3

    (.. in case of |r1|=|r2|<|r3| we would take -a1,k/a0,k=rk3[1+2(R/r3)kcos kt]
    and -a2,k/a0,k=rk3Rk(2cos kt+(R/r3)k] and -a3,k/a0,k = R2krk3

    So..-a3,k/a2,k ~ rk3 as (r3/R)k goes to 0 if k goes to infinity, and -a1,k/a0,k ~ rk3, and -a2,k/a0,k~R2k

    Sample 2: P(x) = x3+4x2-2x-20=0, by taking Q16(z)=z3-8.446003*107z2+1.000553*1016z-6.5536*1020
    , which means R~3.162332 and r3~1.999931,(the correct figures are r3=2 and R=3.162278)
    Based on the knowledge of r3=1.999931 and that r1+r2+r3=-4 follows...that cos t = 0.5R-1(-5.999931)
    So, we obtain (as sign is not determined) t = 2.819756 resp. t = 0.3218366 ... a verification would tell us,
    that the first solution is the correct one., So r1=-3+i and r2=-3-i.

    One friendly hint here: in case of cubic polynomial with only real solutions, all polynomials Q2k(x)
    must have positive roots, which means their coefficients must have the sign pattern ( + - + - ), which
    is not the case in our example, consequently our cubic polynomial must have complex solutions.
    It is by definition clear, that polynomials with real solutions having the same absolute value cannot be
    handled by this method but those are just exceptional cases.

  • Back to MAIN