From aa6c516ec17fb44dff85b1abf3a1b05d58d3bc01 Mon Sep 17 00:00:00 2001 From: Gael Guennebaud Date: Mon, 16 Feb 2015 13:19:05 +0100 Subject: Fix many long to int conversion warnings: - fix usage of Index (API) versus StorageIndex (when multiple indexes are stored) - use StorageIndex(val) when the input has already been check - use internal::convert_index(val) when val is potentially unsafe (directly comes from user input) --- Eigen/src/SparseLU/SparseLU.h | 4 ++-- Eigen/src/SparseLU/SparseLUImpl.h | 2 +- Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h | 6 ++--- Eigen/src/SparseLU/SparseLU_Utils.h | 2 +- Eigen/src/SparseLU/SparseLU_column_bmod.h | 2 +- Eigen/src/SparseLU/SparseLU_column_dfs.h | 22 +++++++++--------- Eigen/src/SparseLU/SparseLU_copy_to_ucol.h | 2 +- Eigen/src/SparseLU/SparseLU_heap_relax_snode.h | 15 ++++++------ Eigen/src/SparseLU/SparseLU_panel_dfs.h | 32 +++++++++++++------------- Eigen/src/SparseLU/SparseLU_pivotL.h | 4 ++-- Eigen/src/SparseLU/SparseLU_pruneL.h | 2 +- Eigen/src/SparseLU/SparseLU_relax_snode.h | 8 +++---- 12 files changed, 50 insertions(+), 51 deletions(-) (limited to 'Eigen/src/SparseLU') diff --git a/Eigen/src/SparseLU/SparseLU.h b/Eigen/src/SparseLU/SparseLU.h index f60e4ac9d..1e448f2ab 100644 --- a/Eigen/src/SparseLU/SparseLU.h +++ b/Eigen/src/SparseLU/SparseLU.h @@ -397,7 +397,7 @@ void SparseLU::analyzePattern(const MatrixType& mat) if (!m_symmetricmode) { IndexVector post, iwork; // Post order etree - internal::treePostorder(m_mat.cols(), m_etree, post); + internal::treePostorder(StorageIndex(m_mat.cols()), m_etree, post); // Renumber etree in postorder @@ -479,7 +479,7 @@ void SparseLU::factorize(const MatrixType& matrix) else { //FIXME This should not be needed if the empty permutation is handled transparently m_perm_c.resize(matrix.cols()); - for(Index i = 0; i < matrix.cols(); ++i) m_perm_c.indices()(i) = i; + for(StorageIndex i = 0; i < matrix.cols(); ++i) m_perm_c.indices()(i) = i; } Index m = m_mat.rows(); diff --git a/Eigen/src/SparseLU/SparseLUImpl.h b/Eigen/src/SparseLU/SparseLUImpl.h index e735fd5c8..731d1652c 100644 --- a/Eigen/src/SparseLU/SparseLUImpl.h +++ b/Eigen/src/SparseLU/SparseLUImpl.h @@ -40,7 +40,7 @@ class SparseLUImpl Index snode_bmod (const Index jcol, const Index fsupc, ScalarVector& dense, GlobalLU_t& glu); Index pivotL(const Index jcol, const RealScalar& diagpivotthresh, IndexVector& perm_r, IndexVector& iperm_c, Index& pivrow, GlobalLU_t& glu); template - void dfs_kernel(const Index jj, IndexVector& perm_r, + void dfs_kernel(const StorageIndex jj, IndexVector& perm_r, Index& nseg, IndexVector& panel_lsub, IndexVector& segrep, Ref repfnz_col, IndexVector& xprune, Ref marker, IndexVector& parent, IndexVector& xplore, GlobalLU_t& glu, Index& nextl_col, Index krow, Traits& traits); diff --git a/Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h b/Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h index f7ffc2d9c..b37b93cf1 100644 --- a/Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h +++ b/Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h @@ -178,11 +178,11 @@ class MappedSuperNodalMatrix * \brief InnerIterator class to iterate over nonzero values of the current column in the supernodal matrix L * */ -template -class MappedSuperNodalMatrix::InnerIterator +template +class MappedSuperNodalMatrix::InnerIterator { public: - InnerIterator(const MappedSuperNodalMatrix& mat, Eigen::Index outer) + InnerIterator(const MappedSuperNodalMatrix& mat, Index outer) : m_matrix(mat), m_outer(outer), m_supno(mat.colToSup()[outer]), diff --git a/Eigen/src/SparseLU/SparseLU_Utils.h b/Eigen/src/SparseLU/SparseLU_Utils.h index b48157d9f..9e3dab44d 100644 --- a/Eigen/src/SparseLU/SparseLU_Utils.h +++ b/Eigen/src/SparseLU/SparseLU_Utils.h @@ -53,7 +53,7 @@ void SparseLUImpl::fixupL(const Index n, const IndexVector& { Index fsupc, i, j, k, jstart; - Index nextl = 0; + StorageIndex nextl = 0; Index nsuper = (glu.supno)(n); // For each supernode diff --git a/Eigen/src/SparseLU/SparseLU_column_bmod.h b/Eigen/src/SparseLU/SparseLU_column_bmod.h index bda01dcb3..be190997d 100644 --- a/Eigen/src/SparseLU/SparseLU_column_bmod.h +++ b/Eigen/src/SparseLU/SparseLU_column_bmod.h @@ -138,7 +138,7 @@ Index SparseLUImpl::column_bmod(const Index jcol, const Ind glu.lusup.segment(nextlu,offset).setZero(); nextlu += offset; } - glu.xlusup(jcol + 1) = nextlu; // close L\U(*,jcol); + glu.xlusup(jcol + 1) = StorageIndex(nextlu); // close L\U(*,jcol); /* For more updates within the panel (also within the current supernode), * should start from the first column of the panel, or the first column diff --git a/Eigen/src/SparseLU/SparseLU_column_dfs.h b/Eigen/src/SparseLU/SparseLU_column_dfs.h index 17c9e6adb..c98b30e32 100644 --- a/Eigen/src/SparseLU/SparseLU_column_dfs.h +++ b/Eigen/src/SparseLU/SparseLU_column_dfs.h @@ -112,13 +112,13 @@ Index SparseLUImpl::column_dfs(const Index m, const Index j // krow was visited before, go to the next nonz; if (kmark == jcol) continue; - dfs_kernel(jcol, perm_r, nseg, glu.lsub, segrep, repfnz, xprune, marker2, parent, + dfs_kernel(StorageIndex(jcol), perm_r, nseg, glu.lsub, segrep, repfnz, xprune, marker2, parent, xplore, glu, nextl, krow, traits); } // for each nonzero ... - Index fsupc, jptr, jm1ptr, ito, ifrom, istop; - Index nsuper = glu.supno(jcol); - Index jcolp1 = jcol + 1; + Index fsupc; + StorageIndex nsuper = glu.supno(jcol); + StorageIndex jcolp1 = StorageIndex(jcol) + 1; Index jcolm1 = jcol - 1; // check to see if j belongs in the same supernode as j-1 @@ -129,8 +129,8 @@ Index SparseLUImpl::column_dfs(const Index m, const Index j else { fsupc = glu.xsup(nsuper); - jptr = glu.xlsub(jcol); // Not yet compressed - jm1ptr = glu.xlsub(jcolm1); + StorageIndex jptr = glu.xlsub(jcol); // Not yet compressed + StorageIndex jm1ptr = glu.xlsub(jcolm1); // Use supernodes of type T2 : see SuperLU paper if ( (nextl-jptr != jptr-jm1ptr-1) ) jsuper = emptyIdxLU; @@ -148,13 +148,13 @@ Index SparseLUImpl::column_dfs(const Index m, const Index j { // starts a new supernode if ( (fsupc < jcolm1-1) ) { // >= 3 columns in nsuper - ito = glu.xlsub(fsupc+1); + StorageIndex ito = glu.xlsub(fsupc+1); glu.xlsub(jcolm1) = ito; - istop = ito + jptr - jm1ptr; + StorageIndex istop = ito + jptr - jm1ptr; xprune(jcolm1) = istop; // intialize xprune(jcol-1) glu.xlsub(jcol) = istop; - for (ifrom = jm1ptr; ifrom < nextl; ++ifrom, ++ito) + for (StorageIndex ifrom = jm1ptr; ifrom < nextl; ++ifrom, ++ito) glu.lsub(ito) = glu.lsub(ifrom); nextl = ito; // = istop + length(jcol) } @@ -166,8 +166,8 @@ Index SparseLUImpl::column_dfs(const Index m, const Index j // Tidy up the pointers before exit glu.xsup(nsuper+1) = jcolp1; glu.supno(jcolp1) = nsuper; - xprune(jcol) = nextl; // Intialize upper bound for pruning - glu.xlsub(jcolp1) = nextl; + xprune(jcol) = StorageIndex(nextl); // Intialize upper bound for pruning + glu.xlsub(jcolp1) = StorageIndex(nextl); return 0; } diff --git a/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h b/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h index bf237951d..c32d8d8b1 100644 --- a/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h +++ b/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h @@ -56,7 +56,7 @@ Index SparseLUImpl::copy_to_ucol(const Index jcol, const In // For each nonzero supernode segment of U[*,j] in topological order Index k = nseg - 1, i; - Index nextu = glu.xusub(jcol); + StorageIndex nextu = glu.xusub(jcol); Index kfnz, isub, segsize; Index new_next,irow; Index fsupc, mem; diff --git a/Eigen/src/SparseLU/SparseLU_heap_relax_snode.h b/Eigen/src/SparseLU/SparseLU_heap_relax_snode.h index 4092f842f..6f75d500e 100644 --- a/Eigen/src/SparseLU/SparseLU_heap_relax_snode.h +++ b/Eigen/src/SparseLU/SparseLU_heap_relax_snode.h @@ -48,15 +48,14 @@ void SparseLUImpl::heap_relax_snode (const Index n, IndexVe // The etree may not be postordered, but its heap ordered IndexVector post; - internal::treePostorder(n, et, post); // Post order etree + internal::treePostorder(StorageIndex(n), et, post); // Post order etree IndexVector inv_post(n+1); - Index i; - for (i = 0; i < n+1; ++i) inv_post(post(i)) = i; // inv_post = post.inverse()??? + for (StorageIndex i = 0; i < n+1; ++i) inv_post(post(i)) = i; // inv_post = post.inverse()??? // Renumber etree in postorder IndexVector iwork(n); IndexVector et_save(n+1); - for (i = 0; i < n; ++i) + for (Index i = 0; i < n; ++i) { iwork(post(i)) = post(et(i)); } @@ -78,7 +77,7 @@ void SparseLUImpl::heap_relax_snode (const Index n, IndexVe StorageIndex k; Index nsuper_et_post = 0; // Number of relaxed snodes in postordered etree Index nsuper_et = 0; // Number of relaxed snodes in the original etree - Index l; + StorageIndex l; for (j = 0; j < n; ) { parent = et(j); @@ -90,8 +89,8 @@ void SparseLUImpl::heap_relax_snode (const Index n, IndexVe } // Found a supernode in postordered etree, j is the last column ++nsuper_et_post; - k = n; - for (i = snode_start; i <= j; ++i) + k = StorageIndex(n); + for (Index i = snode_start; i <= j; ++i) k = (std::min)(k, inv_post(i)); l = inv_post(j); if ( (l - k) == (j - snode_start) ) // Same number of columns in the snode @@ -102,7 +101,7 @@ void SparseLUImpl::heap_relax_snode (const Index n, IndexVe } else { - for (i = snode_start; i <= j; ++i) + for (Index i = snode_start; i <= j; ++i) { l = inv_post(i); if (descendants(i) == 0) diff --git a/Eigen/src/SparseLU/SparseLU_panel_dfs.h b/Eigen/src/SparseLU/SparseLU_panel_dfs.h index f4a908ee5..155df7336 100644 --- a/Eigen/src/SparseLU/SparseLU_panel_dfs.h +++ b/Eigen/src/SparseLU/SparseLU_panel_dfs.h @@ -41,7 +41,7 @@ struct panel_dfs_traits panel_dfs_traits(Index jcol, StorageIndex* marker) : m_jcol(jcol), m_marker(marker) {} - bool update_segrep(Index krep, Index jj) + bool update_segrep(Index krep, StorageIndex jj) { if(m_marker[krep] template -void SparseLUImpl::dfs_kernel(const Index jj, IndexVector& perm_r, +void SparseLUImpl::dfs_kernel(const StorageIndex jj, IndexVector& perm_r, Index& nseg, IndexVector& panel_lsub, IndexVector& segrep, Ref repfnz_col, IndexVector& xprune, Ref marker, IndexVector& parent, IndexVector& xplore, GlobalLU_t& glu, @@ -67,14 +67,14 @@ void SparseLUImpl::dfs_kernel(const Index jj, IndexVector& ) { - Index kmark = marker(krow); + StorageIndex kmark = marker(krow); // For each unmarked krow of jj marker(krow) = jj; - Index kperm = perm_r(krow); + StorageIndex kperm = perm_r(krow); if (kperm == emptyIdxLU ) { // krow is in L : place it in structure of L(*, jj) - panel_lsub(nextl_col++) = krow; // krow is indexed into A + panel_lsub(nextl_col++) = StorageIndex(krow); // krow is indexed into A traits.mem_expand(panel_lsub, nextl_col, kmark); } @@ -83,9 +83,9 @@ void SparseLUImpl::dfs_kernel(const Index jj, IndexVector& // krow is in U : if its supernode-representative krep // has been explored, update repfnz(*) // krep = supernode representative of the current row - Index krep = glu.xsup(glu.supno(kperm)+1) - 1; + StorageIndex krep = glu.xsup(glu.supno(kperm)+1) - 1; // First nonzero element in the current column: - Index myfnz = repfnz_col(krep); + StorageIndex myfnz = repfnz_col(krep); if (myfnz != emptyIdxLU ) { @@ -96,26 +96,26 @@ void SparseLUImpl::dfs_kernel(const Index jj, IndexVector& else { // Otherwise, perform dfs starting at krep - Index oldrep = emptyIdxLU; + StorageIndex oldrep = emptyIdxLU; parent(krep) = oldrep; repfnz_col(krep) = kperm; - Index xdfs = glu.xlsub(krep); + StorageIndex xdfs = glu.xlsub(krep); Index maxdfs = xprune(krep); - Index kpar; + StorageIndex kpar; do { // For each unmarked kchild of krep while (xdfs < maxdfs) { - Index kchild = glu.lsub(xdfs); + StorageIndex kchild = glu.lsub(xdfs); xdfs++; - Index chmark = marker(kchild); + StorageIndex chmark = marker(kchild); if (chmark != jj ) { marker(kchild) = jj; - Index chperm = perm_r(kchild); + StorageIndex chperm = perm_r(kchild); if (chperm == emptyIdxLU) { @@ -128,7 +128,7 @@ void SparseLUImpl::dfs_kernel(const Index jj, IndexVector& // case kchild is in U : // chrep = its supernode-rep. If its rep has been explored, // update its repfnz(*) - Index chrep = glu.xsup(glu.supno(chperm)+1) - 1; + StorageIndex chrep = glu.xsup(glu.supno(chperm)+1) - 1; myfnz = repfnz_col(chrep); if (myfnz != emptyIdxLU) @@ -227,7 +227,7 @@ void SparseLUImpl::panel_dfs(const Index m, const Index w, panel_dfs_traits traits(jcol, marker1.data()); // For each column in the panel - for (Index jj = jcol; jj < jcol + w; jj++) + for (StorageIndex jj = StorageIndex(jcol); jj < jcol + w; jj++) { nextl_col = (jj - jcol) * m; @@ -241,7 +241,7 @@ void SparseLUImpl::panel_dfs(const Index m, const Index w, Index krow = it.row(); dense_col(krow) = it.value(); - Index kmark = marker(krow); + StorageIndex kmark = marker(krow); if (kmark == jj) continue; // krow visited before, go to the next nonzero diff --git a/Eigen/src/SparseLU/SparseLU_pivotL.h b/Eigen/src/SparseLU/SparseLU_pivotL.h index 01f5ba4e9..562128b69 100644 --- a/Eigen/src/SparseLU/SparseLU_pivotL.h +++ b/Eigen/src/SparseLU/SparseLU_pivotL.h @@ -89,7 +89,7 @@ Index SparseLUImpl::pivotL(const Index jcol, const RealScal // Test for singularity if ( pivmax == 0.0 ) { pivrow = lsub_ptr[pivptr]; - perm_r(pivrow) = jcol; + perm_r(pivrow) = StorageIndex(jcol); return (jcol+1); } @@ -110,7 +110,7 @@ Index SparseLUImpl::pivotL(const Index jcol, const RealScal } // Record pivot row - perm_r(pivrow) = jcol; + perm_r(pivrow) = StorageIndex(jcol); // Interchange row subscripts if (pivptr != nsupc ) { diff --git a/Eigen/src/SparseLU/SparseLU_pruneL.h b/Eigen/src/SparseLU/SparseLU_pruneL.h index 13133fcc2..ad32fed5e 100644 --- a/Eigen/src/SparseLU/SparseLU_pruneL.h +++ b/Eigen/src/SparseLU/SparseLU_pruneL.h @@ -124,7 +124,7 @@ void SparseLUImpl::pruneL(const Index jcol, const IndexVect } } // end while - xprune(irep) = kmin; //Pruning + xprune(irep) = StorageIndex(kmin); //Pruning } // end if do_prune } // end pruning } // End for each U-segment diff --git a/Eigen/src/SparseLU/SparseLU_relax_snode.h b/Eigen/src/SparseLU/SparseLU_relax_snode.h index 21c182d56..c408d01b4 100644 --- a/Eigen/src/SparseLU/SparseLU_relax_snode.h +++ b/Eigen/src/SparseLU/SparseLU_relax_snode.h @@ -48,10 +48,10 @@ void SparseLUImpl::relax_snode (const Index n, IndexVector& { // compute the number of descendants of each node in the etree - Index j, parent; + Index parent; relax_end.setConstant(emptyIdxLU); descendants.setZero(); - for (j = 0; j < n; j++) + for (Index j = 0; j < n; j++) { parent = et(j); if (parent != n) // not the dummy root @@ -59,7 +59,7 @@ void SparseLUImpl::relax_snode (const Index n, IndexVector& } // Identify the relaxed supernodes by postorder traversal of the etree Index snode_start; // beginning of a snode - for (j = 0; j < n; ) + for (Index j = 0; j < n; ) { parent = et(j); snode_start = j; @@ -69,7 +69,7 @@ void SparseLUImpl::relax_snode (const Index n, IndexVector& parent = et(j); } // Found a supernode in postordered etree, j is the last column - relax_end(snode_start) = j; // Record last column + relax_end(snode_start) = StorageIndex(j); // Record last column j++; // Search for a new leaf while (descendants(j) != 0 && j < n) j++; -- cgit v1.2.3