OpenMP — Распараллеливание метода Крамера решения СЛАУ

      Комментарии к записи OpenMP — Распараллеливание метода Крамера решения СЛАУ отключены

Главная Форумы Программирование Параллельное программирование Использование OpenMP OpenMP — Распараллеливание метода Крамера решения СЛАУ

В этой теме 0 ответов, 1 участник, последнее обновление  Васильев Владимир Сергеевич 2 нед., 5 дн. назад.

  • Автор
    Сообщения
  • #3885

    Ранее мы рассматривали алгоритм Крамера решения СЛАУ, а также его реализацию на языке С++. Если коротко, то алгоритм сводится к вычислению главного определителя матрицы (Δ), а затем N дополнительных определителей Δi (для матриц, полученных заменой i-того столбца на столбец свободных членов). После проделанных вычислений i-тый корень можно получить по формуле xi = Δi/Δ.

    Для оценки времени работы программ на матрицах разного размера я написал следующий код:

    #include "matrix.h"
    #include <iostream>
    #include <vector>
    #include "omp.h"
    
    int main() {
      for (int n = 50; n < 350; n += 50) {
        std::vector<std::vector<double> > matrix(n);
    
        for (int i = 0; i < n; ++i) {
          matrix[i].resize(n);
          for (int j = 0; j < n; ++j)
            matrix[i][j] = rand();
        }
    
        std::vector<double> column(n);
        for (int j = 0; j < n; ++j)
          column[j] = rand();
    
        double start_time = omp_get_wtime();
    
        std::vector<double> solution = cramer_solving(matrix, column, n);
    
        std::cout << n << " " << omp_get_wtime() - start_time << std::endl;
      }
    }

    Что в этом алгоритме может быть выполнено параллельно?
    У нас вычисляется N определителей и именно в этой части алгоритма программа находится большую часть времени — очевидно, нет смысла распараллеливать вычисление xi для i=0..(N-1), имеющее линейную асимптотическую сложность. Очевидно, можно распараллелить алгоритм вычисления определителя, однако он сводится к сведению матрицы к треугольному виду и вычислению произведения элементов диагонали — т.е. целесообразно параллельно выполнять только триангуляцию. В статье «Распараллеливание алгоритма триангуляции матрицы с OpenMP» рассмотрено три варианта параллельных алгоритмов, возьмем лучший из них.

    Блок-схема алгоритма триангуляции (часть алгоритма, которую целесообразно распараллеливать):

    Другой подход к распараллеливанию алгоритма может заключаться в независимом вычислении Δi в разных потоках. Во время вычисления потоки модифицируют матрицу, переставляя в ней i-тый столбец со столбцом свободных членов, поэтому у каждого потока должна быть своя локальная копия матрицы:

    template <typename T>
    std::vector<T> cramer_solving(std::vector<std::vector<T> > &matrix,
                                  std::vector<T> &free_term_column, int n) {
      T mainDet = gauss_determinant(matrix, n);
      if (std::abs(mainDet) < 0.0001)
        throw NoSolution();
      std::vector<T> solution(n);
      #pragma omp parallel
      {
        std::vector<std::vector<T> > private_matrix = matrix;
      #pragma omp for
        for (int i = 0; i < n; ++i) {
          swap_column(private_matrix, free_term_column, i, n);
          solution[i] = gauss_determinant(private_matrix, n) / mainDet;
          swap_column(private_matrix, free_term_column, i, n);
        }
      }
      return solution;
    }

    Результаты работы обоих вариантов параллельных программ, а также последовательной программы сведены в таблицу:

    размерность50100150200250300
    последовательная программа0.02713780.402571.994686.2426515.230632.3143
    распараллеливание cramer_solving0.02270530.2140061.056493.269727.9462516.9398
    параллельный triangulation по столбцам0.02359080.2317681.089093.465448.409316.9586

    Вычисления производились на двухъядерном компьютере: Intel® Pentium(R) CPU G2130 @ 3.20GHz × 2 с 5.8 Гб оперативной памяти и операционной системой openSUSE Leap 42.1 (x86_64) 64-бит.

    Видно, что при маленьком размере матрицы наиболее эффективен последовательный алгоритм. Это связано с издержками на создание потока и локальных переменных. При больших размерах системы уравнений оба параллельных алгоритма показывают схожие результаты с ускорением близким к 90% (это неплохо для двух ядер). Этот пример хорошо показывает, что в параллельном программировании не всегда надо писать новый параллельный алгоритм (как мы сделали переписав cramer_solving) — иногда схожего эффекта можно получить проанализировав структуру алгоритма и подставив параллельную версию используемой в нем функции (это бывает значительно проще).

Для ответа в этой теме необходимо авторизоваться.