使用QR分解 求一元四次方程的根

        在求特征值的时候，通过QR迭代后就是一个拟上三角矩阵，但不一定是上三角矩阵。

        在一定条件下，由QR算法生成的序列{Ak}收敛为Schur分块上三角形，对角块按特征值的模从大到小排列。但有特殊情况，当收敛结果为Schur分块上三角形时，序列{Ak}的对角块以上的元素以及2阶块的元素不一定收敛。需要继续求得该二阶矩阵的特征值，补充到原矩阵的特征值中即可。

        其matlab代码如下：

% a1*x^4+a2*x^3+a3*x^2+a4*x+a5 = 0

A = [-a2/a1, -a3/a1, -a4/a1, -a5/a1;
     1,0,0,0;
     0,1,0,0;
     0,0,1,0];
i = 0;
while i<1
    i = i+1;
    B1 = A  ;
    r1 = sqrt(B1(1,1)^2+B1(2,1)^2);
    c1 = B1(1,1)/r1;
    s1 = B1(2,1)/r1 ;   
    R21 = [c1,s1,0,0;
           -s1,c1,0,0;
           0,0,1,0;
           0,0,0,1] ; 
    B2 = R21*B1;
    r2 = sqrt(B2(2,2)^2+B2(3,2)^2);
    c2 = B2(2,2)/r2;
    s2 = B2(3,2)/r2  ;
    R32 = [1,0,0,0;
           0,c2,s2,0;
           0,-s2,c2,0;
           0,0,0,1]  ;
    B3 = R32*B2 ;
    r3 = sqrt(B3(3,3)^2+B3(4,3)^2);
    c3 = B3(3,3)/r3;
    s3 = B3(4,3)/r3   ; 
    R43 = [1,0,0,0;
           0,1,0,0;
           0,0,c3,s3;
           0,0,-s3,c3];
    B4 = R43*B3  ;  
    Q = R21'*R32'*R43';
    R = R43*R32*R21*B1;
    A1 = R*Q;
    norm(A-Q*R);
    A = A1 ;
end

g1 = A(1,1);
g2 = A(1,2);
g3 = A(2,1);
g4 = A(2,2);
g5 = A(3,3);
g6 = A(3,4);
g7 = A(4,3);
g8 = A(4,4);
m = (g1+g4)/2;
n = g1*g4-g2*g3;
m1 = (g5+g8)/2;
n1 = g5*g8-g6*g7;

x1 = m+sqrt(m^2-n)
x2 = m-sqrt(m^2-n)
x3 = m1+sqrt(m1^2-n1)
x4 = m1-sqrt(m1^2-n1)