HarmonicSmoothing code uses CG for non-symmetric matrix
Eigen backend fails for HarmonicSmoothing test because HarmonicSmoothing code uses CG for non-symmetric matrix.
Comments (11)
-
-
reporter It fails because Eigen only uses the lower or upper part of the matrix.
Using CG for non-symmetric matrices is unsafe because backends may exploit symmetry for performance reasons.
-
Pushed some fixes to branch:jan/feature-424-eigen-backend.
I'd rather say that such backends exploit wrong assumptions. CG is in this case mathematically correct and PETSc respects it. The same trick is actually employed in
navier-stokes
demo; it breaks also with Eigen CG. Eigen should haveWhole
apart fromUpper
andLower
. -
I was looking today if I can propose a fix to Eigen. Little problem is that I cannot compile DOLFIN with Eigen dev version
In file included from /usr/users/blechta/fenics/eigen/Eigen/SparseCore:37:0, from /usr/users/blechta/fenics/eigen/Eigen/Sparse:19, from /usr/users/blechta/fenics/fenics-dev-dbg/src/dolfin/dolfin/la/EigenMatrix.h:31, from /usr/users/blechta/fenics/fenics-dev-dbg/src/dolfin/dolfin/la/EigenMatrix.cpp:18: /usr/users/blechta/fenics/eigen/Eigen/src/SparseCore/SparseCwiseUnaryOp.h: In member function 'Derived& Eigen::SparseMatrixBase<Derived>::operator*=(const Scalar&) [with Derived = Eigen::Block<Eigen::SparseMatrix<double, 1>, 1, -0x00000000000000001, true>, Eigen::SparseMatrixBase<Derived>::Scalar = double]': /usr/users/blechta/fenics/fenics-dev-dbg/src/dolfin/dolfin/la/EigenMatrix.cpp:229:24: instantiated from here /usr/users/blechta/fenics/eigen/Eigen/src/SparseCore/SparseCwiseUnaryOp.h:169:42: error: no type named 'InnerIterator' in 'class Eigen::Block<Eigen::SparseMatrix<double, 1>, 1, -0x00000000000000001, true>' /usr/users/blechta/fenics/eigen/Eigen/src/SparseCore/SparseCwiseUnaryOp.h:169:42: error: no type named 'InnerIterator' in 'class Eigen::Block<Eigen::SparseMatrix<double, 1>, 1, -0x00000000000000001, true>' /usr/users/blechta/fenics/eigen/Eigen/src/SparseCore/SparseCwiseUnaryOp.h:169:42: error: no type named 'InnerIterator' in 'class Eigen::Block<Eigen::SparseMatrix<double, 1>, 1, -0x00000000000000001, true>' /usr/users/blechta/fenics/eigen/Eigen/src/SparseCore/SparseCwiseUnaryOp.h:169:42: error: no type named 'InnerIterator' in 'class Eigen::Block<Eigen::SparseMatrix<double, 1>, 1, -0x00000000000000001, true>'
@garth-wells, Do you know whether this "non-symmetric CG" trick works with TPetra? Do
HarmonicSmoothing
unit test andnavier-stokes
demo pass? -
reporter @blechta I don't know if it works with Epetra. I don't like to the 'trick' - it's bound to cause problems sooner or later.
I don't know what's up with the Eigen dev version. @chris_richardson might know.
-
Fix proposed upstream. DOLFIN fixes in jan/feature-424-eigen-backend.
-
@blechta - OK, let's make a note in the code, so this can be adopted when the fix is widespread.
-
- changed title to HarmonicSmoothing code uses CG for non-symmetric matrix
-
So far Eigen guys seem to be ignorant of mathematics and being reluctant to merge the fix in.
-
reporter - changed status to resolved
Switched to bicgstab for now.
-
- removed milestone
Removing milestone: 1.6 (automated comment)
- Log in to comment
CG is justifiable here because right-hand side is in the image of symmetric part of the operator and initial approximation is in the kernel of non-symmetric part of the operator. This may fail at least from two reasons:
I'll look into it.