One thing I would like to suggest.

We can provide multiple ways to compute the sparse LU
1. scipy 
2. sage original implementation in src.sage.matrix.matrix2.LU 
<https://github.com/sagemath/sage/blob/acbe15dcd87085d4fa177705ec01107b53a026ef/src/sage/matrix/matrix2.pyx#L13160>
 (Note 
- link 
<https://github.com/sagemath/sage/blob/acbe15dcd87085d4fa177705ec01107b53a026ef/src/sage/matrix/matrix2.pyx#L13249C1-L13254C51>
)
3. convert to dense then factor

It will be up to user to choose based on the complexity.
Is it fine?

On Wednesday, February 28, 2024 at 4:30:51 PM UTC+5:30 Animesh Shree wrote:

> Thank you for reminding
> I went through. 
> We need to Decompose  A11 only and rest can be calculated via taking 
> inverse of L11 or U11.
> Here A11 is square matrix and we can use scipy to decompose square 
> matrices.
> Am I correct?
>
> New and only problem that I see is the returned LU decomposition of 
> scipy's splu is calculated by full permutation of row and column as pointed 
> out by *Nils Bruin*. We will be returning row and col permutation 
> array/matrix separately instead of single row permutation which sage 
> usage generally for plu decomposition.
> User will have to manage row and col permutations. 
> or else
> We can return handler function for reconstruction of matrix from  L, U & 
> p={perm_r, 
> perm_c}
> or
> We can leave that to user
> User will have to permute its input data according to perm_c (like : 
> perm_c * input) before using the perm_r^(-1) * L * U
> as perm_r^(-1) * L * U is PLU decomposition of Original_matrix*perm_c
>
> https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.SuperLU.html
> >>> A = Pr^(-1) *L*U * Pc^(-1) # as told by *Nils Bruin*
> or
> scipy's splu will not do.
>
> On Tuesday, February 27, 2024 at 11:57:02 PM UTC+5:30 Dima Pasechnik wrote:
>
>>
>>
>> On 27 February 2024 17:25:51 GMT, 'Animesh Shree' via sage-devel <
>> sage-...@googlegroups.com> wrote: 
>> >This works good if input is square and I also checked on your idea of 
>> >padding zeros for non square matrices. 
>> >I am currently concerned about the permutation matrix and L, U in case 
>> of 
>> >padded 0s. Because if we pad then how will they affect the outputs, so 
>> that 
>> >we can extract p,l,u for unpadded matrix. 
>>
>> please read details I wrote on how to deal with the non-square case. 
>> There is no padding needed. 
>>
>>
>> > 
>> >On Tuesday, February 27, 2024 at 10:03:25 PM UTC+5:30 Dima Pasechnik 
>> wrote: 
>> > 
>> >> 
>> >> 
>> >> On 27 February 2024 15:34:20 GMT, 'Animesh Shree' via sage-devel < 
>> >> sage-...@googlegroups.com> wrote: 
>> >> >I tried scipy which uses superLU. We get the result but there is 
>> little 
>> >> bit 
>> >> >of issue. 
>> >> > 
>> >> > 
>> >> >--For Dense-- 
>> >> >The dense matrix factorization gives this output using permutation 
>> matrix 
>> >> >sage: a = Matrix(RDF, [[1, 0],[2, 1]], sparse=True) 
>> >> >sage: a 
>> >> >[1.0 0.0] 
>> >> >[2.0 1.0] 
>> >> >sage: p,l,u = a.dense_matrix().LU() 
>> >> >sage: p 
>> >> >[0.0 1.0] 
>> >> >[1.0 0.0] 
>> >> >sage: l 
>> >> >[1.0 0.0] 
>> >> >[0.5 1.0] 
>> >> >sage: u 
>> >> >[ 2.0 1.0] 
>> >> >[ 0.0 -0.5] 
>> >> > 
>> >> 
>> >> you'd probably want to convert the permutation matrix into a 
>> permutation. 
>> >> 
>> >> 
>> >> >--For Sparse-- 
>> >> >But the scipy LU decomposition uses permutations which involves 
>> taking 
>> >> >transpose, also the output permutations are represented as array. 
>> >> 
>> >> It is very normal to represent permutations as arrays. 
>> >> One can reconstruct the permutation matrix from such an array 
>> trivially 
>> >> (IIRC, Sage even has a function for it) 
>> >> 
>> >> I am not sure what you mean by "taking transpose". 
>> >> 
>> >> >sage: p,l,u = a.LU(force=True) 
>> >> >sage: p 
>> >> >{'perm_r': [1, 0], 'perm_c': [1, 0]} 
>> >> >sage: l 
>> >> >[1.0 0.0] 
>> >> >[0.0 1.0] 
>> >> >sage: u 
>> >> >[1.0 2.0] 
>> >> >[0.0 1.0] 
>> >> > 
>> >> > 
>> >> >Shall I continue with this? 
>> >> 
>> >> sure, you are quite close to getting it all done it seems. 
>> >> 
>> >> 
>> >> >On Tuesday, February 6, 2024 at 11:29:07 PM UTC+5:30 Dima Pasechnik 
>> wrote: 
>> >> > 
>> >> >> Non-square case for LU is in fact easy. Note that if you have A=LU 
>> as 
>> >> >> a block matrix 
>> >> >> A11 A12 
>> >> >> A21 A22 
>> >> >> 
>> >> >> then its LU-factors L and U are 
>> >> >> L11 0 and U11 U12 
>> >> >> L21 L22 0 U22 
>> >> >> 
>> >> >> and A11=L11 U11, A12=L11 U12, A21=L21 U11, A22=L21 U12+L22 U22 
>> >> >> 
>> >> >> Assume that A11 is square and full rank (else one may apply 
>> >> >> permutations of rows and columns in the usual way). while A21=0 and 
>> >> >> A22=0. Then one can take L21=0, L22=U22=0, while A12=L11 U12 
>> >> >> implies U12=L11^-1 A12. 
>> >> >> That is, we can first compute LU-decomposition of a square matrix 
>> A11, 
>> >> >> and then compute U12 from it and A. 
>> >> >> 
>> >> >> Similarly, if instead A12=0 and A22=0, then we can take U12=0, 
>> >> >> L22=U22=0, and A21=L21 U11, 
>> >> >> i.e. L21=A21 U11^-1, and again we compute LU-decomposition of A11, 
>> and 
>> >> >> then L21=A21 U11^-1. 
>> >> >> 
>> >> >> ---------------- 
>> >> >> 
>> >> >> Note that in some cases one cannot get LU, but instead must go for 
>> an 
>> >> >> PLU,with P a permutation matrix. 
>> >> >> For non-square matrices this seems a bit more complicated, but, 
>> well, 
>> >> >> still doable. 
>> >> >> 
>> >> >> HTH 
>> >> >> Dima 
>> >> >> 
>> >> >> 
>> >> >> 
>> >> >> 
>> >> >> On Mon, Feb 5, 2024 at 6:00 PM Nils Bruin <nbr...@sfu.ca> wrote: 
>> >> >> > 
>> >> >> > On Monday 5 February 2024 at 02:31:04 UTC-8 Dima Pasechnik wrote: 
>> >> >> > 
>> >> >> > 
>> >> >> > it is the matter of adding extra zero rows or columns to the 
>> matrix 
>> >> you 
>> >> >> want to decompose. This could be a quick fix. 
>> >> >> > 
>> >> >> > (in reference to computing LU decompositions of non-square 
>> matrices) 
>> >> -- 
>> >> >> in a numerical setting, adding extra zero rows/columns may not be 
>> such 
>> >> an 
>> >> >> attractive option: if previously you know you had a maximal rank 
>> >> matrix, 
>> >> >> you have now ruined it by the padding. It's worth checking the 
>> >> >> documentation and literature if padding is appropriate/desirable 
>> for 
>> >> the 
>> >> >> target algorithm/implementation. 
>> >> >> > 
>> >> >> > -- 
>> >> >> > You received this message because you are subscribed to the 
>> Google 
>> >> >> Groups "sage-devel" group. 
>> >> >> > To unsubscribe from this group and stop receiving emails from it, 
>> >> send 
>> >> >> an email to sage-devel+...@googlegroups.com. 
>> >> >> > To view this discussion on the web visit 
>> >> >> 
>> >> 
>> https://groups.google.com/d/msgid/sage-devel/622a01e0-9197-40c5-beda-92729c4e4a32n%40googlegroups.com
>>  
>> >> >> . 
>> >> >> 
>> >> > 
>> >> 
>> > 
>>
>

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/9aa9deef-5fd9-4515-ae73-6e92213dcb59n%40googlegroups.com.

Reply via email to