Transcript Effective polynomial

Efficiant polynomial interpolation
algorithms
Federico Tessmann, Ferienakademie TUM
Overview

Introduction to Vandermonde Matrices and its
utilities

Univariate Interpolation

Multivariate Interpolation
Federico Tessmann, Ferienakademie TUM
Properties of Vandermonde Matrices

Easy to ensure that they are non-singular

Systems of linear equations whose
coefficients form Vandermonde matrices are
easy to solve exactly
Federico Tessmann, Ferienakademie TUM
The Vandermonde Matrix
1 k1

1 k 2
Vn  
 
1 k
n

 k 

 k 


 
n 1 
 kn 
Federico Tessmann, Ferienakademie TUM
n 1
1
n 1
2
Generalized Vandermonde
 k1e1
 e
1
k
 2
Vn  


 k e1
 n
where
e2
1
e2
2
k
k

e2
kn
 k 

 k 

  
en 
 kn 
e1  e2  en
Federico Tessmann, Ferienakademie TUM
en
1
en
2
Determinant of a Vandermonde
det Vn 
1 k1
1 k2
 k1n 1
 k 2n 1
 
1 kn


 k nn 1
1
0

0
1 k 2  k1



 k 2n 1  k 2n  2 k1


1 k n  k1
 k nn 1  k nn  2 k1
Federico Tessmann, Ferienakademie TUM
Determinant of a Vandermonde
1
0
1 k 2  k1


0
k 2n 1  k 2n  2 k1




k nn 1  k nn  2 k1

1 k n  k1
1 j  n
det Vn  det Vn1  (k j  k1 )  det Vn2
1 j n
det Vn 
 det Vn 1  (k j  k1 )
 (k
1i  j n
Federico Tessmann, Ferienakademie TUM
 (k
2 j n
j
j
 k2 )  (k j  k1 )
 ki )
1 j n
Determinant of a Vandermonde
det Vn 
 (k
j
 ki )
1i  j n
The Vandermonde matrix
is non-singular  the ki are distinct
Federico Tessmann, Ferienakademie TUM
The previous result can not be applyed for generalized Vandermonde
matrices

Example
e
1
e
2
1 k
det V2 
 k 2e  k1e
1 k
wich is 0 also when
 k2

 k1
e

  1

Federico Tessmann, Ferienakademie TUM
Non-singularity of generalized
Vandermonde matrices

Proposition 1:
If the ki are distinct positiv real numbers
=> the matrix is non-zero
Federico Tessmann, Ferienakademie TUM
The inverse of a Vandermonde matrix
1 k1

1 k 2
 

1 k
n

 k1n 1   p11
 
n 1
 k 2   p21
*


  

 k nn 1   pn1
p12
p22

pn 2
p1n   P1 (k1 ) P2 (k1 )
 
 p2 n   P1 (k 2 ) P2 (k 2 )


 


 
 pnn   P1 (k n ) P2 (k n )

 Pn (k1 ) 

 Pn (k 2 ) 

 

 Pn (k n ) 
Pj (ki )  p1 j  p2 j ki  p3 j ki    pnj ki
2
Federico Tessmann, Ferienakademie TUM
n 1
The inverse of a Vandermonde matrix
 P1 ( k1 )

 P1 ( k 2 )
En  


 P (k )
 1 n
P2 ( k1 )

P2 ( k 2 ) 


P2 ( k n ) 
Z  ki
Pj ( Z )  
i  j k j  ki
1i  n
Federico Tessmann, Ferienakademie TUM
Pn ( k1 ) 

Pn ( k 2 ) 
 

Pn ( k n ) 
Solving a Vandermonde system of
equations
x1  k1 x2  k1 x3    k1
2
n 1
x1  k 2 x2  k 2 x3    k 2
2
xn  w1
n 1
xn  w2





2
n 1
x1  k n x2  k n x3    k n xn  wn
Federico Tessmann, Ferienakademie TUM
Solving a Vandermonde system of
equations

T
x  ( x1 ,, xn )

T
w  (w1,, wn )
 
Vn  x  w
Vn
1

1 
Vn  x  Vn  w


1 
x  Vn  w
Federico Tessmann, Ferienakademie TUM
Solving a Vandermonde system of
equations
 x1   p11
  
 x2   p21
   
  
x   p
 n   n1
p12

p22



pn 2

p1n   w1 
  
p2 n   w1 
 



  
pnn   wn 
xi  pi1w1  pi 2 w2   pin wn
Pj (Z )  p1 j   pnj Z n1

xi  coef (P1 (Z ), Z i1 )w1   coef (Pn (Z ), Z i1 )wn
Federico Tessmann, Ferienakademie TUM
The algorithm to solve the system
Q( Z ) 
 (Z  k )
i
1i  n
Q( Z )
Qi ( Z ) 
Z  ki
Qi ( Z )
Pi ( Z ) 
Qi (ki )
Federico Tessmann, Ferienakademie TUM
The algorithm to solve the system
The computation of the xi is arranged as follows:
 wn  coef ( Pn , Z 0 ) 
 x1   w1  coef ( P1 , Z 0 ) 



  
1
1
 wn  coef ( Pn , Z ) 
 x2   w1  coef ( P1 , Z ) 
 


  





  
 x   w  coef ( P , Z n 1 ) 
 w  coef ( P , Z n 1 ) 
 n  1
n
1

 n

Calculate each vector and add it to the accumulating X
Federico Tessmann, Ferienakademie TUM
Analysis of the algorithm
1. By calculating the vectors one after the other we only
need to compute one Pi(Z) at the time
2. Each Pi(Z) only needs O(n) time and since we have n
polinoms to compute, the complexity is O(n2) and the
space needed is O(n)
3. Because the inverse of the transposed matrix is the
transpose of the inverse of the matrix, the algorithm
only need a little adjustment to solve a transposed
Vandermonde system of equations
On the Appendix there is an example of this alorithm taken from Zippel
Federico Tessmann, Ferienakademie TUM
Univariate Interpolation

Lagrange Interpolation

Newton Interpolation

Abstract Interpolation
Federico Tessmann, Ferienakademie TUM
Lagrange Interpolation
Giving are a set of distinct evaluation points with its
correspondating functional values
z0 , z1 ,, zn
P( zi )  wi
The goal is to find the polinome
P(Z )  p0  p1Z   pn Z n
Federico Tessmann, Ferienakademie TUM
Lagrange Interpolation
p0  p1 z0  p2 z0    pn z0  w0
2
n
p0  p1 z1  p2 z1    pn z1  w1
2
n





2
n
p0  p1 z n  p2 z n    pn z n  wn
This is a Vandermonde system where
pi  coef ( P0 (Z ), Z i )w0   coef ( Pn (Z ), Z i )wn
Federico Tessmann, Ferienakademie TUM
Lagrange Interpolation
P( Z )  coef ( P0 ( Z ), Z 0 ) w0    coef ( Pn ( Z ), Z 0 ) wn 
(coef ( P0 ( Z ), Z 1 ) w0    coef ( Pn ( Z ), Z 1 ) wn ) Z 

(coef ( P0 ( Z ), Z n ) w0    coef ( Pn ( Z ), Z n ) wn ) Z n

P(Z )  P0 (Z )w0   Pn (Z )wn
Federico Tessmann, Ferienakademie TUM
Lagrange Interpolation
P(Z )  P0 (Z )w0   Pn (Z )wn



Z  zi
P( Z )    
0 j  n  0i  n z j  zi
 i j
Federico Tessmann, Ferienakademie TUM


w j


Newton Interpolation
f(a)=f(x)(mod (x-a))
Federico Tessmann, Ferienakademie TUM
The Chinese remainder algorithm over Z
x  k1 (moda1 )
x  k2 (moda2 )

1
1
x  k1  (k2  k1 )  (a moda2 )  a1 (moda1  a2 )
a11 moda2  gcd(a1 , a2 )  1
Federico Tessmann, Ferienakademie TUM
Chinese remainder with polinoms
When given
z1 ,, zn
and
f zi   wi
Then we change it to the following situation:
Given
pi  x  zi
f ( X )  wi (mod pi )
Compute
f ( x)  an x   a0 (mod(p1  pn ))
n
Federico Tessmann, Ferienakademie TUM
Newton Interpolation algorithm
1.
Let f(x)=0, q(x)=1
2.
Loop for n times doing following:
1.
2.
f(x)=f(x)+q(ki)-1q(x)(wi-f(ki))
q(x)=(x-ki)q(x)
Federico Tessmann, Ferienakademie TUM
Newton´s interpolation formula
Let f i  f (ki )
Newton´s interpolation formula claims that there
exist constants  i such that
f ( x)  0  1 ( x  k1 )  2 ( x  k1 )(x  k2 )  
In fact, 0  f1
and
1 is the solution of
f 2  f1  1 (k2  k1 )
Federico Tessmann, Ferienakademie TUM
Newton´s interpolation formula
Then
f ( 2) ( x)  f1 
f 2  f1
( X  k1 )
k 2  k1
And more generally
fl  fl 1  l (kl  k1 )(kl  kl 1 )
Solving the  i gives
f (l ) ( x)  f l 1 ( x)  ( f l  f l 1 (kl )) 
Federico Tessmann, Ferienakademie TUM
( x  k1 ) ( x  kl 1 )
(kl  k1 ) (kl  kl 1 )
Multivariate Interpolation

Dense Interpolation

Probabilistic Sparse Interpolation

Deterministic Sparse Interpolation without
degree bounds
Federico Tessmann, Ferienakademie TUM
Multivariate dense Interpolation
We are given a black box  P with a degree bound „d“ for
the polinom P(xi,..,xn)
So we can assume that P has the form
P( X1,, X n )  P0 ( X 2 ,, X n ) X1d   Pd ( X 2 ,, X n )
 P ( x10 , x20 ,, xn 0 ) 


  P( X 1 , x20 , , xn 0 )
 P ( x1d , x20 ,, xn 0 )
Federico Tessmann, Ferienakademie TUM
Multivariate dense Interpolation
So we get the values of Pi ( x20 ,, X n0 )
which are the coeficients found by interpolating P
on X1
By doing this procedure we compute recursively
P(X1,...,Xk,x(k+1)0,...,xn0)
Federico Tessmann, Ferienakademie TUM
Multivariate dense Interpolation
 P ( x10 , x20 ,  , xn 0 ) 




  P ( X 1 , x20 ,  , xn 0 ) 

 P ( x1d , x20 ,  , xn 0 )



  P ( X 1 , X 2 , x30 ,  , xn 0 )
 P ( x10 , x2 d ,  , xn 0 ) 




  P ( X 1 , x 2 d ,  , xn 0 ) 

 P ( x1d , x2 d ,  , xn 0 )


Federico Tessmann, Ferienakademie TUM
The complexity of the dense
interpolation
Let I(d) be the complexity of interpolating d+1 values to produce a
univariate plynomial of degree „d“ and Nk the complexity for the first k
variables
N1  I ( d )
N 2  (d  1) N1  (d  1) I (d )

N k  (d  1) N k 1  (d  1) k 1 I (d )

N k  k (d  1)
k 1
I (d )
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation

Formal Presentation

Example

Analysis
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation
Assume we want to dermine P(X1,..., Xn) which is an element of L[X]
where L is a field of cardinal q and the degree of each Xi is bounded
by „d“ and there are no more than T non-zero monomials
 e
 e
P( X1,, X n )  c1 X 1   cT X T
Def :


skelP  (e1 ,, eT )

 
 

skelk P  e1 , , ek  e  skelP, e  e1 ,, ek , 
skelP( X 1 , , X k , xk 1 , , xn )  skelk P( X 1 ,, X n )
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation
Def:
( x1 ,, xn )
is a precise evaluation point if:
k , skelP( X 1 ,, X k , xk 1 ,, xn )  skelk P( X 1 ,, X n )
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation
The probability by wich ( x1 ,, xn ) is an imprecise evaluation point:
For each k we can write
 e
 e
1k
P( X1,, X n )  c1k ( X k 1,, X n ) X   cTk ( X k 1,, X n ) X Tk
It is an imprecise evaluation point if one of the cik = 0 (n  k )dT
q
And the probability that this happends is no more than
(n  1)dT
dT n(n  1)dT


q
q
2q
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation

y  ( y1 ,, yk 1 )
Given is a k-1 tuple
The
probability
that


 ei  e j  
y  y , ei  e j
is 0 if we are we are working on a field of
d
characteristic 0 or at least (2  3 pk 1 )
When working on a field of q elements the
probability is bounded by
Federico Tessmann, Ferienakademie TUM
T (T  1)d
2q
Probabilistic Sparse Interpolation
So the following probability is then one that underlines the
Probabilistic Sparse Interpolation
T (T  1)d n(n  1)dT dT (t  2n 2 )



2q
2q
2q
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation
Assume we want to dermine P(X1,..., Xn) which is an element of L[X]
where L is a field of cardinal q and the degree of each Xi is bounded
by „d“ and there are no more than
T non-zero monomials


e1
eT
1
n
1
T


P( X ,, X )  c X   c X
As in the dense interpolation we Interpolate
P( X1,, X n )  P0 ( X 2 ,, X n ) X1d   Pd ( X 2 ,, X n )
 P ( x10 , x20 ,, xn 0 ) 


  P( X 1 , x20 , , xn 0 )
 P ( x1d , x20 ,, xn 0 )
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation
At the kth stage the first computation gives us:


P( X1,, X k 1, xk 0 ,, xn0, )  p10 X e1   pT 0 X eT
We then assume that


skelk P  e1 ,, ek 
The probability of that being the right skeleton is 1 
(n  k  1)(n  k )dT
2q
We then pick a (k-1) tuple

y  ( y1 ,, yk 1 )
And we set up the following transposed Vandermonde system of linear
ecuations
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation
P(1, ,1, xkj )  p1 j    pTj
 e1
 eT
P( y1 , , yk 1 , xkj )  p1 j y    pTj y
 2 e1
 2 eT
2
2
P( y1 , , yk 1 , xkj )  p1 j y    pTj y

P( y , , y
T
1
T
k 1
 Te1
 TeT
, xkj )  p1 j y    pTj y
So each of the pij can be computed using O(n2) and we can avoid
computing the other interpolations
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation
The probability that the Vandermonde system of
equation is non-singular is bounded by
T (T  1)d
1
2q
Federico Tessmann, Ferienakademie TUM
Probabilistic Sparse Interpolation
So we get

e1
P( X1,, X k 1, xkj )  p1 j X   pTj X
for each k
Then we solve pi  Pi ( X k )
trough the dense interpolation
We then expand it and we get

eT

e0
P( X1,, X k , x( k 1)0 ,, xn0 )  p10 X   pT 0 X
And we are ready to compute the (k+1)th stage
Federico Tessmann, Ferienakademie TUM

eT
Probabilistic Sparse Interpolation
Example
Lets assume we are given a Black Box representing the following
polinom
P( X , Y , Z )  X 5 Z 2  X 5 Z  XY 4  XYZ5  Y 5 Z
 P ( x0 , y0 , z0 )



5

  c0 X  c1 X  c2 

 P ( x5 , y0 , z0 ) 

 P ( x0 , y1 , z0 ) 


5

  c3 X  c4 X  c5   d X 5  (d Y 4  d Y ) X  d Y 5
0
1
2
3
 P ( x2 , y1 , z0 )




 P ( x0 , y5 , z0 )


5

  c15 X  c16 X  c17 

 P ( x2 , y5 , z0 )

Federico Tessmann, Ferienakademie TUM
Deterministic Sparse Interpolation
without degree bounds
Given are a bound on the number of non-zero
terms „T“ and the number of variables „n“
We want
 to compute


P( X )  c1 X   cT X
e1
eT
By choosing a distinct
 e prime for each Xi then the
quantities mi  X 1 will all be distinct.
j
j
j
j
v j  P( p )
Let p  ( x1 ,, xn )
j
j
Then we get: c1m1   cT mT  v j
Federico Tessmann, Ferienakademie TUM
Deterministic Sparse Interpolation
without degree bounds
c1    cT  v0
1
1
c1m1    cT mT  v1

t 1
t 1
c1m1    cT mT  vt 1
The rank of the system of equations is exactly the number of non-zero
monomials in P
This could be easily done by taking the first T equations and computing their
rank which requires O(T3)
Federico Tessmann, Ferienakademie TUM
Deterministic Sparse Interpolation
without degree bounds
Let
Q(Z )  (Z  m1 )(Z  mt )  qt Z t   q0
and consider
i
i
c
Q
(
m
)

(
c
m



c
m
 j j  1 1
t t )qi
1 j t
so
1i t
0  vt  vt 1qt 1   v0q0
i k
i k
c
m
Q
(
m
)

(
c
m



c
m
 1 1
consider also  j j
j
t t )qi
k
1 j t
Federico Tessmann, Ferienakademie TUM
1i t
Deterministic Sparse Interpolation
without degree bounds
Then we get the following Toeplitz system of linear
equations
v0 q0    vt 1qt 1  vt
v1q0    vt qt 1  vt 1

vt 1q0    v2t  2 qt 1  v2t 1
Federico Tessmann, Ferienakademie TUM
Deterministic Sparse Interpolation
without degree bounds
v1  vt 1 
 v0


vt 
 v1 v2 

 

 



v

 t 1 vt  v2t  2 
1

1   c1
 1

 
m2  mt   0
 m1

 





 t 1


t 1
t 1  
m
m

m
2
t
 1
 0
0
c2

0
0  1 m1
 
0  1 m2


 


 
 ct  1 mt


 m1t 1 

t 1
 m2 

 
 mtt 1 
So the system is non-singular if the mi are distinct
Federico Tessmann, Ferienakademie TUM
Deterministic Sparse Interpolation
without degree bounds
v0 q0    vt 1qt 1  vt
v1q0    vt qt 1  vt 1

vt 1q0    v2t  2 qt 1  v2t 1
•So the system can be solved by Gaussian elimination O(t3)
•So then we get Q(Z)
•In order to find the mi we just need to find the zeroes of Q which are in fact
positive integers making this procedure much easier
•Knowing the mi, the ei can easily be determined by factoring each of the mi,
which is in fact very easy because the possible divisors are the first n primes
allready known
•By knowing the mi it is allso easy to compute the ci just by solving the
Vandermonde system, formed by the first t equations
Federico Tessmann, Ferienakademie TUM
The End
Questions?
Bibliography: Richard Zippel
Federico Tessmann, Ferienakademie TUM
Appendix
pseudo-codes for the interpolation algorithms
Federico Tessmann, Ferienakademie TUM
Code for Vandermonde Matrices
SolveVanderm onde(k[n], w[n]){
Q( Z ) : ( Z  k[1]) ( Z  k[n]);
x[] : newArray[n]  [0;...;0];
do{
for i : 1 to n
Qi ( Z ) : PolyQuotient(Q( Z ), Z  k[i ]);
Pi ( Z ) : Qi ( Z ) / Qi (k[i ]);
for
j : 1 to n
do{
x[ j ] : x[ j ]  w[i ]  coef ( Pi ( Z ), Z j 1 );
}
}
return( x);
Federico Tessmann, Ferienakademie TUM
Code for Lagrange Interpolation
LagrangeInterpolation(k[n], w[n]){
Q( Z ) : ( Z  k[1]) ( Z  k[n]);
f : newPolinomial  0;
for i : 1 to n
do{
Qi ( Z ) : PolyQuotient(Q( Z ), Z  k[i ]);
Pi ( Z ) : Qi ( Z ) / Qi (k[i ]);
f : f  w[i ]  Pi ( Z );
}
return( f );
Federico Tessmann, Ferienakademie TUM
Code for Newton Interpolation
NewtonInterpolation(k[n], w[n]){
q : newPolinomial  1;
f : newPolinomial  0;
for i : 1 to n
do{
f : f  q(k[i ])1  q( X )  ( w[i ]  f (k[i ]));
q : ( X  k[i ])q;
}
return( f );
Federico Tessmann, Ferienakademie TUM