GCC Code Coverage Report

Directory: ./
File: tasks/smyshlaev_a_sle_cg_seq/omp/src/ops_omp.cpp
Date: 2026-04-02 17:12:27
Exec Total Coverage
Lines: 57 92 62.0%
Functions: 8 11 72.7%
Branches: 34 80 42.5%
Line Branch Exec Source
1 #include "smyshlaev_a_sle_cg_seq/omp/include/ops_omp.hpp"
2
3 #include <cmath>
4 #include <cstddef>
5 #include <vector>
6
7 #include "smyshlaev_a_sle_cg_seq/common/include/common.hpp"
8 #include "util/include/util.hpp"
9
10 namespace smyshlaev_a_sle_cg_seq {
11
12 namespace {
13
14 ✗ double ComputeDotProduct(const std::vector<double> &v1, const std::vector<double> &v2, int num_threads) {
15 double result = 0.0;
16 ✗ int n = static_cast<int>(v1.size());
17 ✗ #pragma omp parallel for default(none) shared(n, v1, v2) num_threads(num_threads) reduction(+ : result)
18 for (int i = 0; i < n; ++i) {
19 result += v1[i] * v2[i];
20 }
21 ✗ return result;
22 }
23 void ComputeAp(const std::vector<double> &matrix, const std::vector<double> &p, std::vector<double> &ap, int n,
24 int num_threads) {
25 ✗ #pragma omp parallel for default(none) shared(n, matrix, p, ap) num_threads(num_threads)
26 for (int i = 0; i < n; ++i) {
27 double sum = 0.0;
28 for (int j = 0; j < n; ++j) {
29 sum += matrix[(i * n) + j] * p[j];
30 }
31 ap[i] = sum;
32 }
33 }
34 ✗ double UpdateResultAndResidual(std::vector<double> &result, std::vector<double> &r, const std::vector<double> &p,
35 const std::vector<double> &ap, double alpha, int num_threads) {
36 double rs_new = 0.0;
37 ✗ int n = static_cast<int>(result.size());
38 ✗ #pragma omp parallel for default(none) shared(n, result, r, p, ap, alpha) num_threads(num_threads) reduction(+ : rs_new)
39 for (int i = 0; i < n; ++i) {
40 result[i] += alpha * p[i];
41 r[i] -= alpha * ap[i];
42 rs_new += r[i] * r[i];
43 }
44 ✗ return rs_new;
45 }
46
47 void UpdateP(std::vector<double> &p, const std::vector<double> &r, double beta, int num_threads) {
48 ✗ int n = static_cast<int>(p.size());
49 ✗ #pragma omp parallel for default(none) shared(n, p, r, beta) num_threads(num_threads)
50 for (int i = 0; i < n; ++i) {
51 p[i] = r[i] + (beta * p[i]);
52 }
53 }
54
55 double ComputeDotProduct(const std::vector<double> &v1, const std::vector<double> &v2) {
56 double result = 0.0;
57 96 int n = static_cast<int>(v1.size());
58
2/2
✓ Branch 0 taken 252 times.
✓ Branch 1 taken 96 times.
 348 for (int i = 0; i < n; ++i) {
59 252 result += v1[i] * v2[i];
60 }
61 return result;
62 }
63 96 void ComputeAp(const std::vector<double> &matrix, const std::vector<double> &p, std::vector<double> &ap, int n) {
64
2/2
✓ Branch 0 taken 252 times.
✓ Branch 1 taken 96 times.
 348 for (int i = 0; i < n; ++i) {
65 double sum = 0.0;
66
2/2
✓ Branch 0 taken 708 times.
✓ Branch 1 taken 252 times.
 960 for (int j = 0; j < n; ++j) {
67 708 sum += matrix[(i * n) + j] * p[j];
68 }
69 252 ap[i] = sum;
70 }
71 96 }
72 96 double UpdateResultAndResidual(std::vector<double> &result, std::vector<double> &r, const std::vector<double> &p,
73 const std::vector<double> &ap, double alpha) {
74 double rs_new = 0.0;
75 96 int n = static_cast<int>(result.size());
76
2/2
✓ Branch 0 taken 252 times.
✓ Branch 1 taken 96 times.
 348 for (int i = 0; i < n; ++i) {
77 252 result[i] += alpha * p[i];
78 252 r[i] -= alpha * ap[i];
79 252 rs_new += r[i] * r[i];
80 }
81 96 return rs_new;
82 }
83
84 void UpdateP(std::vector<double> &p, const std::vector<double> &r, double beta) {
85 int n = static_cast<int>(p.size());
86
2/2
✓ Branch 0 taken 120 times.
✓ Branch 1 taken 48 times.
 168 for (int i = 0; i < n; ++i) {
87 120 p[i] = r[i] + (beta * p[i]);
88 }
89 }
90
91 } // namespace
92
93
1/2
✓ Branch 1 taken 48 times.
✗ Branch 2 not taken.
 48 SmyshlaevASleCgTaskOMP::SmyshlaevASleCgTaskOMP(const InType &in) {
94 SetTypeOfTask(GetStaticTypeOfTask());
95 GetInput() = in;
96 48 }
97
98
1/2
✓ Branch 0 taken 48 times.
✗ Branch 1 not taken.
 48 bool SmyshlaevASleCgTaskOMP::ValidationImpl() {
99 const auto &a = GetInput().A;
100 const auto &b = GetInput().b;
101
2/4
✓ Branch 0 taken 48 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 48 times.
✗ Branch 3 not taken.
 48 if (a.empty() || b.empty()) {
102 return false;
103 }
104
1/2
✓ Branch 0 taken 48 times.
✗ Branch 1 not taken.
 48 if (a.size() != b.size()) {
105 return false;
106 }
107
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if (a.size() != a[0].size()) {
108 ✗ return false;
109 }
110 return true;
111 }
112
113 48 bool SmyshlaevASleCgTaskOMP::PreProcessingImpl() {
114 const auto &a = GetInput().A;
115 size_t n = a.size();
116 48 flat_A_.resize(n * n);
117
2/2
✓ Branch 0 taken 132 times.
✓ Branch 1 taken 48 times.
 180 for (size_t i = 0; i < n; ++i) {
118
2/2
✓ Branch 0 taken 396 times.
✓ Branch 1 taken 132 times.
 528 for (size_t j = 0; j < n; ++j) {
119 396 flat_A_[(i * n) + j] = a[i][j];
120 }
121 }
122 48 return true;
123 }
124
125 48 bool SmyshlaevASleCgTaskOMP::RunSequential() {
126 48 const auto &b = GetInput().b;
127 48 int n = static_cast<int>(b.size());
128 48 std::vector<double> r = b;
129
1/2
✓ Branch 1 taken 48 times.
✗ Branch 2 not taken.
 48 std::vector<double> p = r;
130
1/4
✓ Branch 1 taken 48 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 48 std::vector<double> ap(n, 0.0);
131
1/4
✓ Branch 1 taken 48 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 48 std::vector<double> result(n, 0.0);
132
133 double rs_old = 0.0;
134
2/2
✓ Branch 0 taken 132 times.
✓ Branch 1 taken 48 times.
 180 for (int i = 0; i < n; ++i) {
135 132 rs_old += r[i] * r[i];
136 }
137
138 const double epsilon = 1e-9;
139
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if (std::sqrt(rs_old) < epsilon) {
140 ✗ GetOutput() = result;
141 return true;
142 }
143
144 48 const int max_iterations = n * 2;
145
1/2
✓ Branch 0 taken 96 times.
✗ Branch 1 not taken.
 96 for (int iter = 0; iter < max_iterations; ++iter) {
146 96 ComputeAp(flat_A_, p, ap, n);
147 double p_ap = ComputeDotProduct(p, ap);
148
1/2
✓ Branch 0 taken 96 times.
✗ Branch 1 not taken.
 96 if (std::abs(p_ap) < 1e-15) {
149 break;
150 }
151
152 96 double alpha = rs_old / p_ap;
153 96 double rs_new = UpdateResultAndResidual(result, r, p, ap, alpha);
154
2/2
✓ Branch 0 taken 48 times.
✓ Branch 1 taken 48 times.
 96 if (std::sqrt(rs_new) < epsilon) {
155 break;
156 }
157
158 48 double beta = rs_new / rs_old;
159 UpdateP(p, r, beta);
160 rs_old = rs_new;
161 }
162
163
1/2
✓ Branch 1 taken 48 times.
✗ Branch 2 not taken.
 48 GetOutput() = result;
164 return true;
165 }
166
167 ✗ bool SmyshlaevASleCgTaskOMP::RunParallel(int num_threads) {
168 ✗ const auto &b = GetInput().b;
169 ✗ int n = static_cast<int>(b.size());
170 ✗ std::vector<double> r = b;
171 ✗ std::vector<double> p = r;
172 ✗ std::vector<double> ap(n, 0.0);
173 ✗ std::vector<double> result(n, 0.0);
174
175 double rs_old = 0.0;
176
177 ✗ #pragma omp parallel for default(none) shared(n, r) num_threads(num_threads) reduction(+ : rs_old)
178 for (int i = 0; i < n; ++i) {
179 rs_old += r[i] * r[i];
180 }
181
182 const double epsilon = 1e-9;
183 ✗ if (std::sqrt(rs_old) < epsilon) {
184 ✗ GetOutput() = result;
185 return true;
186 }
187
188 ✗ const int max_iterations = n * 2;
189 ✗ for (int iter = 0; iter < max_iterations; ++iter) {
190 ✗ ComputeAp(flat_A_, p, ap, n, num_threads);
191 ✗ double p_ap = ComputeDotProduct(p, ap, num_threads);
192 ✗ if (std::abs(p_ap) < 1e-15) {
193 break;
194 }
195
196 ✗ double alpha = rs_old / p_ap;
197 ✗ double rs_new = UpdateResultAndResidual(result, r, p, ap, alpha, num_threads);
198 ✗ if (std::sqrt(rs_new) < epsilon) {
199 break;
200 }
201
202 ✗ double beta = rs_new / rs_old;
203 UpdateP(p, r, beta, num_threads);
204 rs_old = rs_new;
205 }
206
207 ✗ GetOutput() = result;
208 return true;
209 }
210
211
1/2
✓ Branch 0 taken 48 times.
✗ Branch 1 not taken.
 48 bool SmyshlaevASleCgTaskOMP::RunImpl() {
212 const auto &b = GetInput().b;
213
1/2
✓ Branch 0 taken 48 times.
✗ Branch 1 not taken.
 48 if (b.empty()) {
214 return true;
215 }
216 48 int n = static_cast<int>(b.size());
217
218 48 int num_threads = ppc::util::GetNumThreads();
219
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if (n < num_threads * 100) {
220 num_threads = 1;
221 }
222 ✗ if (num_threads == 1) {
223 48 return RunSequential();
224 }
225
226 ✗ return RunParallel(num_threads);
227 }
228
229 48 bool SmyshlaevASleCgTaskOMP::PostProcessingImpl() {
230 48 return true;
231 }
232
233 } // namespace smyshlaev_a_sle_cg_seq
234