30 #ifndef SPARSELU_PIVOTL_H 31 #define SPARSELU_PIVOTL_H 59 template <
typename Scalar,
typename Index>
60 Index SparseLUImpl<Scalar,Index>::pivotL(
const Index jcol,
const RealScalar& diagpivotthresh, IndexVector& perm_r, IndexVector& iperm_c, Index& pivrow, GlobalLU_t& glu)
63 Index fsupc = (glu.xsup)((glu.supno)(jcol));
64 Index nsupc = jcol - fsupc;
65 Index lptr = glu.xlsub(fsupc);
66 Index nsupr = glu.xlsub(fsupc+1) - lptr;
67 Index lda = glu.xlusup(fsupc+1) - glu.xlusup(fsupc);
68 Scalar* lu_sup_ptr = &(glu.lusup.data()[glu.xlusup(fsupc)]);
69 Scalar* lu_col_ptr = &(glu.lusup.data()[glu.xlusup(jcol)]);
70 Index* lsub_ptr = &(glu.lsub.data()[lptr]);
73 Index diagind = iperm_c(jcol);
74 RealScalar pivmax = 0.0;
76 Index diag = emptyIdxLU;
78 Index isub, icol, itemp, k;
79 for (isub = nsupc; isub < nsupr; ++isub) {
81 rtemp = abs(lu_col_ptr[isub]);
86 if (lsub_ptr[isub] == diagind) diag = isub;
90 if ( pivmax == 0.0 ) {
91 pivrow = lsub_ptr[pivptr];
92 perm_r(pivrow) = jcol;
96 RealScalar thresh = diagpivotthresh * pivmax;
106 rtemp = abs(lu_col_ptr[diag]);
107 if (rtemp != 0.0 && rtemp >= thresh) pivptr = diag;
109 pivrow = lsub_ptr[pivptr];
113 perm_r(pivrow) = jcol;
115 if (pivptr != nsupc )
117 std::swap( lsub_ptr[pivptr], lsub_ptr[nsupc] );
120 for (icol = 0; icol <= nsupc; icol++)
122 itemp = pivptr + icol * lda;
123 std::swap(lu_sup_ptr[itemp], lu_sup_ptr[nsupc + icol * lda]);
127 Scalar temp = Scalar(1.0) / lu_col_ptr[nsupc];
128 for (k = nsupc+1; k < nsupr; k++)
129 lu_col_ptr[k] *= temp;
136 #endif // SPARSELU_PIVOTL_H
Definition: Eigen_Colamd.h:54