Note
Soit
un polynôme que l’on souhaite évaluer en \(\alpha\).
Le schéma d’évaluation de Hörner consiste à écrire
ce qui permet de calculer la valeur \(P(\alpha)\) à l’aide de \(n\) multiplications et additions.
Remarquer la similitude avec l’exponentiation rapide.
Horner(P,a) qui évalue le polynôme \(P\) en \(a\).
La tester sur le polynôme formel P:=sum(a[k]*X^(8-k),k,0,8).On définit précisément la suite de Hörner
de sorte que \(u_n(\alpha)=P(\alpha)\).
Montrer que l’on a
Si l’on aligne les puissances de \(X\) dans le membre de droite (en mettant \(P(\alpha)=u_n(\alpha)\) comme terme d’indice \(k=n\) dans la première somme), on obtient
ce qui vaut bien \(\sum_{k=0}^n a_k X^{n-k}=P(X)\) d’après la récurrence des \(u_k\).
En déduire une fonction horner_divrem(P,a) qui renvoie le quotient et le
reste de la division euclidienne de \(P(X)\) par \((X-a)\).
D’après la question précédente, le schéma de Hörner calcule les coefficients du quotient. Il suffit de les conserver.
horner_divrem(P,a):={
local n,k,u,Q;
n := degree(P);
Q := 0;
u := coeff(P,n);
for(k:=n-1; k>=0; k--) {
Q := x*Q + u;
u := a*u + coeff(P,k);
}
return( Q,u );
}
}
Vérifier les résultats en utilisant la fonction xcas quorem.
horner_divrem(X^8-3,5)
quorem(X^8-3,X-5)
On considère une suite de points \((x_0,\dots x_n)\) et une suite de valeurs \((y_0,\dots y_n)\).
Donner une valeur à n, et utiliser seq et rand pour créer deux
listes lx et ly de \(n+1\) éléments aléatoires dans \([-10,10]\).
n := 10;
lx := [ seq(rand(-10,10),k=0..n) ];
ly := [ seq(rand(-10,10),k=0..n) ];
On rappelle que l’interpolateur de Lagrange de \((y)\) aux points \((x)\) est l’unique polynôme \(P(X)\) de degré inférieur ou égal à \(n\) tel que \(P(x_i)=y_i\).
On a la formule explicite
Écrire une fonction LagrangeNaif(lx,ly) qui calcule
\(P(X)\) à l’aide de la formule explicite.
LagrangeNaif(lx,ly):={
local n,P,k,j;
n := size(lx);
for(k:=0;k<n;k++) {
local(Lk:=1);
for(j:=0;j<n;j++) {
if(j==k) continue;
Lk *= (x-lx[j])/(lx[k]-lx[j]);
}
P += ly[k]*Lk;
}
return(P);
}
Quelle est la complexité de l’algorithme en fonction de \(n\) ?
L’algorithme effectue une double boucle contenant la mutiplication d’un polynôme de taille \(n\) par un monôme, soit \(O(n^3)\) opérations.
Notons \(V(x_0,\dots x_n)\) la matrice de passage de la base \((L_0(X),\dots L_n(X))\) des interpolateurs de Lagrange aux points \((x_i)\), à la base canonique \((1,X,\dots X^n)\).
(la matrice de passage de \(\mathcal B\) à \(\mathcal B'\) est la matrice de l’application identité \((E,\mathcal B') \to (E,\mathcal B)\)).
Définir cette matrice, et déterminer l’instruction xcas permettant de la créer rapidement.
L’utiliser pour calculer les coefficients de l’interpolateur de Lagrange \(P(X)\) dans la base canonique.
En déduire une fonction LagrangeVdm(lx,ly). Quelle est sa complexité ?
On sait que \(P(X)\) a pour coordonnées \((y_0,\dots y_n)\) dans la base des interpolateurs \(L_i\). Pour repasser dans la base canonique, il suffit de multiplier ce vecteur par l’inverse de \(V(x_0,\dots x_n)\).
D’où la fonction suivante, dans laquelle on a transposé l’équation et on renverse la liste à la fin pour avoir l’écriture usuelle en partant du coefficient de \(X^n\):
LagrangeVdm(lx,ly) := {
return revlist(ly*(vandermonde(lx)^-1));
}
La complexité est à nouveau de \(O(n^3)\) pour l’inversion matricielle.
Pour économiser les calculs, on propose la stratégie suivante :
calculer d’abord le produit complet
puis calculer chacun des \(L_i\) à l’aide d’une division selon le schéma de Hörner.
Au fur et à mesure de l’étape précédente, on peut construire \(P(X)\).
Écrire une fonction LagrangeHorner(X,Y) selon ce principe.
En partant de \(A(X)\), \(L_i(X)\) vaut
On commence par modifier la fonction horner_divrem pour
qu’elle calcule à la fois \(A_i(X)\) et \(A_i(x_i)\) (qui vaut
simplement \(u_{n-1}\)).
horner_divrem_A(P,a):={
local n,k,u,Q;
n := degree(P);
Q := 0;
u := 0;
for(k:=n; k>0; k--) {
u := a*u + coeff(P,k);
Q := x*Q + u;
}
return( Q,u );
}
LagrangeHorner(lx,ly) := {
local n,A,P,k;
n := size(lx);
A := product(x-lx[k],k,0,n-1);
P := 0;
for(k:=0;k<n;k++) {
local Ak,ak;
Ak,ak := horner_divrem_A(A,lx[k]);
P += ly[k]*Ak/ak;
}
return(P);
}
Quelle est la complexité obtenue ?
quadratique
On essaie de résoudre le problème par induction en introduisant la suite de polynômes \(P_k\) tels que \(P_k\) soit l’interpolateur de Lagrange des \((y_0,\dots y_k)\) aux points \((x_0,\dots x_k)\).
On peut écrire
où \(\lambda_{x_0,\dots x_k}\) est le coefficient dominant de \(P_k\) et peut se calculer par récurrence par la formule
en partant des \(\lambda_{x_i}=y_i\).
Écrire une fonction lambda(ly) qui renvoie la suite des
\(\lambda_{x_0,\dots x_k}\) pour \(0\leq k<n\).
En déduire une fonction LagrangeDivise de calcul de l’interpolant
de Lagrange \(P(X)\).
lx := [-3,-1,0,2,3];
ly := [2,3,1,-1,2];
difdiv(lx,ly) := {
n := size(lx)-1;
//fig :=
seq(point(lx[k],ly[k]),k=0..n);
L := ly;
for(d:=0;d<=n;d++) {
seq(plot(L[k],x,lx[k]-0.5,lx[k+d]+0.5),k,0,n);
//fig
for(k:=n-d;k>0;k--) {
// L[k] interpolates x_k,...x_{k+d}
// L[k-1] interpolates x_{k-1},... x_{k+d-1}
L[k-1] := (L[k]*(x-lx[k-1])-L[k-1]*(x-lx[k+d]))/(lx[k+d]-lx[k-1]);
}
}
}:;
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 | lx = [-3,-1,0,2,3]
ly = [1.7,3,2.5,0.5,1.5]
def difdiv(lx,ly):
Pts = point(zip(lx,ly))
n = len(lx)
L = ly[:]
G = []
for d in xrange(n+1):
F = Pts + sum( plot(L[k],x,lx[k]-0.5,lx[k+d]+0.5)
for k in xrange(n-d) )
G.append(F)
for k in xrange(n-d-1):
L[k] = (L[k]*(x-lk[k+d+1]) - L[k+1]*(x-lk[k])/(lx[k]-lx[k+d+1])
return G
animate(difdiv(lx,ly))
|