Skip to content

Commit 16f2803

Browse files
jessdtatejessdtate
authored
Merge pull request #1682 from moritzdannhauer/AppendMatrixForSparse
AppendMatrix_reimplementation
2 parents c4b9102 + 777803d commit 16f2803

File tree

8 files changed

+602
-58
lines changed

8 files changed

+602
-58
lines changed

‎src/Core/Algorithms/Math/AppendMatrix.cc‎

Lines changed: 88 additions & 36 deletions
Original file line numberDiff line numberDiff line change
@@ -24,68 +24,120 @@
2424
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
2525
FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
2626
DEALINGS IN THE SOFTWARE.
27+
28+
Author : Moritz Dannhauer (reimplementation)
29+
Date: August 2017
2730
*/
2831

2932
#include <Core/Algorithms/Math/AppendMatrix.h>
3033
#include <Core/Datatypes/DenseMatrix.h>
34+
#include <Core/Datatypes/SparseRowMatrix.h>
35+
#include <Core/Datatypes/SparseRowMatrixFromMap.h>
36+
#include <Core/Datatypes/Legacy/Matrix/MatrixOperations.h>
3137
#include <Core/Algorithms/Base/AlgorithmVariableNames.h>
38+
#include <Core/Datatypes/MatrixTypeConversions.h>
39+
#include <Eigen/Sparse>
3240

3341
using namespace SCIRun::Core::Algorithms;
3442
using namespace SCIRun::Core::Algorithms::Math;
3543
using namespace SCIRun::Core::Datatypes;
3644

45+
const AlgorithmInputName AppendMatrixAlgorithm::InputMatrices("InputMatrices");
46+
3747
AppendMatrixAlgorithm::AppendMatrixAlgorithm()
3848
{
3949
addParameter(Variables::RowsOrColumns, 0);
4050
}
4151

42-
AppendMatrixAlgorithm::Outputs AppendMatrixAlgorithm::run(const AppendMatrixAlgorithm::Inputs& input, const AppendMatrixAlgorithm::Parameters& params) const
52+
53+
bool AppendMatrixAlgorithm::check_dimensions(const Matrix& mat1, const Matrix& mat2, const AppendMatrixAlgorithm::Parameters& params) const
4354
{
44-
DenseMatrixConstHandle lhsPtr = input.get<0>();
45-
DenseMatrixConstHandle rhsPtr = input.get<1>();
46-
if (!lhsPtr || !rhsPtr)
47-
return Outputs(); /// @todo: error
55+
auto rows_m1=mat1.nrows(),rows_m2=mat2.nrows(),cols_m1=mat1.ncols(),cols_m2=mat2.ncols();
56+
if (params == ROWS)
57+
{
58+
if(cols_m1!=cols_m2) return false;
59+
} else
60+
if(rows_m1!=rows_m2) return false;
61+
62+
return true;
63+
}
64+
65+
AppendMatrixAlgorithm::Outputs AppendMatrixAlgorithm::ConcatenateMatrices(const MatrixHandle base_matrix, const std::vector<boost::shared_ptr<Matrix>> input_matrices, const AppendMatrixAlgorithm::Parameters& params) const
66+
{
4867

49-
const DenseMatrix& lhs = *lhsPtr;
50-
const DenseMatrix& rhs = *rhsPtr;
68+
if (input_matrices.size()==0)
69+
return base_matrix;
70+
71+
Outputs outputs = run(boost::make_tuple(base_matrix, input_matrices[0]), params);
72+
for(Eigen::Index c=1; c<input_matrices.size(); c++)
73+
{
74+
auto concatenated = run(boost::make_tuple(outputs, input_matrices[c]), params);
75+
outputs=concatenated;
76+
}
77+
78+
return outputs;
79+
}
5180

52-
if (params == ROWS)
81+
AppendMatrixAlgorithm::Outputs AppendMatrixAlgorithm::run(const AppendMatrixAlgorithm::Inputs& input, const AppendMatrixAlgorithm::Parameters& params) const
82+
{
83+
auto lhsPtr = input.get<0>();
84+
auto rhsPtr = input.get<1>();
85+
if (!lhsPtr || !rhsPtr)
86+
error(" At least two matrices are needed to run this module. ");
87+
88+
if (!((matrixIs::sparse(lhsPtr) && matrixIs::sparse(rhsPtr)) || (matrixIs::dense(lhsPtr) && matrixIs::dense(rhsPtr)) || (matrixIs::column(lhsPtr) && matrixIs::column(rhsPtr))))
5389
{
54-
if (lhs.cols() != rhs.cols())
55-
return Outputs(); /// @todo: error
56-
57-
DenseMatrixHandle output(boost::make_shared<DenseMatrix>(lhs.rows() + rhs.rows(), lhs.cols()));
58-
for (int i = 0; i < lhs.rows(); ++i)
59-
for (int j = 0; j < lhs.cols(); ++j)
60-
(*output)(i,j) = lhs(i,j);
61-
for (int i = 0; i < rhs.rows(); ++i)
62-
for (int j = 0; j < rhs.cols(); ++j)
63-
(*output)(i + lhs.rows(), j) = rhs(i,j);
64-
return output;
90+
error(" Mixing of different matrix types as inputs is not supported. ");
91+
return Outputs();
6592
}
66-
else // columns
93+
94+
if (!check_dimensions(*lhsPtr, *rhsPtr, params))
6795
{
68-
if (lhs.rows() != rhs.rows())
69-
return Outputs(); /// @todo: error
70-
71-
DenseMatrixHandle output(boost::make_shared<DenseMatrix>(lhs.rows(), lhs.cols() + rhs.cols()));
72-
for (int i = 0; i < lhs.rows(); ++i)
73-
for (int j = 0; j < lhs.cols(); ++j)
74-
(*output)(i,j) = lhs(i,j);
75-
for (int i = 0; i < rhs.rows(); ++i)
76-
for (int j = 0; j < rhs.cols(); ++j)
77-
(*output)(i, j + lhs.cols()) = rhs(i,j);
78-
return output;
96+
error(" Input matrix dimensions do not match. ");
97+
return Outputs();
7998
}
99+
100+
Eigen::MatrixXd result;
101+
if(matrixIs::dense(lhsPtr) || matrixIs::column(lhsPtr))
102+
{
103+
if(params == ROWS)
104+
result=Eigen::MatrixXd(lhsPtr->nrows()+rhsPtr->nrows(),lhsPtr->ncols());
105+
else
106+
result=Eigen::MatrixXd(lhsPtr->nrows(),lhsPtr->ncols()+rhsPtr->ncols());
107+
108+
if(matrixIs::dense(lhsPtr))
109+
result << *castMatrix::toDense(lhsPtr), *castMatrix::toDense(rhsPtr);
110+
else
111+
result << *castMatrix::toColumn(lhsPtr), *castMatrix::toColumn(rhsPtr);
112+
113+
if (matrixIs::column(lhsPtr) && (result.rows()==1 || result.cols()==1))
114+
return boost::make_shared<DenseColumnMatrix>(result);
115+
116+
return boost::make_shared<DenseMatrix>(result);
117+
} else
118+
if (matrixIs::sparse(lhsPtr))
119+
return SparseRowMatrixFromMap::concatenateSparseMatrices(*castMatrix::toSparse(lhsPtr),*castMatrix::toSparse(rhsPtr),params==ROWS);
120+
else
121+
{
122+
error(" This matrix type is not supported");
123+
}
124+
125+
return Outputs();
80126
}
81127

82128
AlgorithmOutput AppendMatrixAlgorithm::run(const AlgorithmInput& input) const
83129
{
84-
auto lhs = input.get<DenseMatrix>(Variables::FirstMatrix);
85-
auto rhs = input.get<DenseMatrix>(Variables::SecondMatrix);
86-
87-
auto outputs = run(boost::make_tuple(lhs, rhs), Option(get(Variables::RowsOrColumns).toInt()));
88-
130+
auto lhs = input.get<Matrix>(Variables::FirstMatrix);
131+
auto rhs = input.get<Matrix>(Variables::SecondMatrix);
132+
auto input_matrices = input.getList<Matrix>(AppendMatrixAlgorithm::InputMatrices);
133+
134+
Outputs outputs;
135+
Parameters params=Option(get(Variables::RowsOrColumns).toInt());
136+
outputs = run(boost::make_tuple(lhs, rhs), params);
137+
138+
if (input_matrices.size()>0)
139+
outputs = ConcatenateMatrices(outputs, input_matrices, params);
140+
89141
AlgorithmOutput output;
90142
output[Variables::ResultMatrix] = outputs;
91143
return output;

‎src/Core/Algorithms/Math/AppendMatrix.h‎

Lines changed: 6 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -41,17 +41,17 @@ namespace Math {
4141
{
4242
public:
4343
enum Option { ROWS, COLUMNS };
44-
45-
typedef boost::tuple<SCIRun::Core::Datatypes::DenseMatrixConstHandle, SCIRun::Core::Datatypes::DenseMatrixConstHandle> Inputs;
44+
typedef boost::tuple<SCIRun::Core::Datatypes::MatrixHandle, SCIRun::Core::Datatypes::MatrixHandle> Inputs;
4645
typedef Option Parameters;
47-
typedef SCIRun::Core::Datatypes::DenseMatrixHandle Outputs;
48-
46+
typedef SCIRun::Core::Datatypes::MatrixHandle Outputs;
47+
Outputs ConcatenateMatrices(const Datatypes::MatrixHandle base_matrix, const std::vector<boost::shared_ptr<Datatypes::Matrix>> input_matrices, const AppendMatrixAlgorithm::Parameters& params) const;
4948
Outputs run(const Inputs& input, const Parameters& params) const;
50-
49+
bool check_dimensions(const Datatypes::Matrix& mat1, const Datatypes::Matrix& mat2, const Parameters& params) const;
5150
AppendMatrixAlgorithm();
5251
AlgorithmOutput run(const AlgorithmInput& input) const;
52+
static const AlgorithmInputName InputMatrices;
5353
};
5454

5555
}}}}
5656

57-
#endif
57+
#endif

‎src/Core/Algorithms/Math/Tests/AppendMatrixTests.cc‎

Lines changed: 95 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -45,8 +45,8 @@ TEST(AppendMatrixAlgorithmTests, CanAppendRows)
4545
DenseMatrixHandle m2(SCIRun::TestUtils::matrixNonSquare().clone());
4646
auto result = algo.run(AppendMatrixAlgorithm::Inputs(m1, m2), AppendMatrixAlgorithm::ROWS);
4747

48-
EXPECT_EQ(6, result->rows());
49-
EXPECT_EQ(4, result->cols());
48+
EXPECT_EQ(6, result->nrows());
49+
EXPECT_EQ(4, result->ncols());
5050

5151
std::cout << *result << std::endl;
5252
}
@@ -59,8 +59,8 @@ TEST(AppendMatrixAlgorithmTests, CanAppendColumns)
5959
DenseMatrixHandle m2(SCIRun::TestUtils::matrixNonSquare().clone());
6060
auto result = algo.run(AppendMatrixAlgorithm::Inputs(m1, m2), AppendMatrixAlgorithm::COLUMNS);
6161

62-
EXPECT_EQ(3, result->rows());
63-
EXPECT_EQ(8, result->cols());
62+
EXPECT_EQ(3, result->nrows());
63+
EXPECT_EQ(8, result->ncols());
6464

6565
std::cout << *result << std::endl;
6666
}
@@ -80,7 +80,96 @@ TEST(AppendMatrixAlgorithmTests, ReturnsNullWithSizeMismatch)
8080
TEST(AppendMatrixAlgorithmTests, NullInputReturnsDummyValues)
8181
{
8282
AppendMatrixAlgorithm algo;
83-
83+
8484
auto result = algo.run(AppendMatrixAlgorithm::Inputs(), AppendMatrixAlgorithm::ROWS);
8585
EXPECT_FALSE(result);
86-
}
86+
}
87+
88+
TEST(AppendMatrixAlgorithmTests, AppendSquareSparseMatrix)
89+
{
90+
AppendMatrixAlgorithm algo;
91+
int nr_rows=3, nr_cols=nr_rows;
92+
SparseRowMatrix m1(nr_rows,nr_cols);
93+
for(int i=0;i<nr_rows;i++)
94+
m1.insert(i,i) = i+1;
95+
96+
SparseRowMatrix m2(nr_rows,nr_cols);
97+
for(int i=0;i<nr_rows;i++)
98+
m2.insert(i,i) = i+nr_rows;
99+
100+
auto result = algo.run(AppendMatrixAlgorithm::Inputs(boost::make_shared<SparseRowMatrix>(m1), boost::make_shared<SparseRowMatrix>(m2)), AppendMatrixAlgorithm::ROWS);
101+
auto out = boost::dynamic_pointer_cast<SparseRowMatrix>(result);
102+
103+
for (Eigen::Index k = 0; k < out->nrows(); ++k)
104+
{
105+
for (SparseRowMatrix::InnerIterator it(*out, k); it; ++it)
106+
{
107+
if (it.row()<nr_rows)
108+
ASSERT_TRUE(m1.coeffRef(it.row(),it.col())==it.value());
109+
else
110+
ASSERT_TRUE(m2.coeffRef(it.row()-nr_rows,it.col())==it.value());
111+
}
112+
}
113+
114+
result = algo.run(AppendMatrixAlgorithm::Inputs(boost::make_shared<SparseRowMatrix>(m1), boost::make_shared<SparseRowMatrix>(m2)), AppendMatrixAlgorithm::COLUMNS);
115+
out = boost::dynamic_pointer_cast<SparseRowMatrix>(result);
116+
117+
for (Eigen::Index k = 0; k < out->nrows(); ++k)
118+
{
119+
for (SparseRowMatrix::InnerIterator it(*out, k); it; ++it)
120+
{
121+
if (it.col()<nr_cols)
122+
ASSERT_TRUE(m1.coeffRef(it.row(),it.col())==it.value());
123+
else
124+
ASSERT_TRUE(m2.coeffRef(it.row(),it.col()-nr_cols)==it.value());
125+
}
126+
}
127+
128+
}
129+
130+
TEST(AppendMatrixAlgorithmTests, AppendDenseColumnMatrix)
131+
{
132+
AppendMatrixAlgorithm algo;
133+
int nr_comp1=3, nr_comp2=2; /// assumption for the code in this function: the second column vector is shorter
134+
DenseColumnMatrixHandle m1(boost::make_shared<DenseColumnMatrix>(nr_comp1));
135+
DenseColumnMatrixHandle m2(boost::make_shared<DenseColumnMatrix>(nr_comp2));
136+
137+
for (int i=0;i<nr_comp1;i++)
138+
(*m1)(i) = i+1;
139+
140+
for (int i=0;i<nr_comp2;i++)
141+
(*m2)(i) = 2*(i+1);
142+
143+
auto result = algo.run(AppendMatrixAlgorithm::Inputs(m1, m2), AppendMatrixAlgorithm::COLUMNS);
144+
EXPECT_TRUE(result==nullptr);
145+
result = algo.run(AppendMatrixAlgorithm::Inputs(m1, m2), AppendMatrixAlgorithm::ROWS);
146+
147+
auto out = boost::dynamic_pointer_cast<DenseColumnMatrix>(result);
148+
for (int i = 0; i < out->nrows(); i++)
149+
if(i<nr_comp1)
150+
EXPECT_EQ((*m1)(i),(*out)(i));
151+
else
152+
EXPECT_EQ((*m2)(i-nr_comp1),(*out)(i));
153+
154+
}
155+
156+
TEST(AppendMatrixAlgorithmTests, MixingMatrixInputTypes)
157+
{
158+
AppendMatrixAlgorithm algo;
159+
int nr_rows=3, nr_cols=nr_rows;
160+
SparseRowMatrix sparse_mat(nr_rows,nr_cols);
161+
162+
for(int i=0;i<nr_rows;i++)
163+
sparse_mat.insert(i,i) = i;
164+
165+
int count=0;
166+
167+
DenseMatrixHandle dense_mat(boost::make_shared<DenseMatrix>(nr_rows,nr_cols));
168+
for(int i=0;i<nr_rows;i++)
169+
for(int j=0;j<nr_rows;j++)
170+
(*dense_mat)(i,j)=++count;
171+
172+
auto result = algo.run(AppendMatrixAlgorithm::Inputs(dense_mat, boost::make_shared<SparseRowMatrix>(sparse_mat)), AppendMatrixAlgorithm::ROWS);
173+
EXPECT_TRUE(result==nullptr);
174+
}
175+

‎src/Core/Datatypes/SparseRowMatrixFromMap.cc‎

Lines changed: 54 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -153,4 +153,57 @@ SparseRowMatrixHandle SparseRowMatrixFromMap::appendToSparseMatrixSumming(size_t
153153
SparseRowMatrix originalValuesLarger = sparse + empty;
154154
(*mat) += originalValuesLarger;
155155
return mat;
156-
}
156+
}
157+
158+
159+
SparseRowMatrixHandle SparseRowMatrixFromMap::concatenateSparseMatrices(const SparseRowMatrix& mat1, const SparseRowMatrix& mat2, const bool rows)
160+
{
161+
SparseRowMatrixHandle mat;
162+
size_type offset_rows=0,offset_cols=0;
163+
164+
if ( (rows && mat1.ncols() != mat2.ncols()) || (!rows && mat1.nrows() != mat2.nrows()) )
165+
THROW_INVALID_ARGUMENT(" Matrix dimensions do not match! ");
166+
167+
const size_type nnz = mat1.nonZeros() + mat2.nonZeros();
168+
169+
typedef Eigen::Triplet<double> T;
170+
std::vector<T> tripletList;
171+
tripletList.reserve(nnz);
172+
for (size_type k=0; k < mat1.outerSize(); ++k)
173+
{
174+
for (Core::Datatypes::SparseRowMatrix::InnerIterator it(mat1,k); it; ++it)
175+
{
176+
tripletList.push_back(T(it.row(), it.col(), it.value()));
177+
}
178+
}
179+
180+
if (rows)
181+
{
182+
offset_rows=mat1.nrows();
183+
offset_cols=0;
184+
} else
185+
{
186+
offset_rows=0;
187+
offset_cols=mat1.ncols();
188+
}
189+
190+
for (size_type k=0; k < mat2.outerSize(); ++k)
191+
{
192+
for (Core::Datatypes::SparseRowMatrix::InnerIterator it(mat2,k); it; ++it)
193+
{
194+
tripletList.push_back(T(it.row()+offset_rows, it.col()+offset_cols, it.value()));
195+
}
196+
}
197+
198+
if(rows)
199+
{
200+
mat=boost::make_shared<SparseRowMatrix>(mat1.nrows()+mat2.nrows(),mat1.ncols());
201+
} else
202+
{
203+
mat=boost::make_shared<SparseRowMatrix>(mat1.nrows(),mat1.ncols()+mat2.ncols());
204+
}
205+
206+
mat->setFromTriplets(tripletList.begin(), tripletList.end());
207+
208+
return mat;
209+
}

‎src/Core/Datatypes/SparseRowMatrixFromMap.h‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -62,6 +62,7 @@ namespace SCIRun
6262
static SparseRowMatrixHandle make(size_type rows, size_type cols, const SymmetricValues& values);
6363
static SparseRowMatrixHandle appendToSparseMatrix(size_type rows, size_type cols, const SparseRowMatrix& sparse, const Values& additionalValues);
6464
static SparseRowMatrixHandle appendToSparseMatrixSumming(size_type rows, size_type cols, const SparseRowMatrix& sparse, const Values& additionalValues);
65+
static SparseRowMatrixHandle concatenateSparseMatrices(const SparseRowMatrix& mat1, const SparseRowMatrix& mat2, const bool rows);
6566
private:
6667
SparseRowMatrixFromMap();
6768

0 commit comments

Comments
 (0)