3
\$\begingroup\$

I have written an implementation of the Strassen Vinograde Algorithm, but it's running slowly because of recursive creation of static arrays. I know that dynamic arrays would solve this problem, but I'm not allow to use them. So the main idea of this version of Strassen: We use corner elements of each square instead of copying every sub-square. I think something like this we should apply for the arrays that I create recursively. In general, the main problem is that the algorithm is several times slower than usual even at larger values ​​of n.

    #include "pch.h"
    #include<iostream>
    #include<cstdio>
    #include<conio.h>
    #include<cstdlib>
    #include<cmath>
    #include<ctime>
    #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];
                }
            }
        }
    }

It probably may not even start, because of large amount of arrays, but i have launched it and here is the tests:

enter image description here

With N = 128 and 256 naive multiplication takes nearly 10 seconds, when at the same time I'm waiting for Strassen for ~1-5 minutes.

\$\endgroup\$
  • 1
    \$\begingroup\$ Brainstorming here, what do you mean no dynamic arrays? You're probably going to have to pre-declare storage for all the temporaries, using depth as well as row,col boundaries, so the recursive calls can share space in each NxN temporary. It might be ugly, but it will avoid declaring a bunch of extra MAX by MAX temporary arrays for tiny subproblems. \$\endgroup\$ – Kenny Ostrom Dec 16 '18 at 18:24
  • \$\begingroup\$ Maybe we need source storage, but their number for example for 256x256 matrices will be very large. \$\endgroup\$ – Emik Dec 16 '18 at 18:37
  • \$\begingroup\$ And I do not fully understand your idea. \$\endgroup\$ – Emik Dec 16 '18 at 18:43
  • \$\begingroup\$ In addition, we can create dynamic arrays in main function, but it's not allowed to do it recursively. \$\endgroup\$ – Emik Dec 16 '18 at 19:29
  • \$\begingroup\$ For reference, original question: stackoverflow.com/questions/53802871 \$\endgroup\$ – chtz Dec 21 '18 at 10:15
1
\$\begingroup\$

I am not really sure if this question is correctly tagged, as what you wrote seems a lot like plain C-Code. If you really need static storage, you can always just fall back to std::array, which is more or less a wrapper around a plain C-style array.

So you should be able to replace int[256][256] by std::array<std::array<int, 256>, 256>. Note that in most times it is beneficial to have one array and use indexing to access the individual elements, but that is for another day.

Note that std::array does not feature a constructor, so you would still have to manually fill it.

Since recently you can declare variables as constexpr, which is a good replacement for your local variable m.

Now regarding the actual code:

  1. You should not use using namespace std; This is bad practice that will pollute the entire namespace needlessly. It is quite beneficial to get into the habit of writing std:: whenever necessary.

  2. Your original loops in main all run the same length so justfill the elements once and dont interate so often

  3. Regarding the respective functions you should be able to utilize quite a bit of functionality from the algorithm library.

First lets have a look at 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];
            }
        }
    }
}

Here you can use std::inner_product

constexpr int rowSize = 256;
using row = std::array<int, m>;
using mat = std::array<row , m>

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<int>, multiplyRowElement);
        }
    }
    return c;
}

Note that we also directly return the result of the operation rather than passing it as an inout parameter to the function, which is much cleaner.

It seems your otehr functions create access violations if any of the additional parameters is not equal to 0. In that case you iterate over the border of the array.

\$\endgroup\$
  • \$\begingroup\$ I am very grateful for your hints, but before using various functions, I must learn the principles of their work. So, that is why i'm not using it, i'm just learning the language. And what about the main problem? How to do the algorithm without recursive array creation? \$\endgroup\$ – Emik Dec 18 '18 at 22:18
  • \$\begingroup\$ @Emik The C++-standard-library is part of the language (not part of the core language, but part of the C++ standard) \$\endgroup\$ – chtz Dec 21 '18 at 10:14

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.