Invalid Results Using Linear System Solvers Eigen 3
I am porting my MATLAB code to C ++ from Eigen 3 and I decided to build my linear solver rather than call it from a matrix object so I can reuse it. Unfortunately, it does not give the expected results. After several tests, I traced the issue with what looks like a linear solver object, as evidenced by the following relevant snippet of my code:
MatrixXcd M(6,6), Y(6,6), dvalor(6,6);
// Initialization of the matrices...
FullPivHouseholderQR<MatrixXcd> mdivide;
auto result = mldivide.compute(Y).solve(M).eval() *
dvalor * mldivide.compute(M).solve(Y).eval();
std::cout << "result:" << std::endl
<< result << std::endl << std::endl;
std::cout << "using fullPivHouseholderQr:" << std::endl
<< Y.fullPivHouseholderQr().solve(M) * dvalor *
M.fullPivHouseholderQr().solve(Y) << std::endl << std::endl;
std::cout << "using jacobi:" << std::endl
<< Y.jacobiSvd(ComputeThinU | ComputeThinV).solve(M) * dvalor *
M.jacobiSvd(ComputeThinU | ComputeThinV).solve(Y) << std::endl << std::endl;
std::cout << "using inverse:" << std::endl
<< Y.inverse() * M * dvalor * M.inverse() * Y << std::endl << std::endl;
This is the result:
result:
(0.564196,0.00606298) (-0.15282,-0.00179293) (-0.16564,-0.000220726) (-0.16564,-0.000220726) (-0.15282,-0.00179293) (-0.179235,0.00541725)
(2.73563e+184,-8.89017e+185) (1,1.12229e-10) (-6.04297e-11,5.62259e-12) (-6.02256e-11,2.60221e-12) (-5.33448e-11,-3.05694e-11) (-6.55427e-11,8.17723e-12)
(2.11057e+184,-8.81209e+185) (-5.53615e-11,4.54731e-12) (1,1.11384e-10) (-5.99088e-11,9.29495e-12) (-5.51744e-11,1.52693e-12) (-6.21597e-11,-3.11176e-11)
(2.11057e+184,-8.81209e+185) (-0.171319,-0.00113681) (0.134093,-2.05541e-05) (0.865907,2.05542e-05) (0.171319,0.00113681) (-0.140389,-0.000216684)
(2.73563e+184,-8.89017e+185) (0.218904,0.00143522) (-0.171337,3.98397e-05) (0.171337,-3.98397e-05) (0.781096,-0.00143522) (0.179383,0.000262654)
(1.4668e+184,-8.78677e+185) (-5.52732e-11,5.99039e-12) (-5.69133e-11,-3.2223e-11) (-5.97328e-11,9.63533e-12) (-5.51767e-11,4.01574e-12) (1,1.1063e-10)
using fullPivHouseholderQr:
(1,1.1063e-10) (-5.51766e-11,4.01578e-12) (-5.97331e-11,9.63535e-12) (-5.6913e-11,-3.22229e-11) (-5.52733e-11,5.99029e-12) (-6.46967e-11,1.0973e-11)
(-6.54343e-11,6.20279e-12) (1,1.12229e-10) (-6.04297e-11,5.62259e-12) (-6.02256e-11,2.60221e-12) (-5.33448e-11,-3.05694e-11) (-6.55427e-11,8.17723e-12)
(-6.48883e-11,1.07408e-11) (-5.53615e-11,4.54731e-12) (1,1.11384e-10) (-5.99088e-11,9.29495e-12) (-5.51744e-11,1.52693e-12) (-6.21597e-11,-3.11176e-11)
(-6.21592e-11,-3.11175e-11) (-5.51743e-11,1.52697e-12) (-5.99093e-11,9.29482e-12) (1,1.11385e-10) (-5.53614e-11,4.54728e-12) (-6.48881e-11,1.07408e-11)
(-6.55427e-11,8.1773e-12) (-5.33444e-11,-3.05693e-11) (-6.02261e-11,2.60223e-12) (-6.04293e-11,5.62257e-12) (1,1.12229e-10) (-6.54344e-11,6.20272e-12)
(-6.46967e-11,1.0973e-11) (-5.52732e-11,5.99039e-12) (-5.69133e-11,-3.2223e-11) (-5.97328e-11,9.63533e-12) (-5.51767e-11,4.01574e-12) (1,1.1063e-10)
using jacobi:
(1,1.1063e-10) (-5.51761e-11,4.01562e-12) (-5.9733e-11,9.63488e-12) (-5.69134e-11,-3.22232e-11) (-5.5273e-11,5.9903e-12) (-6.46969e-11,1.09735e-11)
(-6.54346e-11,6.20257e-12) (1,1.12228e-10) (-6.04303e-11,5.62232e-12) (-6.0225e-11,2.60224e-12) (-5.33448e-11,-3.05694e-11) (-6.55414e-11,8.17752e-12)
(-6.48874e-11,1.07408e-11) (-5.53617e-11,4.54725e-12) (1,1.11384e-10) (-5.99097e-11,9.29465e-12) (-5.51742e-11,1.52688e-12) (-6.21596e-11,-3.11175e-11)
(-6.21608e-11,-3.11175e-11) (-5.5174e-11,1.52688e-12) (-5.99104e-11,9.29433e-12) (1,1.11384e-10) (-5.53612e-11,4.54728e-12) (-6.48873e-11,1.07413e-11)
(-6.55423e-11,8.17718e-12) (-5.33459e-11,-3.05695e-11) (-6.02258e-11,2.60206e-12) (-6.04294e-11,5.6224e-12) (1,1.12229e-10) (-6.5434e-11,6.20302e-12)
(-6.46966e-11,1.09729e-11) (-5.5273e-11,5.99037e-12) (-5.69134e-11,-3.22233e-11) (-5.97333e-11,9.63543e-12) (-5.51767e-11,4.01582e-12) (1,1.1063e-10)
using inverse:
(1,1.1063e-10) (-5.51768e-11,4.01574e-12) (-5.97329e-11,9.63535e-12) (-5.6914e-11,-3.2223e-11) (-5.52733e-11,5.99031e-12) (-6.46967e-11,1.0973e-11)
(-6.54338e-11,6.20274e-12) (1,1.12229e-10) (-6.04293e-11,5.6225e-12) (-6.02264e-11,2.6022e-12) (-5.33451e-11,-3.05693e-11) (-6.55425e-11,8.1773e-12)
(-6.48877e-11,1.07408e-11) (-5.53615e-11,4.54726e-12) (1,1.11384e-10) (-5.99097e-11,9.29488e-12) (-5.51743e-11,1.52697e-12) (-6.21592e-11,-3.11176e-11)
(-6.21587e-11,-3.11176e-11) (-5.51744e-11,1.52696e-12) (-5.99093e-11,9.29489e-12) (1,1.11384e-10) (-5.53616e-11,4.54727e-12) (-6.4888e-11,1.07408e-11)
(-6.55419e-11,8.17731e-12) (-5.33448e-11,-3.05693e-11) (-6.02257e-11,2.6022e-12) (-6.04301e-11,5.6225e-12) (1,1.12229e-10) (-6.54343e-11,6.20273e-12)
(-6.46964e-11,1.0973e-11) (-5.52732e-11,5.99031e-12) (-5.69131e-11,-3.2223e-11) (-5.97336e-11,9.63535e-12) (-5.51767e-11,4.01574e-12) (1,1.1063e-10)
The first iteration of the algorithm in MATLAB creates an almost identical matrix, which can be seen in C ++ by accessing matrices and calling the solver from the matrix object.
But the result is completely wrong when using the solver object. I have also tried other solvers from my own documentation page with the same result: https://eigen.tuxfamily.org/dox/group__TutorialLinearAlgebra.html .
Am I missing any step or will it call the solver from the matrix, else will it initialize the solver?
The native version is 3.3.4 and I am not using compiler optimizers yet, just C ++ 11:
g++ -I. -std=c++11 main.cpp
The size of the matrices is constant and it will run inside a large loop. Thus, the more memory allocation I can avoid by pulling out the solver, the better.
Thank.
source to share
Be careful with auto
, see common pitfalls , mainly in your case result
is an expression Product<>
referencing dead objects, Replace auto
with MatrixXcd
and solve the problem.
source to share