// dcs.cpp : Defines the entry point for the console application. // // PLAY WITH THESE PARAMETERS TO TEST DIFFERENT MATRIX CONFIGURATIONS #define SCALAR std::complex #define INDEX ptrdiff_t #define ORDER Eigen::ColMajor // [???!!!] uncommenting these lines causes crush //#ifdef EIGEN_DEFAULT_DENSE_INDEX_TYPE //#undef EIGEN_DEFAULT_DENSE_INDEX_TYPE //#endif //#define EIGEN_DEFAULT_DENSE_INDEX_TYPE INDEX #include #include "Eigen/Sparse" #include "Eigen/Dense" #include "dcs.h" #include #include typedef Eigen::Matrix DenseMatrix; typedef Eigen::Matrix DenseMatrixT; typedef Eigen::SparseMatrix SparseMatrix; typedef Eigen::SparseMatrix SparseMatrixT; typedef dcs::HyperSparseMatrix HyperSparseMatrix; typedef dcs::HyperSparseMatrix HyperSparseMatrixT; typedef Eigen::Matrix MatrixXi_r; typedef Eigen::Matrix MatrixXd_r; typedef Eigen::Matrix, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> MatrixXc_r; void tests(INDEX rows, INDEX cols); void ContQ(); void triplets_test(DenseMatrix& mat_d, HyperSparseMatrix& mat_h, INDEX rows, INDEX cols); template bool check(DM& dm, SM& hm, std::string pass_phrase="passed!\n") { DM test = dm; test -= hm; if ((test).cwiseAbs().maxCoeff() < 1.e-10) { { DenseMatrix d_h = hm.diagonal().eval().transpose(); DenseMatrix d_d = dm.diagonal().transpose(); DenseMatrix diff = d_d - d_h; if (diff.cwiseAbs().maxCoeff()>1.e-10) { std::cout << "[ERROR] Matrix diagonal extraction failed: diagonals are different.\n" << std::flush; std::cout << "Should be: \n"< 0 && new_cols > 0 && ((bool)dm.IsRowMajor == (bool)hm.IsRowMajor)) //[???!!!]line 77 raises excepition if matrices have different storage orders { new_rows = std::rand() % (2*new_rows)+1; new_cols = std::rand() % (2*new_cols)+1; dm2.conservativeResizeLike(DenseMatrix::Zero(new_rows, new_cols)); hm2.conservativeResize(new_rows, new_cols); test = dm2; test -= hm2; // error here if matrices have different storage order; if (test.cwiseAbs().maxCoeff() > 1.e-10) { std::cout << "[ERROR] Conservative resize to [" << new_rows << "," << new_cols << "] failed.\n"; std::cout << "Should be: \n" << dm2 << "\n"; std::cout << "Instead, have: \n" << hm2 << "\n"; std::cout << "Difference is:\n" << test << "\n"; ContQ(); } } } std::cout << pass_phrase; return true; } else { std::cout << "failed: expected contents:\n"; std::cout << dm; std::cout << "\nActual contents:\n"; std::cout << hm; std::cout << "\nDifference:\n"; std::cout << test << "\n"; ContQ(); return false; }; }; int main() { INDEX rows = 7; INDEX cols = 7; int passes = 1; std::cout << "enter matrix dimmensions:\n number of rows: "; std::cin >> rows; std::cout << " number of columns: "; std::cin >> cols; std::cout << "enter number of passes:\n"; std::cin >> passes; while (passes-- > 0) { tests(rows, cols); } std::cout << "[END OF TESTS]\n"; ContQ(); } void tests(INDEX rows, INDEX cols) { DenseMatrix mat_d(rows, cols); SparseMatrix mat_s(rows, cols); ////////////////////////////////////////////////////////// { std::cout << "[INFO] Coefficientwise initialization test [coeffRef()]... "; HyperSparseMatrix mat_h(rows, cols); std::srand((unsigned)std::time(NULL)); mat_d.setZero(); for (INDEX i = 0; i < rows; i++) { INDEX r = (INDEX)std::rand() % rows; INDEX c = (INDEX)std::rand() % cols; SCALAR value = (SCALAR)std::rand(); mat_h.coeffRef(r, c) = value; mat_d(r, c) = value; } check(mat_d, mat_h); } ////////////////////////////////////////////////////////// std::cout << "[INFO] Initialization from triplets... "; { HyperSparseMatrix mat_h(rows, cols); triplets_test(mat_d, mat_h, rows, cols); std::cout << "passed.\n"; ////////////////////////////////////////////////////////// std::cout << "[INFO] Submatrix access test... "; HyperSparseMatrix mh = mat_h.bottomRightCorner(std::max((INDEX)1,rows / 2), std::max((INDEX)1,cols / 2)); DenseMatrix md = mat_d.bottomRightCorner(std::max((INDEX)1,rows / 2), std::max((INDEX)1,cols / 2)); check(md, mh); ////////////////////////////////////////////////////////// std::cout << "[INFO] Copy constructor test 1 [same alignment]... "; { HyperSparseMatrix mat_h2(mat_h); check(mat_d, mat_h2); } std::cout << "[INFO] Copy constructor test 2 [different alignment]... "; { HyperSparseMatrixT mat_h2(mat_h); check(mat_d, mat_h2); } ////////////////////////////////////////////////////////// { std::cout << "[INFO] Sparse to dense (same alignment)... "; mat_d = mat_h.toDense(); check(mat_d, mat_h); std::cout << "[INFO] Sparse to dense (different alignment)... "; HyperSparseMatrixT mht(mat_h); mat_d = mht.toDense(); check(mat_d, mat_h); } ////////////////////////////////////////////////////////// { std::cout << "[INFO] Dense to sparse (same alignment)... "; mat_h.setZero(); mat_h = mat_d.sparseView(); check(mat_d, mat_h); std::cout << "[INFO] Dense to sparse (different alignment)... "; HyperSparseMatrixT mht; mht = mat_d.sparseView(); check(mat_d, mht); } ////////////////////////////////////////////////////////// //{ // std::cout << "[INFO] Matrix transpose... "; // mh = mat_h.transpose().eval(); // mat_d.transposeInPlace(); // check(mat_d, mh); // std::cout << "[INFO] Matrix transpose in place... "; // mat_h.transposeInPlace(); // check(mat_d, mat_h); //} // std::cout << "[INFO] InnerIterators visual check... "< T; std::vector tripletList; mat_h.resize(rows, cols); mat_d.resize(rows, cols); tripletList.reserve(rows); mat_d.setZero(); //for (int r = 0; r < 100; r++) { mat_h.setZero(); for (INDEX i = 0; i < rows; i++) { INDEX r = (INDEX)std::rand() % rows; INDEX c = (INDEX)std::rand() % cols; SCALAR value = (SCALAR)std::rand(); tripletList.push_back(T(r, c, value)); mat_d(r, c) += value; } mat_h.setFromTriplets(tripletList.begin(), tripletList.end()); check(mat_d, mat_h, "*"); } }; void ContQ() { char c; std::cout << "Press any character and ENTER to continue... " << std::flush; std::cin>>c; }