Модераторы: Daevaorn
  

Поиск:

Ответ в темуСоздание новой темы Создание опроса
> Умножение матриц Штрассена и OpenMP 
:(
    Опции темы
DDDENISSS
Дата 8.12.2012, 20:34 (ссылка) | (нет голосов) Загрузка ... Загрузка ... Быстрая цитата Цитата


Новичок



Профиль
Группа: Участник
Сообщений: 4
Регистрация: 16.4.2011

Репутация: нет
Всего: нет



Написал алгоритм штрассена умножения матриц. И оптимизирую его распараллеливанием при помощи 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 секунд
Вот весь проект. Там и билды. Можете проверить их скорость. Только размер матрицы вводите кратный двойке.

Это сообщение отредактировал(а) DDDENISSS - 8.12.2012, 20:40

Присоединённый файл ( Кол-во скачиваний: 28 )
Присоединённый файл  Strassen.rar 80,78 Kb
PM MAIL   Вверх
FCM
Дата 2.1.2013, 19:45 (ссылка) |    (голосов:1) Загрузка ... Загрузка ... Быстрая цитата Цитата


Опытный
**


Профиль
Группа: Участник
Сообщений: 461
Регистрация: 30.3.2009

Репутация: нет
Всего: 9



Сделал несколько прогонов твоего примера без твоих билдов:

mingw-builds 4.7.2.x32   / Code::Blocks Release (O2)/ Windows XP 32 / Intel-E8400 (2 однопоточных  ядра)
Код

N = 2048

strassen time = 16234
strassen omp time = 11563
boost strassen omp = 1.40396

normal time = 156547
normal omp time = 120828
boost normal omp = 1.29562

N = 1024

strassen time = 1953
strassen omp time = 1391
boost strassen omp = 1.40403

normal time = 4265
normal omp time = 3735
boost normal omp = 1.1419

N = 1024

strassen time = 1953
strassen omp time = 1406
boost strassen omp = 1.38905

normal time = 4234
normal omp time = 2875
boost normal omp = 1.4727

N = 512

strassen time = 218
strassen omp time = 157
boost strassen omp = 1.38854

normal time = 500
normal omp time = 250
boost normal omp = 2

N = 256

strassen time = 16
strassen omp time = 31
boost strassen omp = 0.516129

normal time = 31
normal omp time = 16
boost normal omp = 1.9375

N = 256

strassen time = 31
strassen omp time = 16
boost strassen omp = 1.9375

normal time = 31
normal omp time = 31
boost normal omp = 1



Возможны нюансы, связанные с  использованием кэш-памяти.

Это сообщение отредактировал(а) FCM - 3.1.2013, 09:32
PM MAIL   Вверх
  
Ответ в темуСоздание новой темы Создание опроса
Правила форума "С++:Общие вопросы"
Earnest Daevaorn

Добро пожаловать!

  • Черновик стандарта C++ (за октябрь 2005) можно скачать с этого сайта. Прямая ссылка на файл черновика(4.4мб).
  • Черновик стандарта C (за сентябрь 2005) можно скачать с этого сайта. Прямая ссылка на файл черновика (3.4мб).
  • Прежде чем задать вопрос, прочтите это и/или это!
  • Здесь хранится весь мировой запас ссылок на документы, связанные с C++ :)
  • Не брезгуйте пользоваться тегами [code=cpp][/code].
  • Пожалуйста, не просите написать за вас программы в этом разделе - для этого существует "Центр Помощи".
  • C++ FAQ

Если Вам понравилась атмосфера форума, заходите к нам чаще! С уважением, Earnest Daevaorn

 
0 Пользователей читают эту тему (0 Гостей и 0 Скрытых Пользователей)
0 Пользователей:
« Предыдущая тема | C/C++: Общие вопросы | Следующая тема »


 




[ Время генерации скрипта: 0.0485 ]   [ Использовано запросов: 22 ]   [ GZIP включён ]


Реклама на сайте     Информационное спонсорство

 
По вопросам размещения рекламы пишите на vladimir(sobaka)vingrad.ru
Отказ от ответственности     Powered by Invision Power Board(R) 1.3 © 2003  IPS, Inc.