我试图划分和征服矩阵乘法,所以我可以并行化它,但我得到一半随机垃圾数和结果的一半0,例如在2x2矩阵“[[15909360,0] [15909360,0]]”。这是我到目前为止基于the algorithm on Wikipedia 但我不知道从哪里去的地方。我还没有使用指针进行分区或线程。这是家庭作业btw。
void partition(const std::vector<std::vector<IntElement> >& m, std::vector<std::vector<IntElement> >& m11, std::vector<std::vector<IntElement> >& m12,
std::vector<std::vector<IntElement> >& m21, std::vector<std::vector<IntElement> >& m22, int n){
for(int i=0;i<n/2;i++)
for(int j=0;j<n/2;j++){
m11[i][j] = m[i][j]; // top left
m12[i][j] = m[i][j + n / 2]; // top right
m21[i][j] = m[i + n / 2][j]; // bottom left
m22[i][j] = m[i + n / 2][j + n / 2]; // bottom right
}
};
void add(std::vector<std::vector<IntElement> >& C, std::vector<std::vector<IntElement> >& T, int n){
if(n==1){
C[0][0] += C[0][0] + T[0][0];
}
else{
std::vector<std::vector<IntElement> > c11(n/2, std::vector<IntElement>(n/2)), c12(n/2, std::vector<IntElement>(n/2)),
c21(n/2, std::vector<IntElement>(n/2)), c22(n/2, std::vector<IntElement>(n/2));
std::vector<std::vector<IntElement> > t11(n/2, std::vector<IntElement>(n/2)), t12(n/2, std::vector<IntElement>(n/2)),
t21(n/2, std::vector<IntElement>(n/2)), t22(n/2, std::vector<IntElement>(n/2));
partition(C, c11, c12, c21, c22, n);
partition(T, t11, t12, t21, t22, n);
add(c11, t11, n/2);
add(c12, t12, n/2);
add(c21, t21, n/2);
add(c22, t22, n/2);
}
};
void multiply(std::vector<std::vector<IntElement> >& C, const std::vector<std::vector<IntElement> >& A,
const std::vector<std::vector<IntElement> >& B, int n){
if(n==1)
C[0][0] += A[0][0] * B[0][0];
else{
std::vector<std::vector<IntElement> > T(n, std::vector<IntElement>(n));
std::vector<std::vector<IntElement> > a11(n/2, std::vector<IntElement>(n/2)), a12(n/2, std::vector<IntElement>(n/2)),
a21(n/2, std::vector<IntElement>(n/2)), a22(n/2, std::vector<IntElement>(n/2));
std::vector<std::vector<IntElement> > b11(n/2, std::vector<IntElement>(n/2)), b12(n/2, std::vector<IntElement>(n/2)),
b21(n/2, std::vector<IntElement>(n/2)), b22(n/2, std::vector<IntElement>(n/2));
std::vector<std::vector<IntElement> > c11(n/2, std::vector<IntElement>(n/2)), c12(n/2, std::vector<IntElement>(n/2)),
c21(n/2, std::vector<IntElement>(n/2)), c22(n/2, std::vector<IntElement>(n/2));
std::vector<std::vector<IntElement> > t11(n/2, std::vector<IntElement>(n/2)), t12(n/2, std::vector<IntElement>(n/2)),
t21(n/2, std::vector<IntElement>(n/2)), t22(n/2, std::vector<IntElement>(n/2));
partition(A, a11, a12, a21, a22, n);
partition(B, b11, b12, b21, b22, n);
partition(C, c11, c12, c21, c22, n);
partition(T, t11, t12, t21, t22, n);
multiply(c11, a11, b11, n/2);
multiply(c12, a11, b12, n/2);
multiply(c21, a21, b11, n/2);
multiply(c22, a21, b12, n/2);
multiply(t11, a12, b21, n/2);
multiply(t12, a12, b22, n/2);
multiply(t21, a22, b21, n/2);
multiply(t22, a22, b22, n/2);
add(C, T, n);
}
return;
};
SquareMatrix& SquareMatrix::operator*=(const SquareMatrix& m){
std::vector<std::vector<IntElement> > C(n, std::vector<IntElement>(n));
multiply(C, elements, m.elements, n);
elements = C;
return *this;
}
SquareMatrix operator*(const SquareMatrix& a, const SquareMatrix& b){
SquareMatrix c = a;
c *= b;
return c;
}
编辑:我改变了C [0] [0] + = C [0] [0] + T [0] [0]; in add()to C [0] [0] + = T [0] [0];此外,我做了一个非分区函数,它基本上反过来并在乘法和添加后将分区放回C和T:
void unpartition(std::vector<std::vector<IntElement> >& m,std::vector<std::vector<IntElement> >& m11, std::vector<std::vector<IntElement> >& m12,
std::vector<std::vector<IntElement> >& m21, std::vector<std::vector<IntElement> >& m22, int n){
for(int i=0;i<n/2;i++)
for(int j=0;j<n/2;j++){
m[i][j] = m11[i][j]; // top left
m[i][j + n / 2] = m12[i][j]; // top right
m[i + n / 2][j] = m21[i][j]; // bottom left
m[i + n / 2][j + n / 2] = m22[i][j]; // bottom right
}
}
在修复了IntElement类的默认构造函数后,我的向量得到了正确的初始化。
您的partition
函数会在源向量中复制数据。在multiply
和add
中,最后一系列调用(multiply(c11, a11, b11, n/2)
,add(c11, t11, n/2)
)将修改这些子矩阵对象。这不会修改原始的C
或T
矩阵。 T
将保持零填充,C
将偶尔从1x1矩阵更新。您需要将结果“取消”分配回适当的矩阵。
add
中的单个元素案例是错误的。 C[0][0] += C[0][0] + T[0][0];
应该是C[0][0] += T[0][0];
或C[0][0] = C[0][0] + T[0][0];
(但不是两者)。
变量T,a11,b11,c11,d11
中的内部向量永远不会被初始化。它们包含在这4个parititons调用中使用的垃圾。
编辑:上述声明是错误的,因为Alan Birtles指出它们是由IntElement
的默认构造函数初始化的。我仍然相信它可能只是int
的typedef,在这种情况下,我的论点是站立的。