🚀 go-pugleaf

Root Finder viii.
by Jon Giffen

The polynomial

a[0]+a[1]t+a[2]t^2+a[3]t^3+ ... + a[n] = 0

may be expressed as the dot product of the two vectors,

T = (t,t^2,t^3,...,t^n) and
N = (a[1],a[2],a[3],...,a[n])   or,

(t,t^2,t^3,..,t^n)*(a[1],a[2],a[3],..,a[n])+a[0]=0  or

T*N+a[0]=0  is the nth degree polynomial

N may be considered the normal to the more general case
of the plane,

a[0] + a[1]x[1] + a[2]x[2] + ... + a[n]x[n] = 0

where x[1]=t  x[2]=t^2 x[3]=t^3 ... x[n]=t^n


Following T up from the origin, R is the variable vector
parallel to N and orthogonal to T-R, as T varies along
its curve. then

R*(T-R)=0

taking the derivative with respect to t,

R'*(T-R)+R*(T'-R')=0

or

2R'*R = R'*T + R*T'

The projection of T' onto N/|N| is R'

T'*N  N
---- ---- = R'
|N|  |N|

define

Q = (T*N)N/|N|^2 = -a[0]N/|N|^2
Q = (-a[0]/|N|^2)(a[1],a[2].a[3],..a[n])
or
Q = (q[1],q[2],q[3],...,q[4])


Since R' is in the same direction as Q, it holds that
since Q and T'-(T'*N)N/|N|^2 are orthogonal,

Q*(T'-(T'*N)N/|N|^2) = 0

where Q is the shortest vector from the origin to the plane.

T = (t,t^2,t^3,...,t^n)
T'= (1,2t,3t^2,..,nt^(n-1) )

T'*Q =q[1]+2q[2]t+3q[3]t^2+....+ nq[n]t^(n-1)

      =(1/t)T*(q[1],2q[2],3q[3],...,nq[n])

T'*N = a[1]+2a[2]t+3a[3]t^2+....+ na[n]t^(n-1)

      =(1/t)T*(a[1],2a[2],3a[3],...,na[n])

Suppose the vectors D,S,U are defined,

D=(  1 ,   2 ,   3 ,...,n    )
S=(a[1],2a[2],3a[3],...,na[n])
U=(  1 ,   1 ,   1 ,...,1    )

Then

T'* U =(1/t)T*D
T'* N =(1/t)T*S

for instance,
  (1,2t,3t^2,..,nt^(n-1) )*(1,1,1,..)
=(1/t)(t,t^2,t^3,...,t^n)*(1,2,3,..)
and

  (1,2t,3t^2,..,nt^(n-1) )*(a[1],a[2],a[3],..)
=(1/t)(t,t^2,t^3,...,t^n)*(a[1],2a[2],3a[3],..)

dividing the below equations with each other,

T'* U =(1/t)T*D
T'* N =(1/t)T*S

(T'*U)(T*S) - (T'*N)(T*D) = 0

factoring out T,

T*((T'*U)S - (T'*N)D) = 0

so T is orthogonal with (T'*U)S - (T'*N)D
or T is orthogonal with (1/t)((T*D)S - (T*S)D)
or T is orthogonal with (T*D)S - (T*S)D

this is true, since
T*((T*D)S - (T*S)D) = 0

Let,

G = (T*D)S - (T*S)D

If m is some ratio of G such that,

(mG-Q)*N = 0

mG*N - Q*N = 0    since Q*N=T*N,

m[(T*D)(S*N)-(T*S)(D*N)]-(T*N) = 0

or

T*[m((S*N)D-(D*N)S )-N] = 0

One solution for m, among others, can be found from

m((S*N)D-(D*N)S - N = 0 or

m((S*N)D-(D*N)S = N   taking the dot product of both
sides with D and solving for m,

             N*D
m = ---------------------
     (S*N)(D*D)-(D*N)(S*D)


Otherwise, T and m((S*N)D-(D*N)S )-N are orthogonal
Notice that since N*((S*N)D-(D*N)S )= 0,
((S*N)D-(D*N)S ) is orthogonal to N.

T is orthogonal with  m((S*N)D-(D*N)S)-N
N is orthogonal with (S*N)D-(D*N)S


OY=m((S*N)D-(D*N)S)-N

Y
*------------ N              T
  \            |            /
    \          |          /
      \        |        /
        \      |      /
          \    |    /
            \  |  /
   Z _________\|/
               |O
OZ=           |
(S*N)D-(D*N)S |
               |
               |
               |
              -N


T is orthogonal with OY
N is orthogonal with OZ


The plane formed by OY and OT is not necessarily
the plane formed by OY and OT, but

OT*OY=0

Suppose OP is orthogonal to ON and OZ then

OP*ON=0  and
OP*OZ=0
OP*OY=0
OP*OT ?=0

OT*OY=0
or
OT*(mOZ-ON)=0

N = (Q/|Q|)|N|

T*(m(|Q|/|N|)OZ - Q)=0


Q            T
*----------*
|         /
|       /
|     /
|   /
| /
*
O

angle(OQT)=pi/2

QT=m(|Q|/|N|)OZ


Hence, OP*OT=0

since angle YOZ = angle TON,

   OY*OZ     T*N            OY*N
-------- = ------ = {1 - (-------)^2 }^(1/2)
|OY||OZ|   |T||N|         |OY||N|

  (OY*OZ)^2            (OY*N)^2
------------  = 1 -  -----------
|OY|^2|OZ|^2         |OY|^2|N|^2

Multiplying both sides of this equation by
|OY|^2|OZ|^2|N|^2,

(OY*OZ)^2|N|^2 = |OY|^2|OZ|^2|N|^2 - (OY*N)^2|OZ|^2  eqn i.

Letting (S*N)D-(D*N)S = C, OZ=C and OY=mC-N and

(OY*OZ)^2=((mC-N)*C)^2|N|^2 = m^2|C|^2|N|^2
|OY|^2=(mC-N)*(mC-N) = m^2|C|^2+|N|^2
|OZ|^2=|C|^2
(OY*N)^2=((mC-N)*N)^2=-|N|^4

Substituting these in eqn i,

(m^2|C|^2|N|^2)|N|^2
=
(m^2|C|^2+|N|^2)|C|^2|N|^2
+
|N|^4|C|^2  or

m^2|N|^2 = m^2(|C|^2+|N|^2) + N^2


N^2 = m^2(|N|^2 - |C|^2 - |N|^2)

        |N|
m = - -----
        |C|


So

OY=m((S*N)D-(D*N)S)-N  or

       |N|
OY= - ---((S*N)D-(D*N)S)-N
       |C|

or
            |N|
OY= - ---------------((S*N)D-(D*N)S) - N
       |(S*N)D-(D*N)S|

or
       |N|
OY= - --- C - N
       |C|


(qOY-Q)*N=0
qOY*N-Q*N=0
q(0-|N|^2) = Q*N
q(-|N|^2)  = -a[0]
q = a[0]/|N|^2


angle ZOY

           OY*OZ      T*N
cos(ZOY)=-------- = ------ = cos(NOT)
          |OY||OZ|   |T||N|

T*N= -a[0] and

          a[0]|OY||OZ|
|T| = - -------------
          |N|(OY*OZ)



          (a[0])^2|OY|^2|OZ|^2
|T|^2 =  ---------------------
             |N|^2(OY*OZ)^2

|OY|^2 = 2|N|^2
(OY*OZ)^2 = (|N||C|)^2

         (a[0])^2 2|N|^2 |C|^2
|T|^2 = --------------------- = 2|Q|^2
             |N|^2 |N|^2 |C|^2


T = Q + {|T|^2-|Q|^2}^(1/2)(-C/|C|)

            |Q|
T = Q +/- ----- C
            |C|

Where
C = (S*N)D-(D*N)S
D=(  1 ,   2 ,   3 ,...,n    )
S=(a[1],2a[2],3a[3],...,na[n])

Alternatively, suppose

T'*(-U)=(1/t)T*D
T'*  N =(1/t)T*S

Then dividing the two equations,

(T'*U)(T*S)+(T'*N)(T*D) =0

then

C = (S*N)D+(D*N)S


D=(  1 ,   2 ,   3 ,...,n    )
S=(a[1],2a[2],3a[3],...,na[n])
C=(S*N)D+(D*N)S


           |Q|
T = Q +/- --- C
           |C|


EXAMPLE:

t^3+2t^2+t-4=0  t=1
a[0]= -4
a[1]= 1 a[2]=2 a[3]=1 N=(1,2,1) |N|^2=6 Q=(2/3)(1,2,1)
|Q|=(2/3)6^(1/2)
D=(1,2,3) S=(1,4,3) S*N=(1,4,3)*(1,2,1)=1+8+3=12
D*N=(1,2,3)*(1,2,1)=1+4+3=8
C=12(1,2,3)+8(1,4,3)=(12+8,24+32,36+24)=(20,56,60)=4(5,14,15)
|C|=(4)446^(1/2)

                             (2/3)6^(1/2)
(t,t^2,t^3)=(2/3)(1,2,1)+/- ------------ (5,14,15)
                               446^(1/2)


(t,t^2,t^3)=(2/3)(1,2,1)+/- 0.07732(5,14,15)
(t,t^2,t^3)=(0.6666,1.3333,0.6666)+/-(0.3866,1.08254,1.1599)

Selecting the first component,

t  = 0.6666+0.3866 = 1.0532 ~ 1

then

       t^2 + 3t + 4
      ------------------------------
t-1/  t^3 + 2t^2 + t - 4
       t^3 -  t^2
       --------------------
             3t^2 + t - 4
             3t^2 -3t
             -------------
                   4t - 4

and the remaining roots are

     -3+/-{9-16|^(1/2)
t = ----------------- = -1.5+/-1.3229i
           2

Dividing t^3 by t^2,

t^3 = 0.6666+1.1599 = 2.2656
-----------------------------
t^2 = 1.3333+1.0825 = 2.4158

or

t = 0.9378 ~ 1

EXAMPLE

t^7+t-10=0  t=?
a[0]=-10 a[1]=1 a[7]=1 N=(1,0,0,0,0,0,1)
|N|^2=2  Q=5(1,0,0,0,0,0,1) |Q|=(5)2^(1/2)
D=(1,2,3,4,5,6,7)  S=(1,0,0,0,0,0,7)
D*N=1+7=8  S*N=1+7=8
C=8(1,2,3,4,5,6,7)+8(1,0,0,0,0,0,7)
  =8(2,2,3,4,5,6,14)
|C|=(8)290^(1/2)

(t,t^2,t^3,..,t^7)
=
                        5
5(1,0,0,0,0,0,1)+/- ----------(2,2,3,4,5,6,14)
                     145^(1/2)

t^7 = 5(1+14/145^(1/2))=10.813 ->t=1.405

t^7+t-10 = 0
(1.405)^7+1.405-10?=0
10.813+1.405-10 = 2.218

Use t=1.405 in Newton's Method

                  2.218
t[n+1]=1.405 - ------------ = 1.364
                7(1.405)^6+1

t^7+t-10
(1.364)^7+1.364-10 = 0.1734 ~ 0


end

Jon Giffen