Написал алгоритм штрассена умножения матриц. И оптимизирую его распараллеливанием при помощи OpenMI.
На моем ноутбуке с Core i3 увеличивается скорость в 1.1-1.3 раза. Не знаю мало или нормально. Но основная проблема в том, что этот прирост скорости только в debug режиме, а в release работает обычно вообще медленнее.
В настройках включил openmi для обоих режимов.
И такие строчки
#pragma omp parallel
cout << omp_get_thread_num() << endl;
показывают работу 2х потоков.
Использую MinGW и Netbeans.
Почему такое происходит?
И возможно еще оптимизировать код?
| Код |
#pragma once
#include <omp.h>
#include "iostream"
struct Matrix {
private:
int size;
double* mat;
public:
Matrix(int size) {
this->size = size;
mat = new double[size * size];
}
~Matrix() {
delete[] mat;
}
static void mul(const Matrix& a, const Matrix& b, Matrix& result) { // проверочный алгоритм
const int N = result.getSize();
for (int x = 0; x < N; x++) {
for (int y = 0; y < N; y++) {
double sum = 0;
for (int k = 0; k < N; k++) {
sum += a.get(x, k) * b.get(k, y);
}
result.set(x, y) = sum;
}
}
}
static void mul_omp(const Matrix& a, const Matrix& b, Matrix& result) { // проверочный алгоритм
const int N = result.getSize();
#pragma omp parallel for default(shared)
for (int x = 0; x < N; x++) {
for (int y = 0; y < N; y++) {
double sum = 0;
for (int k = 0; k < N; k++) {
sum += a.get(x, k) * b.get(k, y);
}
result.set(x, y) = sum;
}
}
}
inline double& set(int x, int y) {
return mat[y * size + x];
}
inline double get(int x, int y) const {
return mat[y * size + x];
}
int getSize() const {
return size;
}
};
|
| Код |
#pragma once
#include "iostream"
#include "Matrix.h"
#include <omp.h>
namespace Strassen_OMP {
void strassen(const Matrix& a, const Matrix& b, Matrix& result);
void add(const Matrix& a, const Matrix& b, Matrix& result);
void sub(const Matrix& a, const Matrix& b, Matrix& result);
void strassen(const Matrix& a, const Matrix& b, Matrix& result) {
const int size = result.getSize();
if (size <= 64) {
Matrix::mul_omp(a, b, result);
return;
}
const int halfSize = size / 2;
Matrix a11 = Matrix(halfSize);
Matrix a12 = Matrix(halfSize);
Matrix a21 = Matrix(halfSize);
Matrix a22 = Matrix(halfSize);
Matrix b11 = Matrix(halfSize);
Matrix b12 = Matrix(halfSize);
Matrix b21 = Matrix(halfSize);
Matrix b22 = Matrix(halfSize);
Matrix r11 = Matrix(halfSize);
Matrix r12 = Matrix(halfSize);
Matrix r21 = Matrix(halfSize);
Matrix r22 = Matrix(halfSize);
Matrix p1 = Matrix(halfSize);
Matrix p2 = Matrix(halfSize);
Matrix p3 = Matrix(halfSize);
Matrix p4 = Matrix(halfSize);
Matrix p5 = Matrix(halfSize);
Matrix p6 = Matrix(halfSize);
Matrix p7 = Matrix(halfSize);
// разделяем каждую матрицу на 4 матрицы
#pragma omp parallel for default(shared)
for (int x = 0; x < halfSize; x++) {
for (int y = 0; y < halfSize; y++) {
a11.set(x, y) = a.get(x, y);
a12.set(x, y) = a.get(x, y + halfSize);
a21.set(x, y) = a.get(x + halfSize, y);
a22.set(x, y) = a.get(x + halfSize, y + halfSize);
b11.set(x, y) = b.get(x, y);
b12.set(x, y) = b.get(x, y + halfSize);
b21.set(x, y) = b.get(x + halfSize, y);
b22.set(x, y) = b.get(x + halfSize, y + halfSize);
}
}
// Calculating p1 to p7:
#pragma omp parallel default(shared)
#pragma omp single
{
#pragma omp task
{
Matrix aResult = Matrix(halfSize);
Matrix bResult = Matrix(halfSize);
add(a11, a22, aResult); // a11 + a22
add(b11, b22, bResult); // b11 + b22
strassen(aResult, bResult, p1); // p1 = (a11+a22) * (b11+b22)
}
#pragma omp task
{
Matrix aResult = Matrix(halfSize);
add(a21, a22, aResult); // a21 + a22
strassen(aResult, b11, p2); // p2 = (a21+a22) * (b11)
}
#pragma omp task
{
Matrix bResult = Matrix(halfSize);
sub(b12, b22, bResult); // b12 - b22
strassen(a11, bResult, p3); // p3 = (a11) * (b12 - b22)
}
#pragma omp task
{
Matrix bResult = Matrix(halfSize);
sub(b21, b11, bResult); // b21 - b11
strassen(a22, bResult, p4); // p4 = (a22) * (b21 - b11)
}
#pragma omp task
{
Matrix aResult = Matrix(halfSize);
add(a11, a12, aResult); // a11 + a12
strassen(aResult, b22, p5); // p5 = (a11+a12) * (b22)
}
#pragma omp task
{
Matrix aResult = Matrix(halfSize);
Matrix bResult = Matrix(halfSize);
sub(a21, a11, aResult); // a21 - a11
add(b11, b12, bResult); // b11 + b12
strassen(aResult, bResult, p6); // p6 = (a21-a11) * (b11+b12)
}
#pragma omp task
{
Matrix aResult = Matrix(halfSize);
Matrix bResult = Matrix(halfSize);
sub(a12, a22, aResult); // a12 - a22
add(b21, b22, bResult); // b21 + b22
strassen(aResult, bResult, p7); // p7 = (a12-a22) * (b21+b22)
}
#pragma omp taskwait
#pragma omp task
{
add(p3, p5, r12); // c12 = p3 + p5
add(p2, p4, r21); // c21 = p2 + p4
}
#pragma omp task
{
Matrix aResult = Matrix(halfSize);
Matrix bResult = Matrix(halfSize);
add(p1, p4, aResult); // p1 + p4
add(aResult, p7, bResult); // p1 + p4 + p7
sub(bResult, p5, r11); // c11 = p1 + p4 - p5 + p7
}
#pragma omp task
{
Matrix aResult = Matrix(halfSize);
Matrix bResult = Matrix(halfSize);
add(p1, p3, aResult); // p1 + p3
add(aResult, p6, bResult); // p1 + p3 + p6
sub(bResult, p2, r22); // c22 = p1 + p3 - p2 + p6
}
#pragma omp taskwait
}
// соединяем матрицы результата
#pragma omp parallel for default(shared)
for (int x = 0; x < halfSize; x++) {
for (int y = 0; y < halfSize; y++) {
result.set(x, y) = r11.get(x, y);
result.set(x, y + halfSize) = r12.get(x, y);
result.set(x + halfSize, y) = r21.get(x, y);
result.set(x + halfSize, y + halfSize) = r22.get(x, y);
}
}
}
void add(const Matrix& a, const Matrix& b, Matrix& result) {
const int size = result.getSize();
for (int x = 0; x < size; x++) {
for (int y = 0; y < size; y++) {
result.set(x, y) = a.get(x, y) + b.get(x, y);
}
}
}
void sub(const Matrix& a, const Matrix& b, Matrix& result) {
const int size = result.getSize();
for (int x = 0; x < size; x++) {
for (int y = 0; y < size; y++) {
result.set(x, y) = a.get(x, y) - b.get(x, y);
}
}
}
}
|
Добавлено через 10 минут и 40 секунд
Вот весь проект. Там и билды. Можете проверить их скорость. Только размер матрицы вводите кратный двойке. |