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
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
| long double *gauss(int size, long double **main_mat, long double *main_b)
{
int i, j, k;
long double **mat = NULL;
long double *b = NULL;
long double *x = NULL;
int *order = NULL;
//делаем копии переданных данных
if(!(mat = malloc(sizeof(long double*) * size))
|| !(b = malloc(sizeof(long double) * size))
|| !(x = malloc(sizeof(long double) * size))
|| !(order = malloc(sizeof(int) * size)))
{
free(mat);
free(b);
free(x);
free(order);
puts("!");
return NULL;
}
for(i = 0; i < size; i++)
{
if(!(mat[i] = malloc(sizeof(long double) * size)))
{
for(; i >=0; i--)
free(mat[i]);
puts("!");
return NULL;
}
memcpy(mat[i], main_mat[i], sizeof(long double) * size);
order[i] = i;
}
memcpy(b, main_b, sizeof(long double) * size);
puts("\nМетод Гаусса");
for(k = 0; k < size; k++)
{
putchar('\n');
long double max = fabsl(mat[k][k]);
int max_i = k, max_j = k;
if(k != size - 1)
{
//находим максимальный каэффициент
for(i = k; i < size; i++)
{
for(j = k; j < size; j++)
{
if(fabsl(mat[i][order[j]]) > max)
{
max = fabsl(mat[i][order[j]]);
max_i = i;
max_j = j;
}
}
}
//он нулевой?..
if(max == 0)
{
puts("Система ворожденная...");
for(i = 0; i < size; i++)
{
free(mat[i]);
}
free(mat);
free(b);
free(x);
free(order);
return NULL;
}
printf("Оптимальное основание [%d][%d] = %Lf\n",
max_i, max_j, mat[max_i][order[max_j]]);
}
else
puts("Итоговая матрица");
//по не обходимости меняем местами строки
if(max_i != k)
{
long double *ptr = mat[k];
mat[k] = mat[max_i];
mat[max_i] = ptr;
long double tmp = b[k];
b[k] = b[max_i];
b[max_i] = tmp;
}
//и/или переставляем столбцы
if(max_j != k)
{
int tmp = order[k];
order[k] = order[max_j];
order[max_j] = tmp;
}
//делим на основание
for(j = k + 1; j < size; j++)
mat[k][order[j]] /= mat[k][order[k]];
b[k] /= mat[k][order[k]];
mat[k][order[k]] = 1;
//вычитаем из оставшихся строк получившуюя с соответствующим коэффициентом
for(i = k + 1; i < size; i++)
{
long double factor = mat[i][order[k]];
for(j = k; j < size; j++)
mat[i][order[j]] -= mat[k][order[j]] * factor;
b[i] -= b[k] * factor;
}
for(i = 0; i < size; i++)
printf(" x%-7d ", order[i]);
puts(" b");
for(i = 0; i < size; i++)
{
for(j = 0; j < size; j++)
printf("% 8Lf ", mat[i][order[j]]);
printf("% 8Lf\n", b[i]);
}
}
//совстевенно, считаем искомый вектор
for(i = size - 1; i >= 0; i--)
{
for(j = i + 1; j < size; j++)
b[i] -= x[order[j]] * mat[i][order[j]];
x[order[i]] = b[i];
}
for(i = 0; i < size; i++)
{
free(mat[i]);
}
free(mat);
free(b);
free(order);
return x;
} |
(
Сообщение было отмечено Памирыч как решение
