michalkielan / MatrixLibs

Small matrix library, written C++11
Mozilla Public License 2.0
1 stars 2 forks source link

LU for every kind of matrix #38

Open michalkielan opened 7 years ago

michalkielan commented 7 years ago

Current behaviour: the LU decomposition provided that all its leading submatrices have non-zero determinant, if any is equal to zero, exception is throwing as expected:

Matrix<float, 3, 3> A =
{ 1, 2, 3,
  2, 4, 5,
  1, 3, 4 };

Output:

Error: Divide by zero exception

In this particular example we can change the second and third row of matrix A, to make them lu decomposition friendly.

swap_rows(A, 1, 2);
auto LU = lu(A);
auto L = LU.first;
auto U = LU.second;
compare(A, L*U, 0.0001f) // return true

Please have a look at the full implementation of this example: https://github.com/michal915/MatrixLibs/blob/lu_every_matrix_tests/tst/lu/lu.cpp#L90

Expected behaviour: make LU decomposition works for every kind matrices

JaroslawWiosna commented 7 years ago

@michal915 have a look at #39


Sorry about goto: :( But I am not sure about changing rows. Are you sure, that after changing rows output L, U would fit equation A = L * U ?

michalkielan commented 7 years ago

Only in this particular example, otherwise, no idea :) i think more general solution would be needed, don't care about goto at the moment

JaroslawWiosna commented 7 years ago

BTW, @michal915 have a look at page 8

https://arxiv.org/pdf/math/0506382v1.pdf

JaroslawWiosna commented 7 years ago

@michal915 In my opinion, when det(A) eqals zero or if submatrices' det equals zero we cannot perform LU decomposition. Have a look at my comment above.

7332dd7 shouldn't be merged, because changing rows makes input matrix different that output matrix. It's like: "Well, I cannot sell you a car, but here is a bike, it is O.K. for you?" And for some people it would be fine, but in general when someone wants a car, he won't buy a bike.

JaroslawWiosna commented 7 years ago

http://stackoverflow.com/a/38156436

If B is a matrix obtained by swapping two rows of A, then

det(B) = -det(A)

    If you swap rows, flip the sign.

@michal915 So I did it in fcf4df0 but it failed. Do you know why?

EDIT: That's right approach. Here is your favourite 'octave test': click

JaroslawWiosna commented 7 years ago

@michal915 I've noticed that LU procedure failed even when input matrix can be decomposed using basing algorithm (without changing rows). 62c9a79