|Re: [eigen] Status of givens QR decomposition|
[ Thread Index |
| More lists.tuxfamily.org/eigen Archives
- To: eigen@xxxxxxxxxxxxxxxxxxx
- Subject: Re: [eigen] Status of givens QR decomposition
- From: Benoit Jacob <jacob.benoit.1@xxxxxxxxx>
- Date: Tue, 27 Apr 2010 08:19:34 -0400
- Dkim-signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=gmail.com; s=gamma; h=domainkey-signature:mime-version:received:received:in-reply-to :references:date:message-id:subject:from:to:content-type :content-transfer-encoding; bh=900NdVPCopQoUJDYlX2c1TEuWjPijUJiZ0r8tw4Bj8s=; b=PAxflJMigJaIvA0x1T98ePz8pWCHrwYAglJMhBE/rGgdT0fhYFUYmjrnuUNJIVK8lK tizEM0+A80s7udB9IPgTmOmz2RILgm5R/jgREI8hG9uffpBOqGqAzs/onAbeBVnTTP1z dFj8rVrFwFdEGl6mY3Xb7f//9bdeYSFz12J9E=
- Domainkey-signature: a=rsa-sha1; c=nofws; d=gmail.com; s=gamma; h=mime-version:in-reply-to:references:date:message-id:subject:from:to :content-type:content-transfer-encoding; b=DjntSMq9/zDzjxv9ZbtsmNdQLErb5w4p095hvh2MWRNfycaM7j66VH4qdIt5YXCkgR Da8GqKbGGkWxiPPd8DwehRoE6++QHmhxy3LBjFz2IryzXmAg1SGMPASgwkK3gbxRzmcV KaLZDXYQT23A5/gjIIfsLESoXYqR43zv6P3iQ=
2010/4/27 Benoit Jacob <jacob.benoit.1@xxxxxxxxx>:
> 2010/4/27 Adolfo Rodríguez Tsouroukdissian <dofo79@xxxxxxxxx>:
>> On Tue, Jan 26, 2010 at 1:52 PM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx>
>>> 2010/1/26 Thomas Capricelli <orzel@xxxxxxxxxxxxxxx>:
>>> > As a following to:
>>> > http://listengine.tuxfamily.org/lists.tuxfamily.org/eigen/2009/09/msg00009.html
>>> > I'd like to say that i have a very important use case for givens-QR with
>>> > big matrices in the NonlinearOptimization module. This is useful in the
>>> > variant of the Levenberg Marquardt algorithm for very big matrices,
>>> > so-called 'storage efficient': the qr decomposition is done with givens
>>> > rotation, presumably because this is the only decomposition that can be done
>>> > incrementally : the matrix to decompose is not known as a whole, but is
>>> > constructed line by line, 'updating' the qr decomposition and deleting the
>>> > line.
>>> > Andrea, Benoit, what is the status of this fork ? As I have write
>>> > permission, I've merged with mainline once or twice, but the fork seems
>>> > abandoned... was there a use case ? not anymore ?
>>> I haven't touched it in a long time. The original use case was small
>>> matrices, where it is faster than householder. Now you mention large
>>> matrices that are not fully known in advance: interesting.
>>> Finally, about the old debate about whether to store the list of
>>> givens permutations, let me mention that in the different setting of
>>> SVD, I have now decide to try this approach too. It's different as I'm
>>> only storing N Givens rotations per sweep, as opposed to N^2/2 in
>>> GivensQR, but perhaps still worth mentioning.
>> Let's revive this thread a little bit :)
>> I'm curious to have a general idea of the solution method you want to
>> implement here (SVD based on Givens rotations, storing only N
> There was a "little problem" with the above: Givens SVD actually
> needed N^2 Givens rotations per sweep, not N.
> So instead, I'm now looking at divide-and-conquer SVD.
> This is the "modern" approach anyway. It's parallelizable, which
> Givens SVD isn't.
> In LAPACK terms, this is SGESDD, as opposed to SGESVD.
> The big reason for me to look at divide-and-conquer SVD is that it
> requires very little storage. Here we're starting from a bidiagonal
> matrix. In order to reduce the problem of size N into two sub-problems
> of size N/2, if I understand well, I need to store 2N scalars. So I
> expect the whole thing to use O(N) storage.
Sorry, I meant O(N log N).
> Of course, if one starts
> from a general rectangular matrix, one has to bidiagonalize it first,
> but that's done already in Eigen, see UpperBidiagonalization, using
> householder transforms.
> Divide-and-conquer is only applicable to narrow enough band matrices.
> For example, it can't be used to compute the QR of a general matrix,
> or the eigendecomposition of a Hessenberg matrix.
>> Is there a paper that you could refer me to?
> Yes, this one is very good:
>> One sweet
>> thing of using Givens rotations (as opposed to say Householder reflections)
>> is that their computation can be easily parallelized across multiple cores.
> ...depends on the order in which the Givens rotations are applied,
> since parallelization is only possible for mutually commuting (AB=BA)
> transformation --- otherwise you get different results depending on
> which one is applied first.
> For example, Givens SVD is not parallelizable at all.
> Givens QR, yes, it's parallelizable.
>> On a separate note, I'd like to ask if there would be interest in having a
>> SVD implementation tailored for the case of using one decomposition instance
>> to solve multiple "similar" problems. By similar I'm referring to solve
>> initially for matrix A, and then for a perturbed matrix A + deltaA, and so
>> on. In such cases, results from previous computations can greatly increase
>> performance. The KDL developers  brought to my attention a paper that
>> details this solution strategy (also based on Givens rotations), which can
>> be found here . They have also a working implementation  that uses
>> eigen matrices, but its internals are not "eigenized".
> That sounds interesting for certain classes of problems. Though it can
> only work if certain assumptions on the matrix A are made. At least,
> it should sufficiently many distinct singular values. Otherwise, bad
> phenomena can happen. For example, let A = Id be the identity matrix.
> Then a SVD of A is just: A = Id x Id x Id. Now take a small
> perturbation of A: it could have absolutely any singular vectors on
> both sides --- just because A is close to Id doesn't imply anything on
> its singular vectors. So the knowledge that A is close to Identity,
> doesn't help at all.
> But yes, that sounds interesting in principle, we'll see about that post-3.0.
>> Finally, I wanted to ask the more general question of what SVD
>> implementations will likely make it to 3.0.
> I am currently writing the new divide-and-conquer SVD. I would like it
> to replace the current SVD, but of course i'll only do that once it's
> clear that it's at least at good.
> I also want to keep JacobiSVD, as it is the only SVD avoiding
> altogether householder bidiagonalization, which is not perfectly
> numerically stable.
> I know I have been talking about that forever, but now I'm really
> doing it, and i'm starting at mozilla in 2 weeks anyway so it must be
> done before.
>>  http://orocos.org/kdl
>>  http://www.engr.colostate.edu/~aam/pdf/journals/5.pdf
>>> > btw, Am I the only one needing givens QR decompozition on the list ?
>>> > anyone ?
>>> It is for sure needed at least by everyone needing QR decomposition of
>>> very small matrices.
>>> If you need to make intrusive changes in Andrea's fork, perhaps it's
>>> easiest for you to fork it, you can then re-merge with Andrea if he
>>> expresses agreement with your changes. I personally don't have time
>>> right now to think about that... :(
>> Adolfo Rodríguez Tsouroukdissian, Ph. D.
>> Robotics engineer
>> PAL ROBOTICS S.L
>> Tel. +34.93.414.53.47
>> AVISO DE CONFIDENCIALIDAD: Este mensaje y sus documentos adjuntos, pueden
>> contener información privilegiada y/o confidencial que está dirigida
>> exclusivamente a su destinatario. Si usted recibe este mensaje y no es el
>> destinatario indicado, o el empleado encargado de su entrega a dicha
>> persona, por favor, notifíquelo inmediatamente y remita el mensaje original
>> a la dirección de correo electrónico indicada. Cualquier copia, uso o
>> distribución no autorizados de esta comunicación queda estrictamente
>> CONFIDENTIALITY NOTICE: This e-mail and the accompanying document(s) may
>> contain confidential information which is privileged and intended only for
>> the individual or entity to whom they are addressed. If you are not the
>> intended recipient, you are hereby notified that any disclosure, copying,
>> distribution or use of this e-mail and/or accompanying document(s) is
>> strictly prohibited. If you have received this e-mail in error, please
>> immediately notify the sender at the above e-mail address.