Собственные значения матрицы

[krasy]

Написала прогу, которая должна выводить собственные значения, но она выходить из цикла раньше времени на первой итерации. Исходная матрица A преобразовывается в разложение вида StS и для того чтобы найти собственные значения этой матрицы используется алгоритм:
A=st*s;
T=s*st
Идет проверка матрицы T на то, чтобы элементы выше главной диагонали были <= заданной точности eps. Если эти числа больше точности, то уже матрицу Т преобразовывают в разложение, то есть она становится на место матрицы A и тд..когда уже проходит условие то собственные значения это элементы главной диагонали.
Так вот, у меня цикл, который проверяет данное условие. Но он считает только один раз и выходит из цикла. несмотря на то, что не проходит условие.
Помогите пожалуйста разобраться...

#include<iostream>
#include<conio.h>
#include<math.h>
#include<cmath>

using namespace std;

int n,r=0,h=1;
double a[100][100],b[100];

void input();
void sobs();
void menu();

int main()
{
 menu();
 getch();
 return 0;    
}


void input()
{ 
cout<<"Vvedite n:";
cin>>n;
cout<<"Vvedite A:"<<endl;
for(int i=1;i<=n;i++)
 for(int j=1;j<=n;j++)
 {cout<<"a["<<i<<","<<j<<"]=";
  cin>>a[i][j];}
cout<<"Vvedite B:"<<endl;
for(int i=1;i<=n;i++)
{cout<<"b["<<i<<"]=";
  cin>>b[i];}
  
cout<<endl;
cout<<"A:"<<endl;
for(int i=1;i<=n;i++){
 for(int j=1;j<=n;j++)
 cout<<a[i][j]<<"  ";
 cout<<endl;}
 
 cout<<"B:"<<endl;
for(int i=1;i<=n;i++)
cout<<b[i]<<" ";
cout<<endl;
}

void sobs()
{
 double s[100][100],p[100][100],t[100][100];
 double st[100][100],stt[100][100];
 double d[100][100],eps;
 
 cout<<"eps: ";
  cin>>eps;
double dd,f,q,l,m,c=0,det=1;  


do
{
 s[1][1]=sqrt(a[1][1]);
 if(a[1][1]>0) d[1][1]=1;
 if(a[1][1]<0) d[1][1]=-1;
 
 for(int j=2;j<=n;j++)
 s[1][j]=a[1][j]/(s[1][1]*d[1][1]);

for(int i=2;i<=n;i++)
{f=0;q=0;
 for(int k=1;k<i;k++)
  {f+=s[k][i]*s[k][i]*d[k][k];
  q+=fabs(s[k][i])*fabs(s[k][i])*d[k][k];}
  
  dd=a[i][i]-f;
  if((dd>0)||(dd==0)) d[i][i]=1;
  if(dd<0) d[i][i]=-1;
  
  s[i][i]=sqrt(fabs(a[i][i]-q));
  
  for(int j=i+1;j<=n;j++)
  {l=0;
  for(int k=1;k<=(i-1);k++)
  l+=s[k][i]*d[k][k]*s[k][j];
  s[i][j]=(a[i][j]-l)/(s[i][i]*d[i][i]);
  }
  }
  
  cout<<endl<<endl;
 
 for(int i=1;i<=n;i++)
 {for(int j=1;j<=n;j++)
  st[i][j]=s[j][i];
 }
 //ïðîèçâåäåíèå ìàòðèö sTs
 
 for(int i=1;i<=n;i++)
 {for(int j=1;j<=n;j++)
    stt[i][j]=st[i][j]*d[i][i];
 }
  
  
  for(int i=1;i<=n;i++)
  {for(int j=1;j<=n;j++)
   {a[i][j]=0;
   for(int k=1;k<=n;k++)
    {
    a[i][j]+=stt[i][k]*s[k][j];
    t[i][j]+=s[i][k]*stt[k][j];
    }
   }
  }
  
  for(int i=1;i<=n;i++)
  {for(int j=1;j<=n;j++)
    a[i][j]=t[i][j];}
    
  for(int i=1;i<=n;i++)
  {for(int j=i+1;j<=n;j++)
   if(a[i][j]<eps) r++;}
}
while(r>0);

 
 for(int i=1;i<=n;i++){
 for(int j=1;j<=n;j++)
 cout<<a[i][j]<<" ";cout<<endl;}
 
 cout<<endl;
  
 for(int i=1;i<=n;i++)
   cout<<"l"<<i<<": "<<a[i][i]<<endl;
 
  /*for(int i=1;i<=n;i++)
  {for(int j=1;j<=n;j++)
   cout<<round(p[i][j]*10)/10.<<"  ";
   cout<<endl;}
   */
   
 
}

void menu()
{
char ch;
int number;
ch=13;
 while(ch!=27)
	{cout<<"\n1. Input;";
	 cout<<"\n2. Sobstv;";
	 cout<<"\n3. Exit;\n";
	 cin>>number;
	switch(number){
		case 1: input(); break;
		case 2: sobs(); break;
		case 3: ch=27; break;
		default: cout<<"\nNet takovo punkta";
	if(ch!=27) ch=_getch();
	}  
}
}

[p51x]

Для начала: массивы нумеруются с 0

[krasy]

p51x , я поменяла 1 на 0...проблема все та же остается, не думаю что от того как я пронумеровала массив зависит то, что у меня вылетает из цикла раньше времени

[p51x]

А <= поменяли на < ?

[krasy]

конечно

[p51x]

  for(int i=1;i<=n;i++)
  {for(int j=1;j<=n;j++)
   {a[i][j]=0;
   for(int k=1;k<=n;k++)
    {
    a[i][j]+=stt[i][k]*s[k][j];
    t[i][j]+=s[i][k]*stt[k][j];
    }
   }
  }
  
  for(int i=1;i<=n;i++)
  {for(int j=1;j<=n;j++)
    a[i][j]=t[i][j];}

Зачем вы первый раз считали a[i][j], если не используете и перезаписываете?

[krasy]

я как-то проверяла по этой формуле правильность разложения и видимо забыла убрать

[p51x]

Вот мы и дошли до условия выхода из цикла. В первом посте вы пишите, что цикл вам надо прервать, если все элементы выше главной диагонали меньше точности. А в коде пишите:

 for(int i=1;i<=n;i++)
  {for(int j=i+1;j<=n;j++)
   if(a[i][j]<eps) r++;}
}
while(r>0);

Если какой-то элемент меньше точности...

[krasy]

Не понимаю..
мне нужно, чтобы каждый элемент матрицы выше главной диагонали был меньше точности...разве я не правильно написала?

[p51x]

А обнуление r где? Чему там элементы массива равны?