Re: [eigen] fast calculation of matrix rank |

[ Thread Index | Date Index | More lists.tuxfamily.org/eigen Archives ]

*To*: eigen@xxxxxxxxxxxxxxxxxxx*Subject*: Re: [eigen] fast calculation of matrix rank*From*: Pieter Eendebak <pieter.eendebak@xxxxxxxxx>*Date*: Tue, 16 Feb 2016 00:55:34 +0100*Dkim-signature*: v=1; a=rsa-sha256; c=relaxed/relaxed; d=gmail.com; s=20120113; h=mime-version:in-reply-to:references:date:message-id:subject:from:to :content-type; bh=xN02P2xyFspNDnQ/GgI06rklH/U2weXYWKN7QoP76TM=; b=Yxc5cMg9iPPBad18xmRV0GFfyU7gIbOxof1BwVisOLo9iChOBPyUcrczK9fG8WYxNL TfuAwEeWz0AUiYxUqPFz7gFSK2ts9fw1D9Dg2EHot8TnbSaUwntXjoqlgLokfZn2nq2j kqlrnAueQoWYcV/wGtcSX+99qzbu6d9pREeANY3H99JW1DKaln2ypAkKLhTN0OvIa4cq 6UXkKottfGTmwpfb+DHgE9v3aT37GEy4TA4++HfusA7bxqS7Dx4nbWxxe5lkMfI4yjTV NUkayCV1vNp5LmnOQDsGaxIqKNuYxHFH0UzpxqNpoBPjCIlSOZj/Pw3NTI6qDe3dqBvj 1ihQ==

Hi Marco,

Thanks for the reply. A_i is the partioned matrix [A | B_i]. I have succeeded in using Eigen to make my rank calculations faster (factor of 2 to3), which was worth the effort for me.

You are right that in general rank(AB_i) <= min(rank(A), rank(B_i)). In my case however, I have used the transformation with Q and P such that A is orthogonal to

[ 0 ]

[ Y_i ]

Then rank(X_i) = rank( [A [ 0; Y_i ] ]) = rank( A) + rank [ 0; Y_i ] ) = rank( A) + rank [Y_i ] )

Taking the transpose of Q and P is indeed faster than the inverse. The reason for taking Q^{-1} and not Q^T in my code is that I only calculate Q^{-1} once (so computation overhead is small) and I was worried that Q^T might lead to larger rounding errors than Q^{-1}. I am not sure whether these rounding errors can occur (for P they cannot since it has integer entries).

With kind regards,

Pieter Eendebak

On Sat, Jan 30, 2016 at 2:39 AM, marco stamazza <giunghi@xxxxxxxxx> wrote:

Pieter,probably you are better off asking this question at http://math.stackexchange..com/questions/tagged/linear-algebraYou mean that A_i is a partitioned matrix [A | B_i] or a product A * B_i?what is the size of the matrices A and B_i?Note that the rank of product and partitioned matrices is rather more complicated than the sum of the ranks you use (http://link.springer.com/chapter/10.1007%2F978-3-642-10473-2_6#page-1)for example rank(AB_i) <= min(rank(A), rank(B_i)), not sure where the additive property you state comes from.In any case, I don't understand why you are computing Q^(-1) A_i P^(-1) i, since Q iand (thin) P are orthogonal only need take the transpose. Furthermore, since AP = QR, Q^(-1) A_i P^(-1) = Q^(-1) [A B_i] P^(-1) = RB_iP' (or similar for partitioned matrix).. So the transformation doesn't help you, on the contrary, R is already available with its rank, so you still need determine the rank of (B_i)..If A_i = [A | B_i] then you could compute the rank of A and the the ranks of only the B_i's with colPivQR. but you still have to deal with the formula in the reference aboveIf A_i is a product A*B_i, a strategy could be computing If A_iA_i' and then get the rank by computing only its eigenvalues or with LU decomposition, which could be faster then colPivQR, but I am not sure.In my opinion with such small nrow and ncol there is no point in inventing algorithms to save few flops, unless you have some properties of the A and B_i matrices you can exploit. .I hope this helps, Marco--On 30 January 2016 at 05:24, Pieter Eendebak <pieter.eendebak@xxxxxxxxx> wrote:Hi,Using Eigen I am trying to calculate the matrix rank of a large number of arrays. So far I have used FullPivLU but this is too slow for my application. The alternative ColPivHouseholderQR is faster, but still not fast enough. What are options for calculating the rank faster? My matrices have size 64x92 and integer elements (+1 or -1).There is one additional fact that might used. All matrices are of the formA_i = [A B_i]with A a common submatrix of the A_i.I can then calculate the rank of A_i in the following way. First perform ColPivHouseholderQR on A:AP = QRand determine the matrix rank of A as r0. Then for each matrix A_i we calculateX_i = Q^(-1) A_i P^(-1)The rank of X_i is equal to the rank of A_i. Next, take the submatrixY_i = X_i[r0:, r0:]Then rank(A_i) = rank(X_i) = rank(A) + rank(Y_i).Instead of calculating the rank of A_i we only have to calculate the rank of Y_i, which is much smaller in size than A_i. I have checked that calculation of the ColPivHouseholderQR decomposition of the Y_i matrices is indeed faster. However, the multiplication Q^(-1) A_i P^(-1) is almost as expensive as the decomposition on the original array X_i. Is there a way I can use the common submatrix in the rank calculations without performing the multiplication?With kind regards,PieterThe matrices share a common block Most of matrices"Anyone who lives within their means suffers from a lack of imagination."

**Messages sorted by:**[ date | thread ]- Prev by Date:
**[eigen] tensor fft** - Next by Date:
**Re: [eigen] tensor fft** - Previous by thread:
**Re: [eigen] tensor fft** - Next by thread:
**[eigen] unsubscribe**

Mail converted by MHonArc 2.6.19+ | http://listengine.tuxfamily.org/ |