Strassen Vinograde算法

Question

我已经编写了Strassen算法的实现，但是由于递归创建静态数组，它运行得很慢。我知道动态数组可以解决这个问题，但我不允许使用它们。所以这个版本的主要思想是:我们使用每个正方形的角元素，而不是复制每个子正方形。我认为类似这样的东西，我们应该应用于我递归创建的数组。通常，主要的问题是，即使在较大的值下，该算法也比通常慢几倍。

……的。

    #include "pch.h"
    #include
    #include
    #include
    #include
    #include
    #include
    #pragma comment(linker, "/STACK:5813242100")
    using namespace std;
    const int sizs = 256;

    void vivod(int matrix[][256], int n);
    void Matrix_Add(int a[][256], int b[][256], int c[][256], int n, int x1, int y1, int x2, int y2); 
    void Matrix_Sub(int a[][256], int b[][256], int c[][256], int n, int x1, int y1, int x2, int y2);
    void Matrix_Multiply(int a[][256], int b[][256], int c[][256], int x1, int y1, int x2, int y2, int n);
    void strassen(int a[][256], int b[][256], int c[][256], int m, int n, int x1, int y1, int x2, int y2);
    void Naive_Multiply(int a[][256], int b[][256], int c[][256], int n);
    
    int main()
    {
        setlocale(LC_ALL, "Russian");
        int n;
        cout << "Enter the N:";
        cin >> n;
        const int m = 256;
        int A[m][m];
        int B[m][m];
        int C[m][m];
        int k[m][m];
        for (int i = 0; i < n; i++)
            for (int j = 0; j < n; j++)
                A[i][j] = 0;
        for (int i = 0; i < n; i++)
            for (int j = 0; j < n; j++)
                B[i][j] = 0;
        for (int i = 0; i < n; i++)
            for (int j = 0; j < n; j++)
                A[i][j] = rand() % 10;
        for (int i = 0; i < n; i++)
            for (int j = 0; j < n; j++)
                B[i][j] = rand() % 10;
        cout << "First Matrix:" << endl;
        //vivod(A, n);
        cout << "Second Matrix:" << endl;
        //vivod(B, n);
        int begin = clock();
        //for (int i =0; i < 100; i++)
        Naive_Multiply(A, B, k, n);
        int end = clock();
        cout << "Naive Multiply tacts: " << end - begin << endl;
        //vivod(k, n);
        int begin2 = clock();
        //for (int i = 0; i < 100; i++)
        strassen(A, B, C, n, n, 0, 0, 0, 0);
        int end2 = clock();
        cout << "Shtrassen tacts: " << end2 - begin2 << endl;
        //vivod(C, n);
        system("pause");
        return 0;
    }

    void strassen(int a[][256], int b[][256], int c[][256], int m, int n, int x1, int y1, int x2, int y2) {
        m = n / 2;
        if (m != 1)
        {
            int s1[sizs][sizs];
            int s2[sizs][sizs];
            int s3[sizs][sizs];
            int s4[sizs][sizs];
            int s5[sizs][sizs];
            int s6[sizs][sizs];
            int s7[sizs][sizs];
            int s8[sizs][sizs];
            int m1[sizs][sizs];
            int m2[sizs][sizs];
            int m3[sizs][sizs];
            int m4[sizs][sizs];
            int m5[sizs][sizs];
            int m6[sizs][sizs];
            int m7[sizs][sizs];
            int t1[sizs][sizs];
            int t2[sizs][sizs];
            int c11[sizs][sizs];
            int c12[sizs][sizs];
            int c21[sizs][sizs];
            int c22[sizs][sizs];
            Matrix_Add(a, a, s1, m, x1 + m, y1, x1 + m, y1 + m);
            Matrix_Sub(s1, a, s2, m, 0, 0, x1, y1);
            Matrix_Sub(a, a, s3, m, x1, y1, x1 + m, y1);
            Matrix_Sub(a, s2, s4, m, x1, y1 + m, 0, 0);
            Matrix_Sub(b, b, s5, m, x2, y2 + m, x2, y2);
            Matrix_Sub(b, s5, s6, m, x2 + m, y2 + m, 0, 0);
            Matrix_Sub(b, b, s7, m, x2 + m, y2 + m, x2, y2 + m);
            Matrix_Sub(s6, b, s8, m, 0, 0, x2 + m, y2);
            strassen(s2, s6, m1, m, m, 0, 0, 0, 0);
            strassen(a, b, m2, m, m, x1, y1, x2, y2);
            strassen(a, b, m3, m, m, x1, y1 + m, x2 + m, y2);
            strassen(s3, s7, m4, m, m, 0, 0, 0, 0);
            strassen(s1, s5, m5, m, m, 0, 0, 0, 0);
            strassen(s4, b, m6, m, m, 0, 0, x2 + m, y2 + m);
            strassen(a, s8, m7, m, m, x1 + m, y1 + m, 0, 0);

            Matrix_Add(m1, m2, t1, m, 0, 0, 0, 0);
            Matrix_Add(t1, m4, t2, m, 0, 0, 0, 0);
            Matrix_Add(m2, m3, c11, m, 0, 0, 0, 0);
            Matrix_Sub(t2, m7, c21, m, 0, 0, 0, 0);
            Matrix_Add(t1, m5, c12, m, 0, 0, 0, 0);
            Matrix_Add(c12, m6, c12, m, 0, 0, 0, 0);
            Matrix_Add(t2, m5, c22, m, 0, 0, 0, 0);
            for (int i = 0; i < n / 2; i++)
            {
                for (int j = 0; j < n / 2; j++)
                {
                    c[i + n - 2 * m][j + n - 2 * m] = c11[i][j];
                    c[i + n - 2 * m][j + n - m] = c12[i][j];
                    c[i + n - m][j + n - 2 * m] = c21[i][j];
                    c[i + n - m][j + n - m] = c22[i][j];
                }
            }
        }
        else
        {
            Matrix_Multiply(a, b, c, x1, y1, x2, y2, n);
        }
    }
    void vivod(int matrix[][256], int n)
    {
        for (int i = 0; i < n; i++)
        {
            for (int j = 0; j < n; j++)
            {
                cout << matrix[i][j] << " ";
            }
            cout << endl;
        }
        cout << endl;
    }
    void Matrix_Add(int a[][256], int b[][256], int c[][256], int n, int x1, int y1, int x2, int y2)
    {
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                c[i][j] = a[i + x1][j + y1] + b[i + x2][j + y2];
            }
        }
    }

    void Matrix_Sub(int a[][256], int b[][256], int c[][256], int n, int x1, int y1, int x2, int y2)
    {
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                c[i][j] = a[i + x1][j + y1] - b[i + x2][j + y2];
            }
        }
    }
    void Matrix_Multiply(int a[][256], int b[][256], int c[][256], int x1, int y1, int x2, int y2, int n)
    {
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                c[i][j] = 0;
                for (int t = 0; t < n; t++) {
                    c[i][j] = c[i][j] + a[x1 + i][y1 + t] * b[x2 + t][y2 + j];
                }
            }
        }
    }
    void Naive_Multiply(int a[][256], int b[][256], int c[][256], int n)
    {
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                c[i][j] = 0;
                for (int t = 0; t < n; t++) {
                    c[i][j] = c[i][j] + a[i][t] * b[t][j];
                }
            }
        }
    }

由于数组数量庞大，它可能甚至无法启动，但我已经启动了它，下面是测试：

在N= 128和256的情况下，简单乘法需要近10秒的时间，而同时我在等待Strassen 1~5分钟。

Answer 1

我不太确定这个问题是否有正确的标记，因为你写的东西看起来很像普通的C代码。如果您确实需要静态存储，那么始终可以回到std::array，它或多或少是一个普通C样式数组的包装器。

因此，您应该能够用std::array, 256>代替D2。注意，在大多数情况下，有一个数组并使用索引来访问各个元素是有益的，但这是另一天的事情。

请注意，std::array不具有构造函数的功能，因此仍然需要手动填充构造函数。

因为最近您可以将变量声明为constexpr，这是对局部变量m的一个很好的替代。

现在关于实际代码：

1. 您不应该使用using namespace std; --这是一种错误的做法，会不必要地污染整个命名空间。只要有必要，养成编写std::的习惯是非常有益的。
2. 你原来的主循环都是一样长的，所以只要填一次元素就行了，不要太频繁地插入。
3. 关于各个函数，您应该能够利用来自algorithm库的相当多的功能。

首先让我们看看Naive_Multiply：

void Naive_Multiply(int a[][256], int b[][256], int c[][256], int n)
{
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < n; j++) {
            c[i][j] = 0;
            for (int t = 0; t < n; t++) {
                c[i][j] = c[i][j] + a[i][t] * b[t][j];
            }
        }
    }
}

在这里您可以使用std::inner_product

constexpr int rowSize = 256;
using row = std::array;
using mat = std::array

mat Naive_Multiply(const mat& a, const mat& b)
{
    mat c;
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < n; j++) {
            auto multiplyRowElement = [j](const int a, const row& b) {
                return a * b[j];
            };
            c[i][j] = std::inner_product(a[i].begin(), a[i].end(), 
                                         b.begin(), b.end(), 0, 
                                         std::plus, multiplyRowElement);
        }
    }
    return c;
}

请注意，我们还直接返回操作的结果，而不是将其作为inout参数传递给函数，这要干净得多。

如果任何附加参数不等于0，则otehr函数似乎会创建访问冲突。在这种情况下，您将遍历数组的边框。