diff --git a/fastlib/sparse/sparse_matrix.h b/fastlib/sparse/sparse_matrix.h index 80e605d757..a0d6d8ac6c 100644 --- a/fastlib/sparse/sparse_matrix.h +++ b/fastlib/sparse/sparse_matrix.h @@ -194,13 +194,20 @@ class SparseMatrix { where i denotes the global row number of A and j denotes the column number */ void ColumnScale(const Vector &vec) { + if (unlikely(!matrix_->Filled())) { + FATAL("You should call EndLoading first...\n"); + } Epetra_Vector temp(View, *map_, (double*)vec.ptr()); matrix_->RightScale(temp); } /** The matrix will be scaled such that A(i,j) = x(i)*A(i,j) * where i denotes the row number of A and j denotes the column number of A. + * You must have called EndLoading() */ void RowScale(const Vector &vec) { + if (unlikely(!matrix_->Filled())) { + FATAL("You should call EndLoading first...\n"); + } Epetra_Vector temp(View, *map_, (double *)vec.ptr()); matrix_->LeftScale(temp); } @@ -219,8 +226,13 @@ class SparseMatrix { /** * Computes the inverse of the sum of absolute values of the rows * of the matrix + * You must have called EndLoading() */ - void InvRowsSums(Vector *result) { + void InvRowSums(Vector *result) { + if (unlikely(!matrix_->Filled())) { + FATAL("You have to call EndLoading before running eigenvalues otherwise " + "it will fail\n"); + } result->Init(dimension_); Epetra_Vector temp(View, *map_, result->ptr()); matrix_->InvRowSums(temp); @@ -228,33 +240,57 @@ class SparseMatrix { /** * Computes the the sum of absolute values of the rows * of the matrix + * You must have called EndLoading() */ - void RowsSums(Vector *result) { + void RowSums(Vector *result) { + if (unlikely(!matrix_->Filled())) { + FATAL("You have to call EndLoading before running eigenvalues otherwise " + "it will fail\n"); + } result->Init(dimension_); Epetra_Vector temp(View, *map_, result->ptr()); matrix_->InvRowSums(temp); + for(index_t i=0; iFilled())) { + FATAL("You have to call EndLoading before running eigenvalues otherwise " + "it will fail\n"); + } result->Init(dimension_); Epetra_Vector temp(View, *map_, result->ptr()); matrix_->InvRowMaxs(temp); } /** Computes the inverse of the sum of absolute values of the columns of the * matrix + * You must have called EndLoading() */ void InvColSums(Vector *result) { + if (unlikely(!matrix_->Filled())) { + FATAL("You have to call EndLoading before running eigenvalues otherwise " + "it will fail\n"); + } result->Init(dimension_); Epetra_Vector temp(View, *map_, result->ptr()); matrix_->InvColSums(temp); } /** * Computes the inv of max of absolute values of the columns of the matrix + * You must have called EndLoading() */ void InvColMaxs(Vector *result) { + if (unlikely(!matrix_->Filled())) { + FATAL("You have to call EndLoading before running eigenvalues otherwise " + "it will fail\n"); + } + result->Init(num_of_columns_); Epetra_Vector temp(View, *map_, result->ptr()); matrix_->InvColMaxs(temp); } @@ -375,6 +411,13 @@ class SparseMatrix { }; +/** + * Sparsem is more like an interface providing basic lagebraic operations + * addition, subtraction multiplicatiion, for sparse matrices. It should + * have been a namespace, but I prefered to make it a class with static + * member functions so that I can declare it as a friend to the SparseMatrix + * class + */ class Sparsem { public: static inline void Add(const SparseMatrix &a, diff --git a/fastlib/sparse/sparse_matrix_impl.h b/fastlib/sparse/sparse_matrix_impl.h index f8c309f753..b5d202ca98 100644 --- a/fastlib/sparse/sparse_matrix_impl.h +++ b/fastlib/sparse/sparse_matrix_impl.h @@ -284,7 +284,11 @@ double SparseMatrix::get(index_t r, index_t c) const { index_t num_of_entries; double *values; index_t *indices; - matrix_->ExtractGlobalRowView(global_row, num_of_entries, values, indices); + if (matrix_->IndicesAreLocal()) { + matrix_->ExtractMyRowView(global_row, num_of_entries, values, indices); + } else { + matrix_->ExtractGlobalRowView(global_row, num_of_entries, values, indices); + } index_t *pos = std::find(indices, indices+num_of_entries, c); if (pos==indices+num_of_entries) { return 0; diff --git a/fastlib/sparse/sparse_matrix_test.cc b/fastlib/sparse/sparse_matrix_test.cc index c2850655a6..a90b75b37d 100644 --- a/fastlib/sparse/sparse_matrix_test.cc +++ b/fastlib/sparse/sparse_matrix_test.cc @@ -40,9 +40,11 @@ class SparseMatrixTest { } } } + void Destruct() { delete smat_; } + void TestInit1() { smat_ = new SparseMatrix(num_of_rows_, num_of_cols_, @@ -71,6 +73,7 @@ class SparseMatrixTest { } NOTIFY("TestInit1 sucess!!\n"); } + void TestInit2() { std::vector rows; std::vector cols; @@ -104,6 +107,7 @@ class SparseMatrixTest { } NOTIFY("TestInit2 success!!\n"); } + void TestInit3() { FILE *fp = fopen("temp.txt", "w"); if (fp==NULL) { @@ -129,10 +133,12 @@ class SparseMatrixTest { } NOTIFY("TestInit3 success!!"); } + void TestCopyConstructor() { smat_ = new SparseMatrix(); NOTIFY("TestCopyConstructor success!!\n"); } + void TestMakeSymmetric() { TestInit1(); smat_->set(2, 3, 1.44); @@ -149,6 +155,103 @@ class SparseMatrixTest { } NOTIFY("Test MakeSymmetric success!!\n"); } + + void TestNegate() { + TestInit1(); + smat_->Negate(); + for(index_t i=0; iget(i,j), -mat_[i][j], + std::numeric_limits::epsilon()); + } + } + NOTIFY("Test Negate success!!"); + } + + void TestColumnScale() { + TestInit1(); + Vector scale; + scale.Init(num_of_cols_); + for(index_t i=0; iEndLoading(); + smat_->ColumnScale(scale); + for(index_t i=0; iget(i,j), j*mat_[i][j], + std::numeric_limits::epsilon()); + } + } + NOTIFY("Test ColumnScale success!!"); + } + + void TestRowScale() { + TestInit1(); + Vector scale; + scale.Init(num_of_cols_); + for(index_t i=0; iEndLoading(); + smat_->RowScale(scale); + for(index_t i=0; iget(i,j), i*mat_[i][j], + std::numeric_limits::epsilon()); + } + } + NOTIFY("Test RowScale success!!"); + } + + void TestRowSums() { + TestInit1(); + Vector row_sums; + smat_->EndLoading(); + smat_->RowSums(&row_sums); + for(index_t i=0; iEndLoading(); + smat_->InvRowSums(&row_sums); + for(index_t i=0; i::epsilon()); + } + NOTIFY("Test InvRowSums success!!"); + } + + void TestInvColMaxs() { + TestInit1(); + Vector col_maxs; + smat_->EndLoading(); + smat_->InvColMaxs(&col_maxs); + for(index_t i=0; i::epsilon()); + } + NOTIFY("Test InvColMaxs success!!"); + + } + void TestEig() { TestInit1(); smat_->EndLoading(); @@ -159,6 +262,7 @@ class SparseMatrixTest { eigvectors.PrintDebug(); NOTIFY("Test Eigenvector success!!\n"); } + void TestLinSolve() { TestInit1(); Vector b,x; @@ -245,6 +349,24 @@ class SparseMatrixTest { TestMakeSymmetric(); Destruct(); Init(); + TestNegate(); + Destruct(); + Init(); + TestColumnScale(); + Destruct(); + Init(); + TestRowScale(); + Destruct(); + Init(); + TestRowSums(); + Destruct(); + Init(); + TestInvRowSums(); + Destruct(); + Init(); + TestInvColMaxs(); + Destruct(); + Init(); TestEig(); Destruct(); Init();