|
| 1 | +#include <cstddef> |
1 | 2 | #include <iostream>
|
2 | 3 | #include <vector>
|
3 |
| -#include <cstring> |
4 | 4 |
|
5 |
| -void thomas(std::vector<double> const a, std::vector<double> const b, std::vector<double> const c, std::vector<double>& x) { |
6 |
| - int size = a.size(); |
7 |
| - double y[size]; |
8 |
| - memset(y, 0, size * sizeof(double)); |
9 |
| - |
10 |
| - y[0] = c[0] / b[0]; |
11 |
| - x[0] = x[0] / b[0]; |
12 |
| - |
13 |
| - for (size_t i = 1; i < size; ++i) { |
14 |
| - double scale = 1.0 / (b[i] - a[i] * y[i - 1]); |
15 |
| - y[i] = c[i] * scale; |
16 |
| - x[i] = (x[i] - a[i] * x[i - 1]) * scale; |
17 |
| - } |
18 |
| - |
19 |
| - for (int i = size - 2; i >= 0; --i) { |
20 |
| - x[i] -= y[i] * x[i + 1]; |
21 |
| - } |
| 5 | +void thomas( |
| 6 | + std::vector<double> const& a, |
| 7 | + std::vector<double> const& b, |
| 8 | + std::vector<double> const& c, |
| 9 | + std::vector<double>& x) { |
| 10 | + auto y = std::vector<double>(a.size(), 0.0); |
| 11 | + |
| 12 | + y[0] = c[0] / b[0]; |
| 13 | + x[0] = x[0] / b[0]; |
| 14 | + |
| 15 | + for (std::size_t i = 1; i < a.size(); ++i) { |
| 16 | + const auto scale = 1.0 / (b[i] - a[i] * y[i - 1]); |
| 17 | + y[i] = c[i] * scale; |
| 18 | + x[i] = (x[i] - a[i] * x[i - 1]) * scale; |
| 19 | + } |
| 20 | + |
| 21 | + for (std::size_t i = a.size() - 2; i < a.size(); --i) { |
| 22 | + x[i] -= y[i] * x[i + 1]; |
| 23 | + } |
22 | 24 | }
|
23 | 25 |
|
24 | 26 | int main() {
|
25 |
| - std::vector<double> a = {0.0, 2.0, 3.0}; |
26 |
| - std::vector<double> b = {1.0, 3.0, 6.0}; |
27 |
| - std::vector<double> c = {4.0, 5.0, 0.0}; |
28 |
| - std::vector<double> x = {7.0, 5.0, 3.0}; |
| 27 | + const std::vector<double> a = {0.0, 2.0, 3.0}; |
| 28 | + const std::vector<double> b = {1.0, 3.0, 6.0}; |
| 29 | + const std::vector<double> c = {4.0, 5.0, 0.0}; |
| 30 | + std::vector<double> x = {7.0, 5.0, 3.0}; |
29 | 31 |
|
30 |
| - std::cout << "The system" << std::endl; |
31 |
| - std::cout << "[1.0 4.0 0.0][x] = [7.0]" << std::endl; |
32 |
| - std::cout << "[2.0 3.0 5.0][y] = [5.0]" << std::endl; |
33 |
| - std::cout << "[0.0 3.0 6.0][z] = [3.0]" << std::endl; |
34 |
| - std::cout << "has the solution" << std::endl; |
| 32 | + std::cout << "The system\n"; |
| 33 | + std::cout << "[1.0 4.0 0.0][x] = [7.0]\n"; |
| 34 | + std::cout << "[2.0 3.0 5.0][y] = [5.0]\n"; |
| 35 | + std::cout << "[0.0 3.0 6.0][z] = [3.0]\n"; |
| 36 | + std::cout << "has the solution:\n"; |
35 | 37 |
|
36 |
| - thomas(a, b, c, x); |
| 38 | + thomas(a, b, c, x); |
37 | 39 |
|
38 |
| - for (size_t i = 0; i < 3; ++i) { |
39 |
| - std::cout << "[" << x[i] << "]" << std::endl; |
40 |
| - } |
| 40 | + for (auto const& val : x) { |
| 41 | + std::cout << "[" << val << "]\n"; |
| 42 | + } |
41 | 43 |
|
42 |
| - return 0; |
| 44 | + return 0; |
43 | 45 | }
|
0 commit comments