3
3
#include < cmath>
4
4
#include < iostream>
5
5
#include < vector>
6
+ #include < iomanip>
6
7
7
8
void gaussian_elimination (std::vector<double >& a, int cols) {
8
9
assert (a.size () % cols == 0 );
@@ -13,16 +14,12 @@ void gaussian_elimination(std::vector<double>& a, int cols) {
13
14
// Main loop going through all columns
14
15
for (int col = 0 ; col < cols - 1 ; ++col) {
15
16
// Step 1: finding the maximum element for each column
16
- int max_index = [&]() {
17
- int res = row;
18
- int max_element = a[col + row * cols];
19
- for (int r = row + 1 ; r < rows; ++r)
20
- if (max_element < std::abs (a[col + r * cols])) {
21
- max_element = std::abs (a[col + r * cols]);
22
- res = r;
23
- }
24
- return res;
25
- }();
17
+ int max_index = row;
18
+ for (int r = row + 1 ; r < rows; ++r) {
19
+ if (a[col + max_index * cols] < std::abs (a[col + r * cols])) {
20
+ max_index = r;
21
+ }
22
+ }
26
23
27
24
// Check to make sure matrix is good!
28
25
if (a[col + max_index * cols] == 0 ) {
@@ -61,10 +58,6 @@ std::vector<double> back_substitution(const std::vector<double>& a, int cols) {
61
58
// Creating the solution Vector
62
59
std::vector<double > soln (rows);
63
60
64
- // initialize the final element
65
- soln[rows - 1 ] =
66
- a[cols - 1 + (rows - 1 ) * cols] / a[cols - 1 - 1 + (rows - 1 ) * cols];
67
-
68
61
for (int i = rows - 1 ; i >= 0 ; --i) {
69
62
auto sum = 0.0 ;
70
63
for (int j = cols - 2 ; j > i; --j) {
@@ -86,12 +79,12 @@ void gauss_jordan_elimination(std::vector<double>& a, int cols) {
86
79
if (a[col + row * cols] != 0 ) {
87
80
88
81
// divide row by pivot and leaving pivot as 1
89
- for (int i = cols - 1 ; i >= static_cast < int >( col) ; --i)
82
+ for (int i = cols - 1 ; i >= col; --i)
90
83
a[i + row * cols] /= a[col + row * cols];
91
84
92
85
// subtract value from above row and set values above pivot to 0
93
- for (int i = 0 ; i < static_cast < int >( row - 1 ) ; ++i)
94
- for (int j = cols - 1 ; j >= static_cast < int >( col) ; --j)
86
+ for (int i = 0 ; i < row; ++i)
87
+ for (int j = cols - 1 ; j >= col; --j)
95
88
a[j + i * cols] -= a[col + i * cols] * a[j + row * cols];
96
89
++row;
97
90
}
@@ -104,16 +97,22 @@ void print_matrix(const std::vector<double>& a, int cols) {
104
97
for (int i = 0 ; i < rows; ++i) {
105
98
std::cout << " [" ;
106
99
for (int j = 0 ; j < cols; ++j) {
107
- std::cout << a[j + i * cols] << " " ;
100
+ std::cout << std::fixed << a[j + i * cols] << " \t " ;
108
101
}
109
102
std::cout << " ]\n " ;
110
103
}
111
104
}
112
105
113
106
int main () {
114
- std::vector<double > a = {2 , 3 , 4 , 6 , 1 , 2 , 3 , 4 , 3 , -4 , 0 , 10 };
107
+ std::vector<double > a = { 2 , 3 , 4 , 6 ,
108
+ 1 , 2 , 3 , 4 ,
109
+ 3 , -4 , 0 , 10 };
115
110
const int cols = 4 ;
116
- assert (a.size () % cols == 0 );
111
+ if (a.size () % cols != 0 )
112
+ {
113
+ std::cout << " The input dimentions are incorrect\n " ;
114
+ return 1 ;
115
+ }
117
116
118
117
gaussian_elimination (a, cols);
119
118
print_matrix (a, cols);
0 commit comments