Back to General and Gameplay Programming

Strassen's Algorithm

General and Gameplay Programming Programming

Started by mojtaba_e May 22, 2006 01:57 PM

7 comments, last by mojtaba_e 17 years, 10 months ago

mojtaba_e

174

Author

May 22, 2006 01:57 PM

Hi friends, I look on the web for Strassen'S algorithm source Code. (For multiplitation two matrixes) But i could't found any source code for this algorithm.[crying] Can you help me to find a source code for this algorithm. Thanks,

--Mojtaba--

kuphryn

210

May 22, 2006 03:22 PM

algorithm?

Kuphryn

rohde

432

May 22, 2006 03:27 PM

Full implementation with lots of info:

strassen.pdf

"We confess our little faults to persuade people that we have no large ones." -Francois de La Rochefoucauld (1613 - 1680). | My blog

mojtaba_e

174

Author

May 22, 2006 03:49 PM

Thanks friends,
I read this PDF file.But this pdf is not complete!
I wanna a complete project(C++ or VB ,...) to perform this Algorithm.
This project is my university project.
I think that this project is difficult!

--Mojtaba--

Bnty_Hntr

187

May 22, 2006 04:06 PM

Strassen's algorithm is a pretty complex thing.
If you need someone to find you the source code, you probably would have a hard time understanding it.

Plus, as you can even see from the paper posted by rohde, Strassen's algorithm only out-performs regular matrix multiplication on some cases, and not always by that much. You need very large sample sizes for Strassen to make a big difference.

So, to implement Strassen for the regular multiplication that takes place in computer graphics and games is kind of an overhaul. I wouldn't suggest it unless you are working with large matrices.

The only other reason you would be asking for Strassen would be homework, and then we couldn't just hand you the answer anyway. I would be willing to post some pseudo-code to get you started, but you will have to translate that into the actual working code on your own.

- Bnty_HntrCome visit Triple Buffer, a great source for DirectX and other tutorials.Check out my Blog.Or, join Triple Buffer by checking out our Help Wanted post.

mojtaba_e

174

Author

May 22, 2006 04:53 PM

OK Bnty_Hntr,
I don't like to Copy this source code and give this copy to my master.
That shames!
I think that about this algorithm for more time.
I can divide matrix to four matrix,but i can't combine values of this recursion.
for example in simple recursion mode we write :
- - - -
|A11|A12| |C11|C12|
|-------| * ..... = |-------|
|A21|A22| |C21|C22|
- - - -

C(11)=A11*B11+A12*B21

In this code we must return a matrix(C11).
My problem is that i can't return this matrix.

--Mojtaba--

Bnty_Hntr

187

May 22, 2006 05:54 PM

You can divide the matrices into 4 parts, that is a start.

I don't quite know what you mean when you say that you can't combine the values of the recursion. You will have to expand on that one.

Here are a couple links to help understand Strassen:
@Wikipedia
@MathWorld

Here, I'll pump out some pseudo code for you to look at.
It isn't too hard of code to understand actually, even though the algorithm isn't technically simple.

// You will need a Matrix type that handles +, -, and *Strassen(Matrix A, Matrix B){	// When the dimension of the matrices goes under a certain threshold	// You will want to multiply by the normal method	// Strassen gets slow when used with low dimensions	if( dim(A, B) < SOME_THRESHOLD )		return A * B	// You said you can handle this next part	A11 = Get A11 from A	A12 = Get A12 from A	.	.	.	B22 = Get B22 from B	// Now we get the 7 different parts, as defined by Strassen.	// We do this recursively, using the Strassen algorithm	// for all multiplications.	M1 = Strassen(A11 + A22, B11 + B22)	M2 = Strassen(A21 + A22, B11)	.	.	.	M7 = Strassen(A12 - A22, B21 + B22)	// Now we calculate the 4 different parts of C, the answer.	C11 = M1 + M4 - M5 + M7	C12 = M3 + M5	C21 = M2 + M4	C22 = M1 - M2 + M3 + M6	// The only thing left to do is to combine C11, C12, C21 and C22	// into the final matrix C.  	// You can do this with a couple simple loops.	return C}

I took this all straight from Wikipedia's description.
It isn't too hard to get some code from there, or from what I gave you.

I hope that all helps.
If you need help actually understanding the algorithm instead of just needing code, feel free to ask.

Good Luck!

- Bnty_HntrCome visit Triple Buffer, a great source for DirectX and other tutorials.Check out my Blog.Or, join Triple Buffer by checking out our Help Wanted post.

NotAYakk

876

May 22, 2006 06:06 PM

A scetch of a design for a matrix class that supports submatrixes:

The "subrect" class encapsulates the information that describes an integer aligned sub-region of a grid.

The "submatrix" class is a pointer-sematic class that gives operator() access to a sub matrix of another matrix. It supports making sub-matrixes of sub-matrixes. Any rectangular window in a parent matrix works. It also supports "skip" matrixes (like every second line/column) for the heck of it.

The "matrix" class actually holds data.

Submatrix's are only valid for as long as their root "matrix" is valid.

Possible improvements:
Actually try to compile the code.
Implement operator*, +, etc.

struct subrect {  int width;  int height;  int linestride;  int elemstride;  int offset_x;  int offset_y;  int offset() {return offset_x + offset_y;}  subrect makesubrect(int width_, int height_) {    subrect retval;    retval.width = width_;    retval.height = height_;    retval.elemstride = 1;    retval.linestride = width;    retval.offset_x = 0;    retval.offset_y = 0;    return retval;  }  subrect(int width_, int height_) {    *this = makesubrect(width_, height_);  }  subrect makesubrect(    int offset_x, int offset_y,    int width_, int height_,    int horz_skip=1, int vert_skip=1  ) {    subrect retval;    retval.linestride = vert_skip;    retval.elemstride = horz_skip;    retval.width = width_;    retval.height = height_;    retval.offset_x = offset_x;    retval.offset_y = offset_y;    return retval;  }  subrect makesubrect( subrect parent, subrect child ) {    subrect retval;    retval.linestride = parent.linestride*child.linestride;    retval.elemstride = parent.elemstride*child.elemstride;    retval.width = child.width;    retval.height = child.height;    retval.offset_x = parent.offset_x+child.offset_x*parent.elemstride;    retval.offset_y = parent.offset_y+child.offset_y*parent.linestride;    return retval;  }  subrect makesubrect(    subrect parent,    int offset_x, int offset_y,    int width_, int height_,    int horz_skip=1, int vert_skip=1  )  {    subrect temp = makesubrect(offset_x, offset_y, width_, height_, horz_skip, vert_skip );    return makesubrect(parent, temp);  }  subrect( subrect parent, subrect child ) {    *this = makesubrect(parent, child);  }};template<typename T>struct submatrix {  subrect rect;  T* data;  T& operator()(int a, int b) {    if (a<0 || a >=rect.width || b<0 || b>=rect.height) {      ASSERT(false);    }    return data[rect.offset()+b*rect.linestride+a*rect.elemstride];  }  T const& operator()(int a, int b) const {    if (a<0 || a >=rect.width || b<0 || b>=rect.height) {      ASSERT(false);    }    return data[offset+height*linestride+width];  }  submatrix( submatrix parent, subrect child ):    rect(parent.rect, child), data(parent.data)  {}  submatrix& operator=(submatrix const& other) {    ASSERT(rect.width = other.rect.width);    ASSERT(rect.height = other.rect.height);    {for (int i = 0; i < rect.width; ++i) {      {for (int j = 0; j < rect.height; ++j) {        (*this)(i,j) = other(i,j);      }}    }}    return *this;  }};template<typename T>struct matrix {  private:    std::vector<T> data_storage;  public:    matrix(int n, int m) {      rect = subrect(n,m);      data_storage.resize(n*m);      data = &data_storage.front();    }    matrix(const matrix& other) {      data_storage = other.data_storage;      rect = other.rect;      data = &data_storage.front();    }    matrix& operator=(const matrix& other) {      submatrix::operator=(other);      return *this;    }};