mirror of
https://gitlab.com/libeigen/eigen.git
synced 2025-08-14 04:35:57 +08:00
Correct bug for triangular solve within supernodes
This commit is contained in:
parent
b5a83867ca
commit
3095e4a5f9
@ -186,7 +186,7 @@ class SparseLU
|
|||||||
// Triangular solve
|
// Triangular solve
|
||||||
Map<const Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > A( &(Lval[luptr]), nsupc, nsupc, OuterStride<>(nsupr) );
|
Map<const Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > A( &(Lval[luptr]), nsupc, nsupc, OuterStride<>(nsupr) );
|
||||||
Map< Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > U (&(X.data()[fsupc]), nsupc, nrhs, OuterStride<>(X.rows()) );
|
Map< Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > U (&(X.data()[fsupc]), nsupc, nrhs, OuterStride<>(X.rows()) );
|
||||||
U = A.template triangularView<Lower>().solve(U);
|
U = A.template triangularView<UnitLower>().solve(U);
|
||||||
|
|
||||||
// Matrix-vector product
|
// Matrix-vector product
|
||||||
new (&A) Map<const Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > ( &(Lval[luptr+nsupc]), nrow, nsupc, OuterStride<>(nsupr) );
|
new (&A) Map<const Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > ( &(Lval[luptr+nsupc]), nrow, nsupc, OuterStride<>(nsupr) );
|
||||||
|
@ -137,7 +137,7 @@ int LU_column_bmod(const int jcol, const int nseg, BlockScalarVector& dense, Sca
|
|||||||
Map<Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > A( &(lusup.data()[luptr]), segsize, segsize, OuterStride<>(nsupr) );
|
Map<Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > A( &(lusup.data()[luptr]), segsize, segsize, OuterStride<>(nsupr) );
|
||||||
VectorBlock<ScalarVector> u(tempv, 0, segsize);
|
VectorBlock<ScalarVector> u(tempv, 0, segsize);
|
||||||
|
|
||||||
u = A.template triangularView<Lower>().solve(u);
|
u = A.template triangularView<UnitLower>().solve(u);
|
||||||
|
|
||||||
// Dense matrix-vector product y <-- A*x
|
// Dense matrix-vector product y <-- A*x
|
||||||
luptr += segsize;
|
luptr += segsize;
|
||||||
@ -213,7 +213,7 @@ int LU_column_bmod(const int jcol, const int nseg, BlockScalarVector& dense, Sca
|
|||||||
ufirst = xlusup(jcol) + d_fsupc;
|
ufirst = xlusup(jcol) + d_fsupc;
|
||||||
Map<Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > A( &(lusup.data()[luptr]), nsupc, nsupc, OuterStride<>(nsupr) );
|
Map<Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > A( &(lusup.data()[luptr]), nsupc, nsupc, OuterStride<>(nsupr) );
|
||||||
VectorBlock<ScalarVector> u(lusup, ufirst, nsupc);
|
VectorBlock<ScalarVector> u(lusup, ufirst, nsupc);
|
||||||
u = A.template triangularView<Lower>().solve(u);
|
u = A.template triangularView<UnitLower>().solve(u);
|
||||||
|
|
||||||
new (&A) Map<Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > ( &(lusup.data()[luptr+nsupc]), nrow, nsupc, OuterStride<>(nsupr) );
|
new (&A) Map<Matrix<Scalar,Dynamic,Dynamic>, 0, OuterStride<> > ( &(lusup.data()[luptr+nsupc]), nrow, nsupc, OuterStride<>(nsupr) );
|
||||||
VectorBlock<ScalarVector> l(lusup, ufirst+nsupc, nrow);
|
VectorBlock<ScalarVector> l(lusup, ufirst+nsupc, nrow);
|
||||||
|
@ -132,8 +132,9 @@ void LU_panel_bmod(const int m, const int w, const int jcol, const int nseg, Sca
|
|||||||
luptr += nsupr * no_zeros + no_zeros;
|
luptr += nsupr * no_zeros + no_zeros;
|
||||||
// triangular solve with Eigen
|
// triangular solve with Eigen
|
||||||
Map<Matrix<Scalar,Dynamic, Dynamic>, 0, OuterStride<> > A( &(lusup.data()[luptr]), segsize, segsize, OuterStride<>(nsupr) );
|
Map<Matrix<Scalar,Dynamic, Dynamic>, 0, OuterStride<> > A( &(lusup.data()[luptr]), segsize, segsize, OuterStride<>(nsupr) );
|
||||||
|
std::cout<< " Matrix \n" << A << std::endl;
|
||||||
VectorBlock<ScalarVector> u(tempv, 0, segsize);
|
VectorBlock<ScalarVector> u(tempv, 0, segsize);
|
||||||
u = A.template triangularView<Lower>().solve(u);
|
u = A.template triangularView<UnitLower>().solve(u);
|
||||||
|
|
||||||
luptr += segsize;
|
luptr += segsize;
|
||||||
// Dense Matrix vector product y <-- A*x;
|
// Dense Matrix vector product y <-- A*x;
|
||||||
@ -164,7 +165,7 @@ void LU_panel_bmod(const int m, const int w, const int jcol, const int nseg, Sca
|
|||||||
l(i) = Scalar(0);
|
l(i) = Scalar(0);
|
||||||
++isub;
|
++isub;
|
||||||
}
|
}
|
||||||
|
std::cout<< jj << " : " << dense_col.transpose() << std::endl;
|
||||||
} // End for each column in the panel
|
} // End for each column in the panel
|
||||||
|
|
||||||
} // End for each updating supernode
|
} // End for each updating supernode
|
||||||
|
@ -76,7 +76,7 @@ int LU_snode_bmod (const int jcol, const int fsupc, ScalarVector& dense, LU_Glob
|
|||||||
// Solve the triangular system for U(fsupc:jcol, jcol) with L(fspuc:jcol, fsupc:jcol)
|
// Solve the triangular system for U(fsupc:jcol, jcol) with L(fspuc:jcol, fsupc:jcol)
|
||||||
Map<Matrix<Scalar,Dynamic,Dynamic>,0,OuterStride<> > A( &(lusup.data()[luptr]), nsupc, nsupc, OuterStride<>(nsupr) );
|
Map<Matrix<Scalar,Dynamic,Dynamic>,0,OuterStride<> > A( &(lusup.data()[luptr]), nsupc, nsupc, OuterStride<>(nsupr) );
|
||||||
VectorBlock<ScalarVector> u(lusup, ufirst, nsupc);
|
VectorBlock<ScalarVector> u(lusup, ufirst, nsupc);
|
||||||
u = A.template triangularView<Lower>().solve(u); // Call the Eigen dense triangular solve interface
|
u = A.template triangularView<UnitLower>().solve(u); // Call the Eigen dense triangular solve interface
|
||||||
|
|
||||||
// Update the trailing part of the column jcol U(jcol:jcol+nrow, jcol) using L(jcol:jcol+nrow, fsupc:jcol) and U(fsupc:jcol)
|
// Update the trailing part of the column jcol U(jcol:jcol+nrow, jcol) using L(jcol:jcol+nrow, fsupc:jcol) and U(fsupc:jcol)
|
||||||
new (&A) Map<Matrix<Scalar,Dynamic,Dynamic>,0,OuterStride<> > ( &(lusup.data()[luptr+nsupc]), nrow, nsupc, OuterStride<>(nsupr) );
|
new (&A) Map<Matrix<Scalar,Dynamic,Dynamic>,0,OuterStride<> > ( &(lusup.data()[luptr+nsupc]), nrow, nsupc, OuterStride<>(nsupr) );
|
||||||
|
Loading…
x
Reference in New Issue
Block a user