sfact

Paramètres

P

matrice carrée de polynômes réels

Description

Si P est une matrice de polynômes vérifiant P(z) / z^m = P(1/z)' * z^m où n=2*m est le degré maximal (nécessairement pair) de P(z) et si P(z)/z^n est définie positive sur le cercle unité alors sfact(P) renvoie une matrice de polynômes F telle que

P(z) = F(z) * F(1/z)' * z^n

Par cette factorisation les racines de det(P) sont séparées en deux ensembles strictement à l'extérieur ou à l'intérieur du cercle unité, les racines (antistables) de det(F) et les racines (stables) de det(F(1/z)'*z^n), respectivement.

Dans le cas scalaire la condition s'exprime plus simplement : P doit être palindromique (ou réciproque), c'est à dire que les coefficients des monômes de degrés k et n-k doivent être égaux.

En temps continu (voir le troisième exemple ci-dessous), la factorisation peut être appliquée à une matrice P(s) verifiant P(s) = P(-s) et P(s) définie positive sur l'axe imaginaire en utilisant la transformation bilinéaire temps discret/temps continu.

La factorisation spectrale a des applications en théorie des systèmes, contrôle optimal, théorie des reseaux et en automatique. Les détails de l'agorithme se trouvent dans l'article cité dans la bibliographie.

Exemples

// Polynôme unique
z = %z;
p = (z -1/2) * (2 - z) // p est palindromique
w = sfact(p);
w * horner(w, 1/z) * z

// Exemple matriciel
z = %z;
// on définit P vérifiant les conditions
// gtild(F1,'d') est égal à horner(F1, 1/z)'*z^3
F1 = [z-1/2, z+1/2, z^2+2; 1, z, -z; z^3+2*z, z, 1/2-z];
P = F1*gtild(F1,'d'); 
F = sfact(P)
roots(det(P))             // racines de det(P)
roots(det(gtild(F,'d')))  // racines stables
roots(det(F))             // racines anti-stables
clean(P-F*gtild(F,'d'))

// Utilisation en temps continu
s = %s;
z = %z;
p = -3*(s+(1+%i))*(s+(1-%i))*(s+0.5)*(s-0.5)*(s-(1+%i))*(s-(1-%i));
p = real(p);

// p(s) = polynôme dans s^2, cherche un f stable tel que p=f(s)*f(-s)
w = horner(p,(1-z)/(1+z));  // transformation bilinéaire w=p((1-z)/(1+z))
wn = w.num;                 // prend le numérateur
fn = sfact(wn);             // Factorisation
f = horner(fn,(1-s)/(s+1)).num; // et transformation inverse
f = f/sqrt(horner(f*horner(f,-s),0));
f = f*sqrt(horner(p,0));        // normalisation
roots(f)    // f est stable
clean(f*horner(f,-s)-p)    // f(s)*f(-s) est égal à p(s)

Voir aussi

• gtild — tilde operation
• fspecg — stable factorization of continuous time dynamical systems

Bibliographie

J. Ježek, V. Kučera, Efficient Algorithm for Spectral Factorization, IFAC Proceedings Volumes, Volume 17, Issue 2, 1984, Pages 257-262, ISSN 1474-6670, https://doi.org/10.1016/S1474-6670(17)60979-0