1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
| double eps = 0.0001;
int Rows = grid1.RowDefinitions.Count;
int m1columns_m1rows = comboBox.SelectedIndex + 2;
double[] Xk = new double[m1columns_m1rows];//Корни
double[] Rk = new double[m1columns_m1rows];//вектор невязки. Отличие полученных свободных членов от требуемых
double[] Sz = new double[m1columns_m1rows];
double[] Zk = new double[m1columns_m1rows];//Вектор спуска
double[] Pk = new double[m1columns_m1rows];
double[] St = new double[m1columns_m1rows];
double[,] t = Transpose(matrix_left);
double[] Sk = new double[m1columns_m1rows];
double alpha, beta, mf;
double Spr, Spr1, Spz;
int i, j;
double tmp,a,b,c;
double[] prev_x = new double[m1columns_m1rows];
/* Вычисляем сумму квадратов элементов вектора right_matrix(правой части)*/
for (mf = 0, i = 0; i < Rows; i++)
{
mf += matrix_right[0, i]* matrix_right[0, i];
}
/* Задаем начальное приближение корней. В Хk хранятся значения корней
* к-й итерации. */
for (i = 0; i < Rows; i++)
{
Xk[i] = 0;
}
/* Задаем начальное значение r0 и z0. */
for (i = 0; i < Rows; i++)
{
for (Sz[i] = 0, j = 0; j < Rows; j++)
{
Sz[i] += matrix_left[i, j] * Xk[j];
}
Rk[i] = matrix_right[0, i] - Sz[i];
Zk[i] = Rk[i];
Sk[i] = Rk[i];
Pk[i] = Rk[i];
}
do
{
/* Вычисляем числитель и знаменатель для коэффициента
* alpha = (pk-1,rk-1)/(Azk-1,sk-1) */
Spz = 0;
Spr = 0;
for (i = 0; i < Rows; i++)
{
for (Sz[i] = 0,St[i] = 0, j = 0; j < Rows; j++)
{
Sz[i] += matrix_left[i, j] * Zk[j];
St[i] += t[i, j] * Sk[j];
}
Spz += Sz[i] * Sk[i];
Spr += Pk[i] * Rk[i];
}
alpha = Spr / Spz; /* alpha */ //величина смещения
/* Вычисляем вектор решения: xk = xk-1+ alpha * zk-1,
вектор невязки: rk = rk-1 - alpha * A * zk-1 и числитель для beta равный (rk,rk) */
Spr1 = 0;
a = 0;
b = 0;
for (i = 0; i < Rows; i++)
{
Xk[i] += alpha * Zk[i];
Rk[i] -= alpha * Sz[i];
Pk[i] -= alpha * St[i];
Spr1 += Pk[i] * Rk[i];
}
/* Вычисляем beta beta = (pk * rk)/(pk-1 * rk-1) (*/
beta = Spr1 / Spr;
/* Вычисляем вектор спуска: zk = rk+ beta * zk-1 */
for (i = 0; i < Rows; i++)
{
Zk[i] = Rk[i] + beta * Zk[i];
Sk[i] = Pk[i] + beta * Sk[i];
}
}
/* Проверяем условие выхода из итерационного цикла */
while (Spr1/mf > eps);
int R;
for (i = 0; i < Rows; i++)
{
R = i + 1;
textBox3.Text += "x[" + (i + 1) + "]= " + Xk[i].ToString() + "\n";
} |