Метод Леверрье-Фаддеева вычисления коэффициентов характеристического полинома

/*функция, предназначенная она для нахождения следа матрицы - это всего-навсего сумма элементов главной диагонали*/

template <class YourOwnFloatType>

YourOwnFloatType Sp(matrix<YourOwnFloatType> &x)

{

YourOwnFloatType res=(YourOwnFloatType)0;

/*предполагается, что ищется след квадратной матрицы, однако соответствующая проверка не выполняется*/

for(long i=0;i<x.getm();i++)

res+=x[i][i];

return res;//возвращаем след

}

cvector Leverrie_Faddeev_eigenval(cmatrix &x)

{

if(x.getn()!=x.getm())

throw xmsg("Для неквадратных матриц собственные значения не определены");

//создаём пока пустой полином

cpolynom p=x.getm()+1;

/*устанавливаем коэффициент при самой старшей степени в единицу*/

p[x.getm()]=1;

cmatrix a=x,b=x;//рабочие матрицы

cmatrix e(x.getm(),x.getm());//нулевая матрица,

for(long i=0;i<e.getm();i++) //плавно переходящая

e[i][i]=1; //в единичную

/*вычисляем коэффициенты характеристического полинома*/

for(long i=0;i<x.getm();i++)

{

p[x.getm()-(i+1)]=-Sp(a)/(i+1);

b=a+p[x.getm()-(i+1)]*e;

a=x*b;

}

return newton(p);/*находим корни модифицированным методом Ньютона*/

}

//Вариант главной функции для тестирования

void main()

{

cmatrix a="test.z";

cout<<a<<endl;

getch();

cvector l=Leverrie_Faddeev_eigenval(a);

cout<<a<<endl<<l<<endl<<endl;

getch();

cmatrix v=eigenvec(a,l);

cout<<v<<endl;

getch();

for(long i=0;i<v.getm();i++)

{

cmatrix x(v.getm(),1);

for(long j=0;j<x.getm();j++)

x[j][0]=v[i][j];

cout<<a*x<<"="<<endl<<l[i]*x<<endl<<endl;

getch();

}

}