sfact
離散時間スペクトル分解
呼び出し手順
F = sfact(P)
パラメータ
- P
実多項式の正方行列.
説明
Pのスペクトル分解 Fを求めます.
P はPの各根が単位円に関する
鏡像となるような多項式行列です.
単位円上に根がある時,問題は特異となります.
sfact(P) は,安定ではない以下のような
多項式行列F(z)を返します.
P = F(z)* F(1/z) *z^n
スカラー多項式の場合,特殊なアルゴリズムが実装されています. アルゴリズムはKuceraの本によるものです.
例
// Simple polynomial example p = (%z -1/2) * (2 - %z) w = sfact(p); w*(horner(w, 1/%z)).num
// matrix example z = %z; 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)) roots(det(gtild(F,'d'))) //The stable roots roots(det(F)) //The antistable roots clean(P-F*gtild(F,'d'))
// Example of continuous time use s = %s; 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) = polynomial in s^2 , looks for stable f such that p=f(s)*f(-s) w = horner(p,(1-s)/(1+s)); // bilinear transform w=p((1-s)/(1+s)) wn = w.num; // take the numerator fn = sfact(wn); f = horner(fn,(1-s)/(s+1)).num; // Factor and back transform f = f/sqrt(horner(f*gtild(f,'c'),0)); f = f*sqrt(horner(p,0)); // normalization roots(f) // f is stable clean(f*gtild(f,'c')-p) // f(s)*f(-s) is p(s)
| Report an issue | ||
| << rowcompr | Polynomials | simp >> |