Нахождение собственных чисел и векторов методом Якоби и степенным методом - файл lab6_nm.doc

Нахождение собственных чисел и векторов методом Якоби и степенным методом
Скачать все файлы (13.7 kb.)

Доступные файлы (2):
n1.cpp
lab6_nm.doc49kb.17.06.2010 22:20скачать

lab6_nm.doc



Завдання
Знайти власні числа методом Якобі та обчислити відповідні власні вектори. Обчислити степеневим методом найбільше та найменше за модулем власні числа матриці на відповідні їм власні вектори. Матриця:

6.7 1.26 0.81 1.18

1.26 3.72 1.3 0.16

0.81 1.3 5.88 2.1

1.18 0.16 2.1 5.66
Текст програми

#include

#include

#include

#include

#include

#include

#define max 10

#define eps 0.00001

using namespace std;
void inp_mf (double a[max][max], int *n, char file[max])

{

//reading matrix from file

int flag = 0, i = 0, j = 0;

{

ifstream f(file);

if (!f)

{

cout << "\n\nFile not found.";

getchar();

getchar ();

exit (1);

}

f >> *n;

if (*n > max || *n < 1)

{

cout << "Matrix dimensions overflow.";

getchar();

exit (1);

};

for (i = 0; i < *n; i++)

for (j = 0; j < *n; j++)

{

f >> a[i][j];

if (fabsl(a[i][j]) > eps) flag = 1;

}

if (flag == 0)

{

cout << "Attention! Zero matrix.";

getchar();

exit (1);

}

for (i = 0; i < *n; i++) f >> a[i][*n];

f.close();

}

}
void inp_mk ()

{

int i = 0, j = 0, n = 0;

double t = 0;

//reading matrix from keyboard and writing to file

{

ofstream f("inp_matr.txt");

if (!f)

{

cout << "File not found and cannot be created.";

getchar();

exit (1);

}

cout << "Enter the dimensions of matrix -> ";

cin >> n;

if (n > max || n < 1)

{

cout << "Matrix dimensions overflow";

getchar();

exit (1);

};

f << n <<"\n";

for (i = 0; i < n; i++)

{

for (j = 0; j < n; j++)

{

cout << "a[" << i + 1 << "][" << j + 1 << "] -> ";

cin >> t;

f << t << ' ';

}

f << "\n";

}

f.close();

}

}
void printm (double a[max][max], int n, ofstream& f)

{

//printing matrix

int i = 0, j = 0;

for (i = 0; i < n; i++)

{

for (j = 0; j < n; j++) f << a[i][j] << " ";

f << "\n";

}

}
void mult_matr (double a[max][max], double b[max][max], int n, double ab[max][max])

{

int i=0, j=0, k;

for ( i=0; i
for ( j=0; j
{ ab[i][j] = 0;

for ( k=0; k
ab[i][j] = ab[i][j]+a[i][k]*b[k][j];

}

}

double normaF (double a[max][max], int n)

{

int i, j = 0;

double s = 0;

for (i = 0; i
for (j = 0; j
s = s + pow(a[i][j],2);

return(pow(s,0.5));

}
double delta (double a[max][max], int n)

{

int i = 0;

double s = 0;

for (i = 0; i
s = s + pow(a[i][i],2);

return(s);

}
double w (double a[max][max], int n)

{

int i, j = 0;

double s = 0;

for (i = 0; i
for (j = i+1; j
s = s + pow(a[i][j],2);

return(s);

}
void JacobiMethod(double a[max][max], int n)

{

ofstream f ("out_matr.txt");

f.precision (7);

f.width (14);

double alpha, beta, gama, dzeta, t, c, s, Q[max][max], Q1[max][max], A[max][max], a1[max][max];

int i, j, k, i_max, j_max = 0;

for (i=0; i
for (j=0; j
{

A[i][j] = a[i][j];

if ( i == j ) Q[i][j] = 1;

else Q[i][j] = 0;

}

f<<" NormaF = "<
int p = 0;

do

{

p++;

f<<"\n Iteration number "<

gama = 0;

for (i = 0; i
for (j = i+1; j
if (fabs(gama) < fabs(a[i][j]))

{ gama = a[i][j];

i_max = i;

j_max = j;

alpha = a[i][i];

beta = a[j][j];

}

dzeta = (beta - alpha)/(2*gama);

if (dzeta <= 0) t = 1/(dzeta - sqrt(dzeta*dzeta+1));

else t = 1/(dzeta + sqrt(dzeta*dzeta+1));

c = 1/sqrt(t*t+1);

s = t*c;

f<<"c = "<
f<<"s = "<
f<<"alpha (a["<
f<<"beta (a["<
f<<"gama (a["<
f<<"dzeta = "<
f<<"t = "<
for ( i=0; i
for ( j=0; j
{

a1[i][j] = a[i][j];

Q1[i][j] = Q[i][j];

}

for ( k=0; k
{

a1[k][i_max] = a[k][i_max]*c - a[k][j_max]*s;

}

for ( k=i_max+1; k<=j_max; k++)

{

a1[i_max][k] = a[i_max][k]*c - a[k][j_max]*s;

}

for ( k=j_max+1; k
{

a1[i_max][k] = a[i_max][k]*c - a[j_max][k]*s;

}

for ( k=0; k
{

a1[k][j_max] = a[k][i_max]*s + a[k][j_max]*c;

}

for ( k=i_max+1; k<=j_max; k++)

{

a1[k][j_max] = a[i_max][k]*s + a[k][j_max]*c;

}

for ( k=j_max+1; k
{

a1[j_max][k] = a[i_max][k]*s + a[j_max][k]*c;

}

a1[i_max][i_max] = a[i_max][i_max]*c*c + a[j_max][j_max]*s*s - 2*c*s*a[i_max][j_max];

a1[j_max][j_max] = a[i_max][i_max]*s*s + a[j_max][j_max]*c*c + 2*c*s*a[i_max][j_max];

a1[i_max][j_max] = 0;

for ( k=0; k
for ( i=k+1; i
a1[i][k] = a1[k][i];

for ( i=0; i
for ( j=0; j
a[i][j] = a1[i][j];

f<<" A"<

printm(a, n, f);

printm(Q, n, f);

for ( j=0; j
{

Q1[j][i_max] = c*Q[j][i_max] - s*Q[j][j_max];

Q1[j][j_max] = s*Q[j][i_max] + c*Q[j][j_max];

}

for ( i=0; i
for ( j=0; j
Q[i][j] = Q1[i][j];

f<<" Q"<

printm(Q, n, f);

f<<" Delta = "<
f<<" NormaF = "<
f<<"W = "<
}

while ( w ( a, n ) >= pow(eps,2) );

mult_matr(A,Q,n,a1);

for ( i=0; i
{

f<<"\n"<
f<<" eigenvector :";

for ( j=0; j
f<
f<<"\n Av-lv : ";

for ( j=0; j
f<
}

}
double norma (double x[max], int n)

{

double S = 0;

for (int i = 0; i < n; i++) S = S + pow (x[i], 2);

return pow (S, 0.5);

}
void stepen (double matr [max][max], int n, double *l, double x[max], int s)

{

double y[max], p[max];

int i, j, k = 0;

for (i = 0; i < n; i++) x[i] = 1;

int flag = 0;

while (flag == 0)

{

flag = 1;

for (i = 0; i < n; i++)

{

y[i] = 0;

for (j = 0; j < n; j++) y[i] = y[i] + x[j]*matr[i][j];

}

*l = norma (y, n)/norma(x, n);

if (x[0]*y[0] < 0) *l = -*l;

for (i = 0; i < n; i++) if (fabsl(x[i]-y[i]/(*l)) > eps) flag = 0;

for (i = 0; i < n; i++) x[i] = y[i]/norma(y, n);

}

}

int main()

{

int n = 0, flag = 0, i = 0, j = 0, k = 0;

double a[max][max], A[max][max];

cout << "Please, enter what do you want to do:\n 1. Read matrix slar from file.\n 2. Enter slar from keyboard.\n - > ";

cin >> flag;

char inp[] = "inp_matr.txt";

if (flag == 1) inp_mf (a, &n, inp);

else

{

inp_mk ();

inp_mf (a, &n, inp);

}

for ( i=0; i
for ( j=0; j
A[i][j] = a[i][j];

JacobiMethod(a, n);

cout << "\n\n Power method \n";

double v1[max], v2[max];

for ( i=0; i
v1[i] = 1;

double lmax = 0;

stepen (A, n, &lmax, v1, 1);

for (i = 0; i < n; i++)

{

v2[i] = 0;

for (j = 0; j < n; j++) v2[i] = v2[i] + A[i][j]*v1[j];

}

cout << "\n Proper number_max = " << lmax;

cout << "\n Eigenvector:";

for (i = 0; i < n; i++) cout << v1[i] << " ";

cout << "\n Residual vector:";

for (i = 0; i < n; i++) cout << v2[i] - v1[i]*lmax << " ";

for (i = 0; i < n; i++)

for (j = 0; j < n; j++)

{

if (i != j) a[i][j] = A[i][j];

else a[i][j] = A[i][j] - lmax;

}

double lmin = 0;

stepen (a, n, &lmin, v1, 2);

for (i = 0; i < n; i++)

{

v2[i] = 0;

for (j = 0; j < n; j++) v2[i] = v2[i] + A[i][j]*v1[j];

}

lmin = norma (v2, n)/norma(v1, n);

if (v1[0]*v2[0] < 0) lmin = -lmin;

cout << "\n\n Proper number_min = " << lmin;

cout << "\n Eigenvector:";

for (i = 0; i < n; i++) cout << v1[i] << " ";

cout << "\n Residual vector:";

for (i = 0; i < n; i++) cout << v2[i] - v1[i]*lmin << " ";

getchar();

}

Результати роботи

Power method
Proper number_max = 9.25747

Eigenvector:0.562094 0.277295 0.571114 0.530081

Residual vector:-4.9032e-005 -1.03899e-005 2.92896e-005 2.58737e-005

Proper number_min = 2.57038

Eigenvector:0.244928 -0.754764 0.476611 -0.378395

Residual vector:8.01378e-006 -3.41725e-005 -2.73742e-005 3.88716e-005
Jacobi method

1. proper number 5.849002

eigenvector :0.7802453 0.1369973 -0.4812273 -0.3753256

Av-lv : 0.003925614 0.005793835 0.006994363 0.002968846

2. proper number 2.570739

eigenvector :-0.2474624 0.7533496 -0.4797225 0.3756235

Av-lv : -0.01795573 -0.009546406 -0.01980769 -0.01848757

3. proper number 9.261516

eigenvector :0.5552177 0.2949443 0.5812657 0.5165958

Av-lv : 0.02983925 -0.09656342 -0.04755405 0.06247804

4. proper number 4.278744

eigenvector :0.1473532 -0.5715802 -0.4476793 0.6716892

Av-lv : 0.06656181 0.03051914 0.07 0.07007236

Учебный текст
© perviydoc.ru
При копировании укажите ссылку.
обратиться к администрации