Форум программистов, компьютерный форум, киберфорум
igorrr37
Войти
Регистрация
Восстановить пароль
Рейтинг: 4.00. Голосов: 1.

Решение СЛАУ методом Гаусса выбором главного элемента по столбцу

Запись от igorrr37 размещена 25.04.2013 в 09:54
Обновил(-а) igorrr37 11.10.2014 в 11:35

C++
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
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
/*
Матрица читается из файла in.txt.
Число столбцов матрицы должно быть на 1 больше чем число строк.
Максимально допустимая размерность матрицы 20х21.
Пример файла in.txt с матрицей 2х3:
 
23   12.2    4
-3   -13     18.2
 
*/
 
#include <iostream>
#include <clocale>
#include <cstdlib>
#include <iomanip>
#include <cstring>
#include <fstream>
#define max_rows 20
#define max_columns (max_rows + 1)
 
// Освобождает память из-под матрицы
void FreeMemory(double** pmatrix, size_t const rows)
{
    for(size_t i = 0; i < rows; ++i)
    {
        delete[] pmatrix[i];
        pmatrix[i] = 0;
    }
    delete[] pmatrix;
    pmatrix = 0;
}
 
// Выделение памяти под матрицу. Возвращает указатель на выделенную память.
double** GetMemory(double**& pmatrix, size_t const rows, size_t const columns)
{
    size_t i = 0;
    pmatrix = new double*[rows];
    for(i = 0; i < rows; ++i)
        pmatrix[i] = 0;
    for(i = 0; i < rows; ++i)
        pmatrix[i] = new double[columns];
    return pmatrix;
}
 
// Ввод матрицы из файла. Возвращает указатель на матрицу в памяти.
double** FillFromFile(char const* filepath, double**& pmatrix, size_t& rows, size_t& columns)
{
    std::ifstream ifs(filepath);
    if(ifs.is_open())
    {
        // чтение матрицы из файла в массив nums
        double nums[max_rows * max_columns];
        int siz;
        for(siz = 0; ifs >> nums[siz]; ++siz)
            ;
        // вычисление размерности считанной матрицы (rows x columns)
        for(rows = 1; rows <= max_rows; ++rows)
        {
            if(siz / rows == rows + 1)
            {
                columns = rows + 1;
                break;
            }
        }
        if(rows != max_columns)
        {
            GetMemory(pmatrix, rows, columns);
            for(size_t i = 0; i < rows; ++i)
            {
                for(size_t j = 0; j < columns; ++j)
                {
                    pmatrix[i][j] = nums[i * columns + j];
                }
            }
        }
        else
            std::cerr << "FillFromFile failed: fail to calculate matrix dimension (rows and columns unknown)\n";
        ifs.close();
    }
    else
        std::cerr << "FillFromFile failed: unable to open input file\n";
    return pmatrix;
}
 
// Вывод матрицы в поток (по-умолчанию - на экран)
std::ostream& MatrixToStream(double* const* const pmatrix, size_t const rows, size_t const columns, std::ostream& ost = std::cout)
{
    ost << "\n-------------------------MatrixToStream-------------------------------\n";
    ost << "Matrix size: " << rows << " x " << columns << "\n";
    for(size_t i = 0; i < rows; ++i)
    {
        for(size_t j = 0; j < columns; ++j)
        {
            ost << std::setw(13) << std::left << pmatrix[i][j];
        }
        ost << "\n\n";
    }
    ost << "\n----------------------------------------------------------------------\n";
    return ost;
}
 
// Умножение массива на число
double* MultiplicationArrayToNumber(double* parray, size_t const columns, double const number)
{
    for(size_t j = 0; j < columns; ++j)
    {
        parray[j] = parray[j] * number;
    }
    return parray;
}
 
// Вычитание массива №2 из массива №1 (они равной длины)
double* SubtractionOfArrays(double* parray1, double* parray2, size_t const columns)
{
    for(size_t j = 0; j < columns; ++j)
    {
        parray1[j] = parray1[j] - parray2[j];
    }
    return parray1;
}
 
// Обмен местами двух массивов в памяти
void SwapArrays(double* parray1, double* parray2, size_t const columns)
{
    double* ptmparray = new double[columns];
    memcpy(ptmparray, parray1, columns * sizeof(double));
    memcpy(parray1, parray2, columns * sizeof(double));
    memcpy(parray2, ptmparray, columns * sizeof(double));
    delete[] ptmparray;
    ptmparray = 0;
}
 
// Поиск строки среди текущей и нижележащих строк, содержащей в заданном столбце максимальный элемент.
//Возвращает индекс найденной строки.
size_t MaxRowIndex(double** pdmatrix, size_t const rows, size_t const columnindex)
{
    size_t maxrowindex = columnindex, i, j;
    for(i = columnindex + 1; i < rows; ++i)
    {
        if(pdmatrix[i][columnindex] > pdmatrix[columnindex][columnindex])
            maxrowindex = i;
    }
    return maxrowindex;
}
 
// Приведение матрицы к верхнетреугольному виду.
 // Возвращает число корней матрицы.
size_t MakeTriangularMatrix(double** pdmatrix, size_t const rows, size_t const columns)
{
    size_t rootcount = rows;
    size_t i = 0, j = 0, k = 0, maxrowindex = 0;
    double* prowcopy = new double[columns];
    for(j = 0; j < rows - 1; ++j)
    {
        maxrowindex = MaxRowIndex(pdmatrix, rows, j);
        if(maxrowindex != j)
            SwapArrays(pdmatrix[maxrowindex], pdmatrix[j], columns);
        if(double(0) != pdmatrix[j][j])
        {
            MultiplicationArrayToNumber(pdmatrix[j], columns, 1 / pdmatrix[j][j]);
            for(i = j + 1; i < rows; ++i)
            {
                memcpy(prowcopy, pdmatrix[j], columns * sizeof(double));
                MultiplicationArrayToNumber(prowcopy, columns, pdmatrix[i][j]);
                SubtractionOfArrays(pdmatrix[i], prowcopy, columns);
            }
        }
    }
    delete[] prowcopy;
    prowcopy = 0;
    // Проверка системы на вырожденность
    for(i = 0; i < rows; ++i)
    {
        bool zerorow = true;
        for(j = 0; j < columns - 1; ++j)
        {
            if(pdmatrix[i][j] != double(0))
            {
                zerorow = false;
                break;
            }
        }
        if(zerorow && pdmatrix[i][columns - 1] != double(0))
        {
            std::cout << "\nCистема НЕСОВМЕСТНА\n\n";
            rootcount = 0;
            break;
        }
        else if(zerorow && pdmatrix[i][columns - 1] == double(0))
        {
            std::cout << "\nCистема имеет БЕСКОНЕЧНОЕ МНОЖЕСТВО РЕШЕНИЙ\n\n";
            rootcount = rows + 5;
            break;
        }
    }
    return rootcount;
}
 
//Вычисляет корни СЛАУ по верхнетреугольной матрице.
//Возвращает указатель на массив корней.
double* GetRoots(double** pdmatrix, size_t const rows, double* pdroots)
{
    int i = 0, j = 0;
    for(i = rows - 1; i >= 0; --i)
    {
        double sum = 0;
        for(j = rows - 1; j > i; --j)
            sum += pdroots[j] * pdmatrix[i][j];
        pdroots[i] = (pdmatrix[i][rows] - sum) / pdmatrix[i][i];
    }
    return pdroots;
}
 
int main()
{
    setlocale(LC_ALL, "rus");
    char const* filepath = "in.txt"; // путь к файлу с матрицей
    size_t rows = 0, columns = 0, i = 0, j = 0;
    double** pdmatrix = 0;
    FillFromFile(filepath, pdmatrix, rows, columns);
    std::cout << "\n\nИсходная матрица:\n";
    MatrixToStream(pdmatrix, rows, columns);
    size_t rootcount = MakeTriangularMatrix(pdmatrix, rows, columns);
    std::cout << "\n\nТреугольная матрица:\n";
    MatrixToStream(pdmatrix, rows, columns);
    if(rows == rootcount)
    {
        double* pdroots = new double[rows];
        GetRoots(pdmatrix, rows, pdroots);
        std::cout << "\n\nКорни системы:\n\n";
        for(i = 0; i < rows; ++i)
            std::cout << "X" << i << ": " << pdroots[i] << '\n';
        delete[] pdroots;
        pdroots = 0;
    }
    else
        std::cout << "Матрица вырождена\n";
    FreeMemory(pdmatrix, rows);
    return 0;
}
Матрица Гильберта десятого порядка для проверки программы
Bash
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
1.000000000000000   0.500000000000000   0.333333333333330   0.250000000000000   0.200000000000000   0.166666666666670   0.142857142857140   0.125000000000000   0.111111111111110   0.100000000000000   2.928968253968250   
0.500000000000000   0.333333333333330   0.250000000000000   0.200000000000000   0.166666666666670   0.142857142857140   0.125000000000000   0.111111111111110   0.100000000000000   0.090909090909090   2.019877344877340
0.333333333333330   0.250000000000000   0.200000000000000   0.166666666666670   0.142857142857140   0.125000000000000   0.111111111111110   0.100000000000000   0.090909090909090   0.083333333333330   1.603210678210670
0.250000000000000   0.200000000000000   0.166666666666670   0.142857142857140   0.125000000000000   0.111111111111110   0.100000000000000   0.090909090909090   0.083333333333330   0.076923076923080   1.346800421800420
0.200000000000000   0.166666666666670   0.142857142857140   0.125000000000000   0.111111111111110   0.100000000000000   0.090909090909090   0.083333333333330   0.076923076923080   0.071428571428570   1.168228993228990
0.166666666666670   0.142857142857140   0.125000000000000   0.111111111111110   0.100000000000000   0.090909090909090   0.083333333333330   0.076923076923080   0.071428571428570   0.066666666666670   1.034895659895660
0.142857142857140   0.125000000000000   0.111111111111110   0.100000000000000   0.090909090909090   0.083333333333330   0.076923076923080   0.071428571428570   0.066666666666670   0.062500000000000   0.930728993228990
0.125000000000000   0.111111111111110   0.100000000000000   0.090909090909090   0.083333333333330   0.076923076923080   0.071428571428570   0.066666666666670   0.062500000000000   0.058823529411770   0.846695379783620
0.111111111111110   0.100000000000000   0.090909090909090   0.083333333333330   0.076923076923080   0.071428571428570   0.066666666666670   0.062500000000000   0.058823529411770   0.055555555555560   0.777250935339180
0.100000000000000   0.090909090909090   0.083333333333330   0.076923076923080   0.071428571428570   0.066666666666670   0.062500000000000   0.058823529411770   0.055555555555560   0.052631578947370   0.718771403175440
 
Эталонное решение x
i=1 Xi = 1.000000000000000
i=2 Xi = 1.000000000000000
i=3 Xi = 1.000000000000000
i=4 Xi = 1.000000000000000
i=5 Xi = 1.000000000000000
i=6 Xi = 1.000000000000000
i=7 Xi = 1.000000000000000
i=8 Xi = 1.000000000000000
i=9 Xi = 1.000000000000000
i=10    Xi = 1.000000000000000
Размещено в Без категории
Показов 4586 Комментарии 0
Всего комментариев 0
Комментарии
 
КиберФорум - форум программистов, компьютерный форум, программирование
Powered by vBulletin
Copyright ©2000 - 2024, CyberForum.ru