Понимание рекурсивного алгоритма для метода Штрассена - PullRequest
0 голосов
/ 09 июня 2018

Итак, я пытаюсь выяснить метод Штрассена для умножения матриц, я использую C ++, но это может быть любой язык.Прямо сейчас это выглядит так:

 typedef vector<long int> ROW;
 typedef vector<ROW> MATRIX;

 void getQuad(const MATRIX& IN, MATRIX& OUT0, MATRIX& OUT1
                    MATRIX& OUT2, MATRIX& OUT3)
 { /*determine quadrants*/ }

 void strassen(const MATRIX& A, const MATRIX& B, MATRIX& C
 {
      if (A.size() == 2 && A[0] == 2) //know that its 2x2, stop
      {
           // Get M1-M7 vars and set MATRIX C with them
      }
      else
      {
           /*
             getQuad(...) returns the quadrants
             ___________
             | X0 | X1 |
             -----------
             | X2 | X3 |
             -----------
           */

        MATRIX A0,A1,A2,A3;
        getQuad(A,A0,A1,A2,A3);

        MATRIX B0,B1,B2,B3;
        getQuad(B,B0,B1,B2,B3);
      }
 }

Я не уверен, куда идти дальше с отдельными квадрантами, то есть, как получить матрицы M1-M7 в этой точке.Я предположил бы, что матрицы M1-M7 (в отличие от примитивных типов данных в базовом случае) будут использоваться так же, как базовый вариант.Я просто не уверен, как будет выглядеть распутывание здесь.

Я знаю, что немного трудно читать чужой код, но, надеюсь, он был понятен.

Я уверен, что мой базовый случай верен, и я уверен, что я разделяюМатрица правильно, мне просто непонятно, куда идти дальше.Возможно, я неправильно написал свой алгоритм.

1 Ответ

0 голосов
/ 09 июня 2018

Полагаю, вы упустили основную мысль алгоритма Штрассена - тот факт, что он рекурсивен псевдокоде алгоритм будет выглядеть примерно так:

MATRIX strassen(const MATRIX&a, const MATRIX&b) {
    int aw = a.width();
    int ah = a.height();
    int bw = b.width();
    int bh = b.height();

    if (aw != bh)
        throw some_exception();

    // Strassen algorithm requires each size to be a power of 2
    int max_size = max(aw, ah, bw);
    int extended_size = next_pow_2(max_size);
    MATRIX aEx = a.extend(extended_size, extended_size);
    MATRIX bEx = a.extend(extended_size, extended_size);
    MATRIX cEx = strassenImpl(aEx, bEx);

    // truncate back from power of 2 to real size
    return cEx.truncate(ah, bw);
}


MATRIX strassenImpl(const MATRIX&A, const MATRIX&B) {
    // if matrix size is relatively small it is faster to do the usual straightforward multiplication
    if (A.size() <= threshold) {
        return usualMultiply(A, B);
    }
    // alternatively threshold is 1 so matrix multiplication is just multiplication of the single values
    //if (A.size() == 1) {
    //    return MATRIX(A[0][0]*B[0][0]);
    //} 
    else {
        MATRIX A11, A12, A21, A22;
        getQuad(A, A11, A12, A21, A22);

        MATRIX B11, B12, B21, B22;
        getQuad(B, B11, B12, B21, B22);

        // recursive calls, note that we don't need to go through the extension step
        // here because if the size is a power of 2, half of the size is also a power of 2
        MATRIX M1 = strassenImpl(A11 + A22, B11 + B22);
        MATRIX M2 = strassenImpl(A21 + A22, B11);
        MATRIX M3 = strassenImpl(A11, B12 - B22);
        MATRIX M4 = strassenImpl(A22, B21 - B11);
        MATRIX M5 = strassenImpl(A11 + A12, B22);
        MATRIX M6 = strassenImpl(A21 - A11, B11 + B12);
        MATRIX M2 = strassenImpl(A12 - A22, B21 + B22);

        MATRIX C11 = M1 + M4 - M5 + M7;
        MATRIX C12 = M3 + M5;           
        MATRIX C21 = M2 + M4;
        MATRIX C22 = M1 - M2 + M3 + M6;

        MATRIX C = buildFromQuads(C11, C12, C21, C22);
        return C;
    }
}
...