Re: [eigen] SVD Bug

[Benoit Jacob <jacob.benoit.1@xxxxxxxxx>]

2010/9/24 Hauke Heibel <hauke.heibel@xxxxxxxxxxxxxx>:
> Just to help me understand what is going on...
>
> Does that mean for a 3x2 the decomposition does actually create
> matrices with dimensions
>
> 3x2 2x2 2x2
>
> and not
>
> 3x3 3x2 2x2 ?

Our SVD creates a matrixU of size 3x3 with the last column filled with
0's (yes that's quite stupid, but again it's not worth breaking
existing behavior in matrixU() which will be deprecated once we have a
cleverU()... ),  a singularValues vector of size 2, and a matrixV of
size 2x2.

>
> I am just asking since I think at an earlier point we argued that for
> USV^T, we always wanted U and V to be square an unitary. I am not sure
> and did not yet test it but I think Eigen 2 actually does that and U
> contains no zeros.

ah! I didn't remember that.

OK, that's a good reason to keep filling U so it's actually unitary.
At least that matrixU() method would remain potentially useful in a
future implementation (having a good mathematical meaning) even if
slow.

To implement that, we need to do Gram-Schmidt to the remaining columns.

of course we need a more efficient SVD with compact U and V anyway....


>
> Let me clarify what you meant with "does not have a solution". In
> fact, what I did in the unit test was sort of stupid because solving a
> system A*x = b where A is overdetermined will of course result in an x
> = pinv(A)*b such that Ax != b but where ||Ax-b||^2 is minimal. So in
> case we can guarantee that our SVD produces a solution that computes
> exactly the x which minimized ||Ax-b||, we can say it is fine though I
> don't know how to unit test such a system besides doing tests based on
> problems for which we know the correct solutions.

Yes, so I was forgetting that in the case of SVD, solve() is
implicitly doing least squares. I agree with all you write above.

Benoit

>
> - Hauke
>
> On Thu, Sep 23, 2010 at 3:58 PM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx> wrote:
>> It turns out that the unit test was failing because it was trying to
>> solve a system that didn't have solution. It took a random rectangular
>> matrix and a random rhs and tried to solve. That was working well on
>> square matrices, but broke on rectangular matrices with rows >= cols.
>>
>> I pushed your patch plus one more changeset fixing the unit test, and
>> fixing the SVD to actually enforce rows >= cols since our code is
>> still relying on that.
>>
>> Also I kept the existing U unchanged in all its stupidity (mxm filled
>> with zeros in the useless area) because it's not worth breaking
>> existing behavior at this point; once SVD is reimplemented the
>> matrixU() method will be deprecated anyway and replaced by something
>> that evaluates only what's needed.
>>
>> In conclusion, problem solved, thanks for your memory error fix.
>>
>> Benoit
>>
>> 2010/9/23 Benoit Jacob <jacob.benoit.1@xxxxxxxxx>:
>>> OK so what happens is that this SVD implementation computes a compact U:
>>>
>>>  m_matU.setZero();
>>>  if (m>=n)
>>>    m_matU.block(0,0,m,n) = A;
>>>  else
>>>    m_matU = A.block(0,0,m,m);
>>>
>>> it doesn't bother extending the remaining rows by orthonormalization.
>>> That actually isn't a bad thing. So the decomposition itself is fine.
>>>  We need to:
>>>  - document that aspect of U when taking the SVD of a rectangular matrix
>>>  - fix the solve method.
>>>
>>> Benoit
>>>
>>> 2010/9/23 Benoit Jacob <jacob.benoit.1@xxxxxxxxx>:
>>>> 2010/9/23 Hauke Heibel <hauke.heibel@xxxxxxxxxxxxxx>:
>>>>> I tried to look into it but failed. I can only say that already the
>>>>> matrix U is strange since its last row is zero and this should not be
>>>>> the case.
>>>>
>>>> Yes, indeed, that's very interesting:
>>>>
>>>> a
>>>> -0.08936    0.971
>>>> -0.5195  -0.1061
>>>>  0.3311  -0.4731
>>>>
>>>> matrixU
>>>> -0.8752  -0.2161        0
>>>> -0.009971   0.9026        0
>>>>  0.4837  -0.3723        0
>>>>
>>>> sigma
>>>>  1.105       0
>>>>     0  0.5873
>>>>     0       0
>>>>
>>>> matrixV.transpose()
>>>>  0.2205  -0.9754
>>>> -0.9754  -0.2205
>>>> This unit test is really very insufficient. It should at least have
>>>> checked that U is unitary.
>>>>
>>>> Benoit
>>>>
>>>>
>>>>>
>>>>> I know that Benoit is still working on rewriting the SVD but at least
>>>>> the bad memory access should be fixed.
>>>>>
>>>>> Cheers,
>>>>> Hauke
>>>>>
>>>>
>>>
>>
>>
>>
>
>
>