| 44 | * and zero coefficients respectively. |
| 45 | */ |
| 46 | template<typename Scalar,int Opt1,int Opt2,typename StorageIndex> void |
| 47 | initSparse(double density, |
| 48 | Matrix<Scalar,Dynamic,Dynamic,Opt1>& refMat, |
| 49 | SparseMatrix<Scalar,Opt2,StorageIndex>& sparseMat, |
| 50 | int flags = 0, |
| 51 | std::vector<Matrix<StorageIndex,2,1> >* zeroCoords = 0, |
| 52 | std::vector<Matrix<StorageIndex,2,1> >* nonzeroCoords = 0) |
| 53 | { |
| 54 | enum { IsRowMajor = SparseMatrix<Scalar,Opt2,StorageIndex>::IsRowMajor }; |
| 55 | sparseMat.setZero(); |
| 56 | //sparseMat.reserve(int(refMat.rows()*refMat.cols()*density)); |
| 57 | int nnz = static_cast<int>((1.5 * density) * static_cast<double>(IsRowMajor ? refMat.cols() : refMat.rows())); |
| 58 | sparseMat.reserve(VectorXi::Constant(IsRowMajor ? refMat.rows() : refMat.cols(), nnz)); |
| 59 | |
| 60 | Index insert_count = 0; |
| 61 | for(Index j=0; j<sparseMat.outerSize(); j++) |
| 62 | { |
| 63 | //sparseMat.startVec(j); |
| 64 | for(Index i=0; i<sparseMat.innerSize(); i++) |
| 65 | { |
| 66 | Index ai(i), aj(j); |
| 67 | if(IsRowMajor) |
| 68 | std::swap(ai,aj); |
| 69 | Scalar v = (internal::random<double>(0,1) < density) ? internal::random<Scalar>() : Scalar(0); |
| 70 | if ((flags&ForceNonZeroDiag) && (i==j)) |
| 71 | { |
| 72 | // FIXME: the following is too conservative |
| 73 | v = internal::random<Scalar>()*Scalar(3.); |
| 74 | v = v*v; |
| 75 | if(numext::real(v)>0) v += Scalar(5); |
| 76 | else v -= Scalar(5); |
| 77 | } |
| 78 | if ((flags & MakeLowerTriangular) && aj>ai) |
| 79 | v = Scalar(0); |
| 80 | else if ((flags & MakeUpperTriangular) && aj<ai) |
| 81 | v = Scalar(0); |
| 82 | |
| 83 | if ((flags&ForceRealDiag) && (i==j)) |
| 84 | v = numext::real(v); |
| 85 | |
| 86 | if (!numext::is_exactly_zero(v)) |
| 87 | { |
| 88 | //sparseMat.insertBackByOuterInner(j,i) = v; |
| 89 | sparseMat.insertByOuterInner(j,i) = v; |
| 90 | ++insert_count; |
| 91 | if (nonzeroCoords) |
| 92 | nonzeroCoords->push_back(Matrix<StorageIndex,2,1> (ai,aj)); |
| 93 | } |
| 94 | else if (zeroCoords) |
| 95 | { |
| 96 | zeroCoords->push_back(Matrix<StorageIndex,2,1> (ai,aj)); |
| 97 | } |
| 98 | refMat(ai,aj) = v; |
| 99 | |
| 100 | // make sure we only insert as many as the sparse matrix supports |
| 101 | if(insert_count == NumTraits<StorageIndex>::highest()) return; |
| 102 | } |
| 103 | } |